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

Merge pull request #1770 from bska/add-tolcrit

Browse Source 
Add Support for TOLCRIT to Sat-Func Field Properties
This commit is contained in:
Bård Skaflestad
2020-06-10 19:19:07 +02:00
committed by GitHub
parent 6a742a20fd f61882ce69
commit 9396d44379
6 changed files with 820 additions and 231 deletions
Download Patch File Download Diff File Expand all files Collapse all files

27
opm/parser/eclipse/EclipseState/Grid/SatfuncPropertyInitializers.hpp
View File

@@ -19,6 +19,7 @@
#ifndef ECLIPSE_SATFUNCPROPERTY_INITIALIZERS_HPP
#define ECLIPSE_SATFUNCPROPERTY_INITIALIZERS_HPP
#include <memory>
#include <string>
#include <vector>
@@ -29,9 +30,35 @@ namespace Opm {
namespace Opm { namespace satfunc {
struct RawTableEndPoints
{
struct {
std::vector<double> gas;
std::vector<double> water;
} connate;
struct {
std::vector<double> oil_in_gas;
std::vector<double> oil_in_water;
std::vector<double> gas;
std::vector<double> water;
} critical;
struct {
std::vector<double> gas;
std::vector<double> water;
} maximum;
};
std::shared_ptr<RawTableEndPoints>
getRawTableEndpoints(const Opm::TableManager& tm,
const Opm::Phases& phases,
const double tolcrit);
std::vector<double> init(const std::string& kewyord,
const TableManager& tables,
const Phases& phases,
const RawTableEndPoints& ep,
const std::vector<double>& cell_depth,
const std::vector<int>& num,
const std::vector<int>& endnum);

2
src/opm/output/eclipse/Tables.cpp
View File

@@ -232,7 +232,7 @@ namespace { namespace SatFunc {
std::transform(begin, end, dest,
[tolcrit](const double kr) -> double
{
return (kr < tolcrit) ? 0.0 : kr;
return (kr > tolcrit) ? kr : 0.0;
});
}

17
src/opm/parser/eclipse/EclipseState/Grid/FieldProps.cpp
View File

@@ -409,6 +409,7 @@ FieldProps::FieldProps(const Deck& deck, const Phases& phases, const EclipseGrid
ny(grid.getNY()),
nz(grid.getNZ()),
m_phases(phases),
m_satfuncctrl(deck),
m_actnum(grid.getACTNUM()),
cell_volume(extract_cell_volume(grid)),
cell_depth(extract_cell_depth(grid)),
@@ -1056,14 +1057,16 @@ void FieldProps::scanEDITSection(const EDITSection& edit_section) {
void FieldProps::init_satfunc(const std::string& keyword, FieldData<double>& satfunc) {
if (this->m_rtep == nullptr)
this->m_rtep = satfunc::getRawTableEndpoints(this->tables, this->m_phases,
this->m_satfuncctrl.minimumRelpermMobilityThreshold());
const auto& endnum = this->get<int>("ENDNUM");
if (keyword[0] == 'I') {
const auto& imbnum = this->get<int>("IMBNUM");
satfunc.default_update(satfunc::init(keyword, this->tables, this->m_phases, this->cell_depth, imbnum, endnum));
} else {
const auto& satnum = this->get<int>("SATNUM");
satfunc.default_update(satfunc::init(keyword, this->tables, this->m_phases, this->cell_depth, satnum, endnum));
}
const auto& satreg = (keyword[0] == 'I')
? this->get<int>("IMBNUM")
: this->get<int>("SATNUM");
satfunc.default_update(satfunc::init(keyword, this->tables, this->m_phases, *this->m_rtep, this->cell_depth, satreg, endnum));
}

4
src/opm/parser/eclipse/EclipseState/Grid/FieldProps.hpp
View File

@@ -19,6 +19,7 @@
#ifndef FIELDPROPS_HPP
#define FIELDPROPS_HPP
#include <memory>
#include <string>
#include <unordered_set>
#include <vector>
@@ -27,6 +28,7 @@
#include <opm/parser/eclipse/Deck/DeckSection.hpp>
#include <opm/parser/eclipse/Units/UnitSystem.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/Box.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/SatfuncPropertyInitializers.hpp>
#include <opm/parser/eclipse/EclipseState/Runspec.hpp>
namespace Opm {
@@ -357,12 +359,14 @@ private:
const UnitSystem unit_system;
std::size_t nx,ny,nz;
Phases m_phases;
SatFuncControls m_satfuncctrl;
std::vector<int> m_actnum;
std::vector<double> cell_volume;
std::vector<double> cell_depth;
const std::string m_default_region;
const EclipseGrid * grid_ptr; // A bit undecided whether to properly use the grid or not ...
const TableManager& tables;
std::shared_ptr<satfunc::RawTableEndPoints> m_rtep;
std::vector<MultregpRecord> multregp;
std::unordered_map<std::string, FieldData<int>> int_data;
std::unordered_map<std::string, FieldData<double>> double_data;

383
src/opm/parser/eclipse/EclipseState/Grid/SatfuncPropertyInitializers.cpp
View File

@@ -27,6 +27,7 @@
#include <opm/parser/eclipse/EclipseState/Tables/SwfnTable.hpp>
#include <opm/parser/eclipse/EclipseState/Tables/SwofTable.hpp>
#include <opm/parser/eclipse/EclipseState/Tables/Tabdims.hpp>
#include <opm/parser/eclipse/EclipseState/Tables/TableColumn.hpp>
#include <opm/parser/eclipse/EclipseState/Tables/TableContainer.hpp>
#include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
@@ -34,34 +35,32 @@
#include <algorithm>
#include <array>
#include <cassert>
#include <exception>
#include <functional>
#include <iterator>
#include <memory>
#include <stdexcept>
#include <string>
#include <type_traits>
#include <utility>
#include <stddef.h>
// Note on deriving critical saturations: All table scanners are implemented
// in terms of std::lower_bound(begin, end, tolcrit, predicate) which returns
// the first position in [begin, end) for which
//
// predicate(*iter, tolcrit)
//
// is false. Using predicate = std::greater<>{} thus determines the first
// position in the sequence for which the elements is less than or equal to
// 'tolcrit'. Similarly, a predicate equivalent to '<=' returns the first
// position for which the elements is strictly greater than 'tolcrit'.
namespace {
struct RawTableEndPoints
{
struct {
std::vector<double> gas;
std::vector<double> water;
} connate;
struct {
std::vector<double> oil_in_gas;
std::vector<double> oil_in_water;
std::vector<double> gas;
std::vector<double> water;
} critical;
struct {
std::vector<double> gas;
std::vector<double> water;
} maximum;
};
using ::Opm::satfunc::RawTableEndPoints;
/*
* See the "Saturation Functions" chapter in the Eclipse Technical
@@ -116,8 +115,6 @@ namespace {
return SatfuncFamily::none;
}
// enum class limit { min, max };
std::vector<double>
findMinWaterSaturation(const Opm::TableManager& tm,
const Opm::Phases& ph)
@@ -260,39 +257,67 @@ namespace {
}
}
/*
* These functions have been ported from an older implementation to instead
* use std::upper_bound and more from <algorithm> to make code -intent-
* clearer. This also made some (maybe intentional) details easier to spot.
* A short discussion:
*
* I don't know if not finding any element larger than 0.0 in the tables
* was ever supposed to happen (or even possible), in which case the vector
* elements remained at their initial value of 0.0. This behaviour has been
* preserved, but is now explicit. The original code was also not clear if
* it was possible to look up columns at index -1 (see critical_water for
* an example), but the new version is explicit about this. Unfortuately
* I'm not familiar enough with the maths or internal structure to make
* more than a guess here, but most of this behaviour should be preserved.
*
*/
template <typename T>
double critical_water( const T& table )
template <typename Predicate>
auto crit_sat_index(const Opm::TableColumn& col,
const double tolcrit,
Predicate&& pred)
{
const auto& col = table.getKrwColumn();
const auto end = col.begin() + table.numRows();
const auto critical = std::upper_bound( col.begin(), end, 0.0 );
const auto index = std::distance( col.begin(), critical );
using SizeT = std::remove_const_t<
std::remove_reference_t<decltype(col.size())>
>;
if( index == 0 || critical == end ) return 0.0;
auto begin = col.begin();
auto pos = std::lower_bound(begin, col.end(), tolcrit,
std::forward<Predicate>(pred));
return table.getSwColumn()[ index - 1 ];
assert ((pos != col.end()) &&
"Detected relative permeability function "
"without immobile state");
return static_cast<SizeT>(std::distance(begin, pos));
}
std::vector< double >
double crit_sat_increasing_KR(const Opm::TableColumn& sat,
const Opm::TableColumn& kr,
const double tolcrit)
{
// First position for which Kr(S) > tolcrit.
const auto i = crit_sat_index(kr, tolcrit,
[](const double kr1, const double kr2)
{
// kr1 <= kr2. Kr2 is 'tolcrit'.
return ! (kr2 < kr1);
});
return sat[i - 1]; // Last saturation for which Kr(S) <= tolcrit
}
double crit_sat_decreasing_KR(const Opm::TableColumn& sat,
const Opm::TableColumn& kr,
const double tolcrit)
{
// First position for which Kr(S) <= tolcrit.
const auto i = crit_sat_index(kr, tolcrit, std::greater<>{});
return sat[i];
}
/// Maximum water saturation for which Krw(Sw) <= tolcrit.
///
/// Expected Table Format:
/// [Sw, Krw(Sw), ...other...]
///
/// Krw increasing.
template <typename T>
double critical_water(const T& table, const double tolcrit)
{
return crit_sat_increasing_KR(table.getSwColumn(),
table.getKrwColumn(), tolcrit);
}
std::vector<double>
findCriticalWater(const Opm::TableManager& tm,
const Opm::Phases& ph)
const Opm::Phases& ph,
const double tolcrit)
{
const auto num_tables = tm.getTabdims().getNumSatTables();
@@ -302,12 +327,14 @@ namespace {
const auto& swofTables = tm.getSwofTables();
const auto& swfnTables = tm.getSwfnTables();
const auto famI = [&swofTables]( int i ) {
return critical_water( swofTables.getTable<Opm::SwofTable>( i ) );
const auto famI = [&swofTables, tolcrit](const int i) -> double
{
return critical_water(swofTables.getTable<Opm::SwofTable>(i), tolcrit);
};
const auto famII = [&swfnTables]( int i ) {
return critical_water( swfnTables.getTable<Opm::SwfnTable>( i ) );
const auto famII = [&swfnTables, tolcrit](const int i) -> double
{
return critical_water(swfnTables.getTable<Opm::SwfnTable>(i), tolcrit);
};
switch( getSaturationFunctionFamily( tm, ph ) ) {
@@ -317,31 +344,40 @@ namespace {
}
}
template< typename T >
double critical_gas( const T& table ) {
const auto& col = table.getKrgColumn();
const auto end = col.begin() + table.numRows();
const auto critical = std::upper_bound( col.begin(), end, 0.0 );
const auto index = std::distance( col.begin(), critical );
if( index == 0 || critical == end ) return 0.0;
return table.getSgColumn()[ index - 1 ];
/// Maximum gas saturation for which Krg(Sg) <= tolcrit.
///
/// Expected Table Format:
/// [Sg, Krg(Sg), ...other...]
///
/// Krg increasing.
template <typename T>
double critical_gas(const T& table, const double tolcrit)
{
return crit_sat_increasing_KR(table.getSgColumn(),
table.getKrgColumn(), tolcrit);
}
double critical_gas( const Opm::SlgofTable& slgofTable ) {
const auto& col = slgofTable.getKrgColumn();
const auto critical = std::upper_bound( col.begin(), col.end(), 0.0 );
const auto index = std::distance( col.begin(), critical );
/// Maximum gas saturation for which Krg(Sg) <= tolcrit.
///
/// Table Format (Sl = So + Swco):
/// [Sl, Krg(Sl), Krog(Sl), Pcgo(Sl)]
///
/// Krg decreasing, Krog increasing, Pcog not increasing.
double critical_gas(const Opm::SlgofTable& slgofTable,
const double tolcrit)
{
const auto sl_at_crit_gas =
crit_sat_decreasing_KR(slgofTable.getSlColumn(),
slgofTable.getKrgColumn(), tolcrit);
if( index == 0 || critical == col.end() ) return 0.0;
return slgofTable.getSlColumn()[ index - 1 ];
// Sg = 1 - Sl
return 1.0 - sl_at_crit_gas;
}
std::vector<double>
findCriticalGas(const Opm::TableManager& tm,
const Opm::Phases& ph)
const Opm::Phases& ph,
const double tolcrit)
{
const auto num_tables = tm.getTabdims().getNumSatTables();
@@ -352,16 +388,19 @@ namespace {
const auto& sgofTables = tm.getSgofTables();
const auto& slgofTables = tm.getSlgofTables();
const auto famI_sgof = [&sgofTables]( int i ) {
return critical_gas( sgofTables.getTable<Opm::SgofTable>( i ) );
const auto famI_sgof = [&sgofTables, tolcrit](const int i) -> double
{
return critical_gas(sgofTables.getTable<Opm::SgofTable>(i), tolcrit);
};
const auto famI_slgof = [&slgofTables]( int i ) {
return critical_gas( slgofTables.getTable<Opm::SlgofTable>( i ) );
const auto famI_slgof = [&slgofTables, tolcrit](const int i) -> double
{
return critical_gas(slgofTables.getTable<Opm::SlgofTable>(i), tolcrit);
};
const auto famII = [&sgfnTables]( int i ) {
return critical_gas( sgfnTables.getTable<Opm::SgfnTable>( i ) );
const auto famII = [&sgfnTables, tolcrit](const int i) -> double
{
return critical_gas(sgfnTables.getTable<Opm::SgfnTable>(i), tolcrit);
};
switch( getSaturationFunctionFamily( tm, ph ) ) {
@@ -382,42 +421,53 @@ namespace {
}
}
double critical_oil_water( const Opm::SwofTable& swofTable ) {
const auto& col = swofTable.getKrowColumn();
/// Maximum oil saturation for which Krow(So) <= tolcrit.
///
/// Table Format:
/// [Sw, Krw(Sw), Krow(Sw), Pcow(Sw)]
///
/// Krw increasing, Krow decreasing, Pcow not increasing.
double critical_oil_water(const Opm::SwofTable& swofTable,
const double tolcrit)
{
const auto sw_at_crit_oil =
crit_sat_decreasing_KR(swofTable.getSwColumn(),
swofTable.getKrowColumn(), tolcrit);
using reverse = std::reverse_iterator< decltype( col.begin() ) >;
auto rbegin = reverse( col.begin() + swofTable.numRows() );
auto rend = reverse( col.begin() );
const auto critical = std::upper_bound( rbegin, rend, 0.0 );
const auto index = std::distance( col.begin(), critical.base() - 1 );
if( critical == rend ) return 0.0;
return 1 - swofTable.getSwColumn()[ index + 1 ];
// So = 1 - Sw
return 1.0 - sw_at_crit_oil;
}
double critical_oil( const Opm::Sof2Table& sof2Table ) {
const auto& col = sof2Table.getKroColumn();
const auto critical = std::upper_bound( col.begin(), col.end(), 0.0 );
const auto index = std::distance( col.begin(), critical );
if( index == 0 || critical == col.end() ) return 0.0;
return sof2Table.getSoColumn()[ index - 1 ];
/// Maximum oil saturation for which Kro(So) <= tolcrit.
///
/// Table Format:
/// [So, Kro(So)]
///
/// Kro increasing.
double critical_oil(const Opm::Sof2Table& sof2Table,
const double tolcrit)
{
return crit_sat_increasing_KR(sof2Table.getSoColumn(),
sof2Table.getKroColumn(), tolcrit);
}
double critical_oil( const Opm::Sof3Table& sof3Table, const Opm::TableColumn& col ) {
const auto critical = std::upper_bound( col.begin(), col.end(), 0.0 );
const auto index = std::distance( col.begin(), critical );
if( index == 0 || critical == col.end() ) return 0.0;
return sof3Table.getSoColumn()[ index - 1 ];
/// Maximum oil saturation for which Kro(So) <= tolcrit.
///
/// Table Format:
/// [So, Krow(So), Krog(So)]
///
/// Krow increasing, Krog increasing.
double critical_oil(const Opm::Sof3Table& sof3Table,
const Opm::TableColumn& col,
const double tolcrit)
{
return crit_sat_increasing_KR(sof3Table.getSoColumn(), col, tolcrit);
}
std::vector<double>
findCriticalOilWater(const Opm::TableManager& tm,
const Opm::Phases& ph)
const Opm::Phases& ph,
const double tolcrit)
{
const auto num_tables = tm.getTabdims().getNumSatTables();
@@ -429,17 +479,20 @@ namespace {
const auto& sof2Tables = tm.getSof2Tables();
const auto& sof3Tables = tm.getSof3Tables();
const auto famI = [&swofTables]( int i ) {
return critical_oil_water( swofTables.getTable<Opm::SwofTable>( i ) );
const auto famI = [&swofTables, tolcrit](const int i) -> double
{
return critical_oil_water(swofTables.getTable<Opm::SwofTable>(i), tolcrit);
};
const auto famII_2p = [&sof2Tables]( int i ) {
return critical_oil( sof2Tables.getTable<Opm::Sof2Table>( i ) );
const auto famII_2p = [&sof2Tables, tolcrit](const int i) -> double
{
return critical_oil(sof2Tables.getTable<Opm::Sof2Table>(i), tolcrit);
};
const auto famII_3p = [&sof3Tables]( int i ) {
const auto& tb = sof3Tables.getTable<Opm::Sof3Table>( i );
return critical_oil( tb, tb.getKrowColumn() );
const auto famII_3p = [&sof3Tables, tolcrit](const int i) -> double
{
const auto& tb = sof3Tables.getTable<Opm::Sof3Table>(i);
return critical_oil(tb, tb.getKrowColumn(), tolcrit);
};
switch( getSaturationFunctionFamily( tm, ph ) ) {
@@ -453,36 +506,41 @@ namespace {
}
}
double critical_oil_gas( const Opm::SgofTable& sgofTable )
/// Maximum oil saturation for which Krog(So) <= tolcrit.
///
/// Table Format:
/// [Sg, Krg(Sg), Krog(Sg), Pcgo(Sg)]
///
/// Krg increasing, Krog decreasing, Pcgo not decreasing.
double critical_oil_gas(const Opm::SgofTable& sgofTable,
const double tolcrit)
{
const auto& col = sgofTable.getKrogColumn();
const auto sg_at_crit_oil =
crit_sat_decreasing_KR(sgofTable.getSgColumn(),
sgofTable.getKrogColumn(), tolcrit);
using reverse = std::reverse_iterator< decltype( col.begin() ) >;
auto rbegin = reverse( col.begin() + sgofTable.numRows() );
auto rend = reverse( col.begin() );
const auto critical = std::upper_bound( rbegin, rend, 0.0 );
if( critical == rend ) {
return 0.0;
}
const auto index = std::distance( col.begin(), critical.base() - 1 );
return 1.0 - sgofTable.getSgColumn()[ index + 1 ];
// So = 1 - Sg
return 1.0 - sg_at_crit_oil;
}
double critical_oil_gas( const Opm::SlgofTable& sgofTable )
/// Maximum oil saturation for which Krog(So) <= tolcrit.
///
/// Table Format (Sl = So + Swco):
/// [Sl, Krg(Sl), Krog(Sl), Pcgo(Sl)]
///
/// Krg decreasing, Krog increasing, Pcgo not increasing.
double critical_oil_gas(const Opm::SlgofTable& slgofTable,
const double tolcrit)
{
const auto& col = sgofTable.getKrogColumn();
const auto critical = std::upper_bound( col.begin(), col.end(), 0.0 );
if (critical == col.end()) {
return 0.0;
}
const auto index = std::distance( col.begin(), critical - 1);
return sgofTable.getSlColumn()[ index ];
return crit_sat_increasing_KR(slgofTable.getSlColumn(),
slgofTable.getKrogColumn(), tolcrit);
}
std::vector<double>
findCriticalOilGas(const Opm::TableManager& tm,
const Opm::Phases& ph,
const std::vector<double>& swco)
const std::vector<double>& swco,
const double tolcrit)
{
const auto num_tables = tm.getTabdims().getNumSatTables();
@@ -495,23 +553,25 @@ namespace {
const auto& sof2Tables = tm.getSof2Tables();
const auto& sof3Tables = tm.getSof3Tables();
const auto famI_sgof = [&sgofTables, &swco](const int i) -> double
const auto famI_sgof = [&sgofTables, &swco, tolcrit](const int i) -> double
{
return critical_oil_gas(sgofTables.getTable<Opm::SgofTable>(i)) - swco[i];
return critical_oil_gas(sgofTables.getTable<Opm::SgofTable>(i), tolcrit) - swco[i];
};
const auto famI_slgof = [&slgofTables, &swco](const int i) -> double
const auto famI_slgof = [&slgofTables, &swco, tolcrit](const int i) -> double
{
return critical_oil_gas(slgofTables.getTable<Opm::SlgofTable>(i)) - swco[i];
return critical_oil_gas(slgofTables.getTable<Opm::SlgofTable>(i), tolcrit) - swco[i];
};
const auto famII_2p = [&sof2Tables]( int i ) {
return critical_oil( sof2Tables.getTable<Opm::Sof2Table>( i ) );
const auto famII_2p = [&sof2Tables, tolcrit](const int i) -> double
{
return critical_oil(sof2Tables.getTable<Opm::Sof2Table>(i), tolcrit);
};
const auto famII_3p = [&sof3Tables]( int i ) {
const auto& tb = sof3Tables.getTable<Opm::Sof3Table>( i );
return critical_oil( tb, tb.getKrogColumn() );
const auto famII_3p = [&sof3Tables, tolcrit](const int i) -> double
{
const auto& tb = sof3Tables.getTable<Opm::Sof3Table>(i);
return critical_oil(tb, tb.getKrogColumn(), tolcrit);
};
switch( getSaturationFunctionFamily( tm, ph ) ) {
@@ -973,26 +1033,6 @@ namespace {
}
}
RawTableEndPoints
getRawTableEndpoints(const Opm::TableManager& tm,
const Opm::Phases& phases)
{
auto ep = RawTableEndPoints{};
ep.connate.gas = findMinGasSaturation(tm, phases);
ep.connate.water = findMinWaterSaturation(tm, phases);
ep.critical.oil_in_gas = findCriticalOilGas(tm, phases, ep.connate.water);
ep.critical.oil_in_water = findCriticalOilWater(tm, phases);
ep.critical.gas = findCriticalGas(tm, phases);
ep.critical.water = findCriticalWater(tm, phases);
ep.maximum.gas = findMaxGasSaturation(tm, phases);
ep.maximum.water = findMaxWaterSaturation(tm, phases);
return ep;
}
double selectValue(const Opm::TableContainer& depthTables,
int tableIdx,
const std::string& columnName,
@@ -1533,10 +1573,33 @@ namespace {
}
} // namespace Anonymous
std::shared_ptr<Opm::satfunc::RawTableEndPoints>
Opm::satfunc::getRawTableEndpoints(const Opm::TableManager& tm,
const Opm::Phases& phases,
const double tolcrit)
{
auto ep = std::make_shared<RawTableEndPoints>();
ep->connate.gas = findMinGasSaturation(tm, phases);
ep->connate.water = findMinWaterSaturation(tm, phases);
ep->critical.oil_in_gas = findCriticalOilGas(tm, phases, ep->connate.water, tolcrit);
ep->critical.oil_in_water = findCriticalOilWater(tm, phases, tolcrit);
ep->critical.gas = findCriticalGas(tm, phases, tolcrit);
ep->critical.water = findCriticalWater(tm, phases, tolcrit);
ep->maximum.gas = findMaxGasSaturation(tm, phases);
ep->maximum.water = findMaxWaterSaturation(tm, phases);
return ep;
}
std::vector<double>
Opm::satfunc::init(const std::string& keyword,
const TableManager& tables,
const Phases& phases,
const RawTableEndPoints& ep,
const std::vector<double>& cell_depth,
const std::vector<int>& num,
const std::vector<int>& endnum)
@@ -1585,7 +1648,5 @@ Opm::satfunc::init(const std::string& keyword,
+ keyword + '\''
};
const auto ep = getRawTableEndpoints(tables, phases);
return func->second(tables, phases, ep, cell_depth, num, endnum);
}

618
tests/parser/FieldPropsTests.cpp
View File

@@ -19,8 +19,10 @@
#include <stdexcept>
#include <iostream>
#include <iomanip>
#include <memory>
#include <numeric>
#include <sstream>
#include <string>
#define BOOST_TEST_MODULE FieldPropsTests
@@ -39,6 +41,7 @@
#include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/FieldPropsManager.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/EclipseGrid.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/SatfuncPropertyInitializers.hpp>
#include <opm/parser/eclipse/EclipseState/Runspec.hpp>
#include "src/opm/parser/eclipse/EclipseState/Grid/FieldProps.hpp"
@@ -588,6 +591,63 @@ MULTZ
BOOST_CHECK_EQUAL( multx[0], 2 );
}
BOOST_AUTO_TEST_CASE(OPERATE) {
std::string deck_string = R"(
GRID
PORO
6*1.0 /
OPERATE
PORO 1 3 1 1 1 1 'MAXLIM' PORO 0.50 /
PORO 1 3 2 2 1 1 'MAXLIM' PORO 0.25 /
/
PERMX
6*1/
PERMY
6*1000/
OPERATE
PERMX 1 3 1 1 1 1 'MINLIM' PERMX 2 /
PERMX 1 3 2 2 1 1 'MINLIM' PERMX 4 /
PERMY 1 3 1 1 1 1 'MAXLIM' PERMY 100 /
PERMY 1 3 2 2 1 1 'MAXLIM' PERMY 200 /
PERMZ 1 3 1 1 1 1 'MULTA' PERMY 2 1000 /
PERMZ 1 3 2 2 1 1 'MULTA' PERMX 3 300 /
/
)";
UnitSystem unit_system(UnitSystem::UnitType::UNIT_TYPE_METRIC);
auto to_si = [&unit_system](double raw_value) { return unit_system.to_si(UnitSystem::measure::permeability, raw_value); };
EclipseGrid grid(3,2,1);
Deck deck = Parser{}.parseString(deck_string);
FieldPropsManager fpm(deck, Phases{true, true, true}, grid, TableManager());
const auto& poro = fpm.get_double("PORO");
BOOST_CHECK_EQUAL(poro[0], 0.50);
BOOST_CHECK_EQUAL(poro[3], 0.25);
const auto& permx = fpm.get_double("PERMX");
BOOST_CHECK_EQUAL(permx[0], to_si(2));
BOOST_CHECK_EQUAL(permx[3], to_si(4));
const auto& permy = fpm.get_double("PERMY");
BOOST_CHECK_EQUAL(permy[0], to_si(100));
BOOST_CHECK_EQUAL(permy[3], to_si(200));
const auto& permz = fpm.get_double("PERMZ");
for (std::size_t i = 0; i < 3; i++) {
BOOST_CHECK_EQUAL(permz[i] , 2*permy[i] + to_si(1000));
BOOST_CHECK_EQUAL(permz[i+3], 3*permx[i+3] + to_si(300));
}
}
namespace {
std::string satfunc_model_setup()
{
@@ -869,6 +929,15 @@ SOF3
/ )" };
}
std::string tolCrit(const double t)
{
std::ostringstream os;
os << "TOLCRIT\n " << std::scientific << std::setw(11) << t << "\n/\n";
return os.str();
}
std::string end()
{
return { R"(
@@ -877,7 +946,370 @@ END
}
} // namespace Anonymous
BOOST_AUTO_TEST_CASE(SatFunc_EndPts_Family_I) {
BOOST_AUTO_TEST_CASE(RawTableEndPoints_Family_I_TolCrit_Zero) {
const auto es = ::Opm::EclipseState {
::Opm::Parser{}.parseString(satfunc_model_setup() + satfunc_family_I() + end())
};
const auto& tm = es.getTableManager();
const auto& ph = es.runspec().phases();
const auto tolcrit = 0.0;
auto rtepPtr = satfunc::getRawTableEndpoints(tm, ph, tolcrit);
// Water end-points
{
const auto swl = rtepPtr->connate .water;
const auto swcr = rtepPtr->critical.water;
const auto swu = rtepPtr->maximum .water;
BOOST_CHECK_CLOSE(swl [0], 0.071004, 1.0e-10);
BOOST_CHECK_CLOSE(swcr[0], 0.071004, 1.0e-10); // == SWL. TOLCRIT = 0.0
BOOST_CHECK_CLOSE(swu [0], 1.0 , 1.0e-10);
}
// Oil end-points
{
const auto sowcr = rtepPtr->critical.oil_in_water;
const auto sogcr = rtepPtr->critical.oil_in_gas;
BOOST_CHECK_CLOSE(sowcr[0], 1.0 - 0.791004, 1.0e-10); // TOLCRIT = 0.0
BOOST_CHECK_CLOSE(sogcr[0], 1.0 - 0.858996 - 0.071004, 1.0e-10); // Include SWL
}
// Gas end-points
{
const auto sgl = rtepPtr->connate .gas;
const auto sgcr = rtepPtr->critical.gas;
const auto sgu = rtepPtr->maximum .gas;
BOOST_CHECK_CLOSE(sgl [0], 0.0, 1.0e-10);
BOOST_CHECK_CLOSE(sgcr[0], 0.03, 1.0e-10);
BOOST_CHECK_CLOSE(sgu [0], 0.928996, 1.0e-10);
}
}
BOOST_AUTO_TEST_CASE(RawTableEndPoints_Family_II_TolCrit_Zero) {
const auto es = ::Opm::EclipseState {
::Opm::Parser{}.parseString(satfunc_model_setup() + satfunc_family_II() + end())
};
const auto& tm = es.getTableManager();
const auto& ph = es.runspec().phases();
const auto tolcrit = 0.0;
auto rtepPtr = satfunc::getRawTableEndpoints(tm, ph, tolcrit);
// Water end-points
{
const auto swl = rtepPtr->connate .water;
const auto swcr = rtepPtr->critical.water;
const auto swu = rtepPtr->maximum .water;
BOOST_CHECK_CLOSE(swl [0], 0.071004, 1.0e-10);
BOOST_CHECK_CLOSE(swcr[0], 0.071004, 1.0e-10); // == SWL. TOLCRIT = 0.0
BOOST_CHECK_CLOSE(swu [0], 1.0 , 1.0e-10);
}
// Oil end-points
{
const auto sowcr = rtepPtr->critical.oil_in_water;
const auto sogcr = rtepPtr->critical.oil_in_gas;
BOOST_CHECK_CLOSE(sowcr[0], 1.0 - 0.791004, 1.0e-10); // TOLCRIT = 0.0
BOOST_CHECK_CLOSE(sogcr[0], 1.0 - 0.858996 - 0.071004, 1.0e-10); // Include SWL
}
// Gas end-points
{
const auto sgl = rtepPtr->connate .gas;
const auto sgcr = rtepPtr->critical.gas;
const auto sgu = rtepPtr->maximum .gas;
BOOST_CHECK_CLOSE(sgl [0], 0.0, 1.0e-10);
BOOST_CHECK_CLOSE(sgcr[0], 0.03, 1.0e-10);
BOOST_CHECK_CLOSE(sgu [0], 0.928996, 1.0e-10);
}
}
BOOST_AUTO_TEST_CASE(RawTableEndPoints_Family_I_TolCrit_Default) {
const auto es = ::Opm::EclipseState {
::Opm::Parser{}.parseString(satfunc_model_setup() + satfunc_family_I() + end())
};
const auto& tm = es.getTableManager();
const auto& ph = es.runspec().phases();
const auto tolcrit = es.runspec().saturationFunctionControls()
.minimumRelpermMobilityThreshold(); // 1.0e-6.
auto rtepPtr = satfunc::getRawTableEndpoints(tm, ph, tolcrit);
// Water end-points
{
const auto swl = rtepPtr->connate .water;
const auto swcr = rtepPtr->critical.water;
const auto swu = rtepPtr->maximum .water;
BOOST_CHECK_CLOSE(swl [0], 0.071004, 1.0e-10);
BOOST_CHECK_CLOSE(swcr[0], 0.091004, 1.0e-10); // Max Sw for which Krw(Sw) <= TOLCRIT
BOOST_CHECK_CLOSE(swu [0], 1.0 , 1.0e-10);
}
// Oil end-points
{
const auto sowcr = rtepPtr->critical.oil_in_water;
const auto sogcr = rtepPtr->critical.oil_in_gas;
BOOST_CHECK_CLOSE(sowcr[0], 0.228996, 1.0e-10); // Max So for which Krow(So) <= TOLCRIT
BOOST_CHECK_CLOSE(sogcr[0], 0.070000, 1.0e-10); // Max So for which Krog(So) <= TOLCRIT
}
// Gas end-points
{
const auto sgl = rtepPtr->connate .gas;
const auto sgcr = rtepPtr->critical.gas;
const auto sgu = rtepPtr->maximum .gas;
BOOST_CHECK_CLOSE(sgl [0], 0.0, 1.0e-10);
BOOST_CHECK_CLOSE(sgcr[0], 0.03, 1.0e-10); // Max Sg for which Krg(Sg) <= TOLCRIT
BOOST_CHECK_CLOSE(sgu [0], 0.928996, 1.0e-10);
}
}
BOOST_AUTO_TEST_CASE(RawTableEndPoints_Family_II_TolCrit_Default) {
const auto es = ::Opm::EclipseState {
::Opm::Parser{}.parseString(satfunc_model_setup() + satfunc_family_II() + end())
};
const auto& tm = es.getTableManager();
const auto& ph = es.runspec().phases();
const auto tolcrit = es.runspec().saturationFunctionControls()
.minimumRelpermMobilityThreshold(); // 1.0e-6.
auto rtepPtr = satfunc::getRawTableEndpoints(tm, ph, tolcrit);
// Water end-points
{
const auto swl = rtepPtr->connate .water;
const auto swcr = rtepPtr->critical.water;
const auto swu = rtepPtr->maximum .water;
BOOST_CHECK_CLOSE(swl [0], 0.071004, 1.0e-10);
BOOST_CHECK_CLOSE(swcr[0], 0.091004, 1.0e-10); // Max Sw for which Krw(Sw) <= TOLCRIT
BOOST_CHECK_CLOSE(swu [0], 1.0 , 1.0e-10);
}
// Oil end-points
{
const auto sowcr = rtepPtr->critical.oil_in_water;
const auto sogcr = rtepPtr->critical.oil_in_gas;
BOOST_CHECK_CLOSE(sowcr[0], 0.228996, 1.0e-10); // Max So for which Krow(So) <= TOLCRIT
BOOST_CHECK_CLOSE(sogcr[0], 0.070000, 1.0e-10); // Max So for which Krog(So) <= TOLCRIT
}
// Gas end-points
{
const auto sgl = rtepPtr->connate .gas;
const auto sgcr = rtepPtr->critical.gas;
const auto sgu = rtepPtr->maximum .gas;
BOOST_CHECK_CLOSE(sgl [0], 0.0, 1.0e-10);
BOOST_CHECK_CLOSE(sgcr[0], 0.03, 1.0e-10); // Max Sg for which Krg(Sg) <= TOLCRIT
BOOST_CHECK_CLOSE(sgu [0], 0.928996, 1.0e-10);
}
}
BOOST_AUTO_TEST_CASE(RawTableEndPoints_Family_I_TolCrit_Large) {
const auto es = ::Opm::EclipseState {
::Opm::Parser{}.parseString(satfunc_model_setup() + satfunc_family_I() + end())
};
const auto& tm = es.getTableManager();
const auto& ph = es.runspec().phases();
const auto tolcrit = 0.01;
auto rtepPtr = satfunc::getRawTableEndpoints(tm, ph, tolcrit);
// Water end-points
{
const auto swl = rtepPtr->connate .water;
const auto swcr = rtepPtr->critical.water;
const auto swu = rtepPtr->maximum .water;
BOOST_CHECK_CLOSE(swl [0], 0.071004, 1.0e-10);
BOOST_CHECK_CLOSE(swcr[0], 0.231004, 1.0e-10); // Max Sw for which Krw(Sw) <= TOLCRIT
BOOST_CHECK_CLOSE(swu [0], 1.0 , 1.0e-10);
}
// Oil end-points
{
const auto sowcr = rtepPtr->critical.oil_in_water;
const auto sogcr = rtepPtr->critical.oil_in_gas;
BOOST_CHECK_CLOSE(sowcr[0], 0.448996, 1.0e-10); // Max So for which Krow(So) <= TOLCRIT
BOOST_CHECK_CLOSE(sogcr[0], 0.238996, 1.0e-10); // Max So for which Krog(So) <= TOLCRIT
}
// Gas end-points
{
const auto sgl = rtepPtr->connate .gas;
const auto sgcr = rtepPtr->critical.gas;
const auto sgu = rtepPtr->maximum .gas;
BOOST_CHECK_CLOSE(sgl [0], 0.0, 1.0e-10);
BOOST_CHECK_CLOSE(sgcr[0], 0.09, 1.0e-10); // Max Sg for which Krg(Sg) <= TOLCRIT
BOOST_CHECK_CLOSE(sgu [0], 0.928996, 1.0e-10);
}
}
BOOST_AUTO_TEST_CASE(RawTableEndPoints_Family_II_TolCrit_Large) {
const auto es = ::Opm::EclipseState {
::Opm::Parser{}.parseString(satfunc_model_setup() + satfunc_family_II() + end())
};
const auto& tm = es.getTableManager();
const auto& ph = es.runspec().phases();
const auto tolcrit = 0.01;
auto rtepPtr = satfunc::getRawTableEndpoints(tm, ph, tolcrit);
// Water end-points
{
const auto swl = rtepPtr->connate .water;
const auto swcr = rtepPtr->critical.water;
const auto swu = rtepPtr->maximum .water;
BOOST_CHECK_CLOSE(swl [0], 0.071004, 1.0e-10);
BOOST_CHECK_CLOSE(swcr[0], 0.231004, 1.0e-10); // Max Sw for which Krw(Sw) <= TOLCRIT
BOOST_CHECK_CLOSE(swu [0], 1.0 , 1.0e-10);
}
// Oil end-points
{
const auto sowcr = rtepPtr->critical.oil_in_water;
const auto sogcr = rtepPtr->critical.oil_in_gas;
BOOST_CHECK_CLOSE(sowcr[0], 0.448996, 1.0e-10); // Max So for which Krow(So) <= TOLCRIT
BOOST_CHECK_CLOSE(sogcr[0], 0.248996, 1.0e-10); // Max So for which Krog(So) <= TOLCRIT
}
// Gas end-points
{
const auto sgl = rtepPtr->connate .gas;
const auto sgcr = rtepPtr->critical.gas;
const auto sgu = rtepPtr->maximum .gas;
BOOST_CHECK_CLOSE(sgl [0], 0.0, 1.0e-10);
BOOST_CHECK_CLOSE(sgcr[0], 0.09, 1.0e-10); // Max Sg for which Krg(Sg) <= TOLCRIT
BOOST_CHECK_CLOSE(sgu [0], 0.928996, 1.0e-10);
}
}
// =====================================================================
BOOST_AUTO_TEST_CASE(SatFunc_EndPts_Family_I_TolCrit_Zero) {
const auto es = ::Opm::EclipseState {
::Opm::Parser{}.parseString(satfunc_model_setup() + tolCrit(0.0) + satfunc_family_I() + end())
};
auto fp = es.fieldProps();
// Water end-points
{
const auto swl = fp.get_double("SWL");
const auto swcr = fp.get_double("SWCR");
const auto swu = fp.get_double("SWU");
BOOST_CHECK_CLOSE(swl [0], 0.071004, 1.0e-10);
BOOST_CHECK_CLOSE(swcr[0], 0.071004, 1.0e-10); // == SWL. TOLCRIT = 0.0
BOOST_CHECK_CLOSE(swu [0], 1.0 , 1.0e-10);
const auto krwr = fp.get_double("KRWR"); // Krw(Sw=1-Sowcr-Sgl) = Krw(Sw=0.791004)
const auto krw = fp.get_double("KRW"); // Krw(Swmax) = Krw(Sw=1)
BOOST_CHECK_CLOSE(krwr[0], 0.911429, 1.0e-10);
BOOST_CHECK_CLOSE(krw [0], 1.0 , 1.0e-10);
const auto pcw = fp.get_double("PCW");
BOOST_CHECK_CLOSE(pcw[0], 7.847999*unit::barsa, 1.0e-10);
}
// Oil end-points
{
const auto sowcr = fp.get_double("SOWCR");
const auto sogcr = fp.get_double("SOGCR");
BOOST_CHECK_CLOSE(sowcr[0], 1.0 - 0.791004, 1.0e-10); // TOLCRIT = 0.0
BOOST_CHECK_CLOSE(sogcr[0], 1.0 - 0.858996 - 0.071004, 1.0e-10); // Include SWL
const auto krorw = fp.get_double("KRORW"); // Krow(So=1-Swcr-Sgl) = Krow(So=0.928996)
const auto krorg = fp.get_double("KRORG"); // Krog(So=1-Sgcr-Swl) = Krog(So=0.898996)
const auto kro = fp.get_double("KRO"); // Krow(So=Somax) = Krog(So=Somax)
BOOST_CHECK_CLOSE(krorw[0], 1.0, 1.0e-10); // TOLCRIT = 0.0
BOOST_CHECK_CLOSE(krorg[0], 0.896942, 1.0e-10);
BOOST_CHECK_CLOSE(kro [0], 1.0, 1.0e-10);
}
// Gas end-points
{
const auto sgl = fp.get_double("SGL");
const auto sgcr = fp.get_double("SGCR");
const auto sgu = fp.get_double("SGU");
BOOST_CHECK_CLOSE(sgl [0], 0.0, 1.0e-10);
BOOST_CHECK_CLOSE(sgcr[0], 0.03, 1.0e-10);
BOOST_CHECK_CLOSE(sgu [0], 0.928996, 1.0e-10);
const auto krgr = fp.get_double("KRGR"); // Krg(Sg=1-Sogcr-Swl) = Krg(Sg=0.858996)
const auto krg = fp.get_double("KRG"); // Krg(Sgmax) = Krg(Sg=0.928996)
BOOST_CHECK_CLOSE(krgr[0], 0.866135, 1.0e-10);
BOOST_CHECK_CLOSE(krg [0], 1.0 , 1.0e-10);
const auto pcg = fp.get_double("PCG");
BOOST_CHECK_CLOSE(pcg[0], 0.0, 1.0e-10);
}
}
BOOST_AUTO_TEST_CASE(SatFunc_EndPts_Family_II_TolCrit_Zero) {
const auto es = ::Opm::EclipseState {
::Opm::Parser{}.parseString(satfunc_model_setup() + tolCrit(0.0) + satfunc_family_II() + end())
};
auto fp = es.fieldProps();
// Water end-points
{
const auto swl = fp.get_double("SWL");
const auto swcr = fp.get_double("SWCR");
const auto swu = fp.get_double("SWU");
BOOST_CHECK_CLOSE(swl [0], 0.071004, 1.0e-10);
BOOST_CHECK_CLOSE(swcr[0], 0.071004, 1.0e-10); // == SWL. TOLCRIT = 0.0
BOOST_CHECK_CLOSE(swu [0], 1.0 , 1.0e-10);
const auto krwr = fp.get_double("KRWR"); // Krw(Sw=1-Sowcr-Sgl) = Krw(Sw=0.791004)
const auto krw = fp.get_double("KRW"); // Krw(Swmax) = Krw(Sw=1)
BOOST_CHECK_CLOSE(krwr[0], 0.911429, 1.0e-10);
BOOST_CHECK_CLOSE(krw [0], 1.0 , 1.0e-10);
const auto pcw = fp.get_double("PCW");
BOOST_CHECK_CLOSE(pcw[0], 7.847999*unit::barsa, 1.0e-10);
}
// Oil end-points
{
const auto sowcr = fp.get_double("SOWCR");
const auto sogcr = fp.get_double("SOGCR");
BOOST_CHECK_CLOSE(sowcr[0], 1.0 - 0.791004, 1.0e-10); // TOLCRIT = 0.0
BOOST_CHECK_CLOSE(sogcr[0], 1.0 - 0.858996 - 0.071004, 1.0e-10); // Include SWL
const auto krorw = fp.get_double("KRORW"); // Krow(So=1-Swcr-Sgl) = Krow(So=0.928996)
const auto krorg = fp.get_double("KRORG"); // Krog(So=1-Sgcr-Swl) = Krog(So=0.898996)
const auto kro = fp.get_double("KRO"); // Krow(So=Somax) = Krog(So=Somax)
BOOST_CHECK_CLOSE(krorw[0], 1.0, 1.0e-10); // TOLCRIT = 0.0
BOOST_CHECK_CLOSE(krorg[0], 0.896942, 1.0e-10);
BOOST_CHECK_CLOSE(kro [0], 1.0, 1.0e-10);
}
// Gas end-points
{
const auto sgl = fp.get_double("SGL");
const auto sgcr = fp.get_double("SGCR");
const auto sgu = fp.get_double("SGU");
BOOST_CHECK_CLOSE(sgl [0], 0.0, 1.0e-10);
BOOST_CHECK_CLOSE(sgcr[0], 0.03, 1.0e-10);
BOOST_CHECK_CLOSE(sgu [0], 0.928996, 1.0e-10);
const auto krgr = fp.get_double("KRGR"); // Krg(Sg=1-Sogcr-Swl) = Krg(Sg=0.858996)
const auto krg = fp.get_double("KRG"); // Krg(Sgmax) = Krg(Sg=0.928996)
BOOST_CHECK_CLOSE(krgr[0], 0.866135, 1.0e-10);
BOOST_CHECK_CLOSE(krg [0], 1.0 , 1.0e-10);
const auto pcg = fp.get_double("PCG");
BOOST_CHECK_CLOSE(pcg[0], 0.0, 1.0e-10);
}
}
BOOST_AUTO_TEST_CASE(SatFunc_EndPts_Family_I_TolCrit_Default) {
// TOLCRIT = 1.0e-6
const auto es = ::Opm::EclipseState {
::Opm::Parser{}.parseString(satfunc_model_setup() + satfunc_family_I() + end())
};
@@ -890,12 +1322,12 @@ BOOST_AUTO_TEST_CASE(SatFunc_EndPts_Family_I) {
const auto swcr = fp.get_double("SWCR");
const auto swu = fp.get_double("SWU");
BOOST_CHECK_CLOSE(swl [0], 0.071004, 1.0e-10);
BOOST_CHECK_CLOSE(swcr[0], 0.071004, 1.0e-10); // == SWL. TOLCRIT not currently implemented.
BOOST_CHECK_CLOSE(swcr[0], 0.091004, 1.0e-10); // Max Sw for which Krw(Sw) <= TOLCRIT
BOOST_CHECK_CLOSE(swu [0], 1.0 , 1.0e-10);
const auto krwr = fp.get_double("KRWR"); // Krw(Sw=1-Sowcr-Sgl) = Krw(Sw=0.791004)
const auto krwr = fp.get_double("KRWR"); // Krw(Sw=1-Sowcr-Sgl) = Krw(Sw=0.771004)
const auto krw = fp.get_double("KRW"); // Krw(Swmax) = Krw(Sw=1)
BOOST_CHECK_CLOSE(krwr[0], 0.911429, 1.0e-10);
BOOST_CHECK_CLOSE(krwr[0], 0.835916, 1.0e-10);
BOOST_CHECK_CLOSE(krw [0], 1.0 , 1.0e-10);
const auto pcw = fp.get_double("PCW");
@@ -906,13 +1338,13 @@ BOOST_AUTO_TEST_CASE(SatFunc_EndPts_Family_I) {
{
const auto sowcr = fp.get_double("SOWCR");
const auto sogcr = fp.get_double("SOGCR");
BOOST_CHECK_CLOSE(sowcr[0], 1.0 - 0.791004, 1.0e-10); // TOLCRIT currently not supported.
BOOST_CHECK_CLOSE(sowcr[0], 1.0 - 0.771004, 1.0e-10); // Max So for which Krow(So) <= TOLCRIT.
BOOST_CHECK_CLOSE(sogcr[0], 1.0 - 0.858996 - 0.071004, 1.0e-10); // Include SWL
const auto krorw = fp.get_double("KRORW"); // Krow(So=1-Swcr-Sgl) = Krow(So=0.928996)
const auto krorw = fp.get_double("KRORW"); // Krow(So=1-Swcr-Sgl) = Krow(So=0.908996)
const auto krorg = fp.get_double("KRORG"); // Krog(So=1-Sgcr-Swl) = Krog(So=0.898996)
const auto kro = fp.get_double("KRO"); // Krow(So=Somax) = Krog(So=Somax)
BOOST_CHECK_CLOSE(krorw[0], 1.0, 1.0e-10); // TOLCRIT currently not supported.
BOOST_CHECK_CLOSE(krorw[0], 0.882459, 1.0e-10);
BOOST_CHECK_CLOSE(krorg[0], 0.896942, 1.0e-10);
BOOST_CHECK_CLOSE(kro [0], 1.0, 1.0e-10);
}
@@ -936,7 +1368,8 @@ BOOST_AUTO_TEST_CASE(SatFunc_EndPts_Family_I) {
}
}
BOOST_AUTO_TEST_CASE(SatFunc_EndPts_Family_II) {
BOOST_AUTO_TEST_CASE(SatFunc_EndPts_Family_II_TolCrit_Default) {
// TOLCRIT = 1.0e-6
const auto es = ::Opm::EclipseState {
::Opm::Parser{}.parseString(satfunc_model_setup() + satfunc_family_II() + end())
};
@@ -949,12 +1382,12 @@ BOOST_AUTO_TEST_CASE(SatFunc_EndPts_Family_II) {
const auto swcr = fp.get_double("SWCR");
const auto swu = fp.get_double("SWU");
BOOST_CHECK_CLOSE(swl [0], 0.071004, 1.0e-10);
BOOST_CHECK_CLOSE(swcr[0], 0.071004, 1.0e-10); // == SWL. TOLCRIT not currently implemented.
BOOST_CHECK_CLOSE(swcr[0], 0.091004, 1.0e-10); // Max Sw for which Krw(Sw) <= TOLCRIT
BOOST_CHECK_CLOSE(swu [0], 1.0 , 1.0e-10);
const auto krwr = fp.get_double("KRWR"); // Krw(Sw=1-Sowcr-Sgl) = Krw(Sw=0.791004)
const auto krwr = fp.get_double("KRWR"); // Krw(Sw=1-Sowcr-Sgl) = Krw(Sw=0.771004)
const auto krw = fp.get_double("KRW"); // Krw(Swmax) = Krw(Sw=1)
BOOST_CHECK_CLOSE(krwr[0], 0.911429, 1.0e-10);
BOOST_CHECK_CLOSE(krwr[0], 0.835916, 1.0e-10);
BOOST_CHECK_CLOSE(krw [0], 1.0 , 1.0e-10);
const auto pcw = fp.get_double("PCW");
@@ -965,21 +1398,21 @@ BOOST_AUTO_TEST_CASE(SatFunc_EndPts_Family_II) {
{
const auto sowcr = fp.get_double("SOWCR");
const auto sogcr = fp.get_double("SOGCR");
BOOST_CHECK_CLOSE(sowcr[0], 1.0 - 0.791004, 1.0e-10); // TOLCRIT currently not supported.
BOOST_CHECK_CLOSE(sogcr[0], 1.0 - 0.858996 - 0.071004, 1.0e-10); // Include SWL
BOOST_CHECK_CLOSE(sowcr[0], 0.228996, 1.0e-10); // Max So for which Krow(So) <= TOLCRIT
BOOST_CHECK_CLOSE(sogcr[0], 0.070000, 1.0e-10); // Max So for which Krow(So) <= TOLCRIT
const auto krorw = fp.get_double("KRORW"); // Krow(So=1-Swcr-Sgl) = Krow(So=0.928996)
const auto krorw = fp.get_double("KRORW"); // Krow(So=1-Swcr-Sgl) = Krow(So=0.908996)
const auto krorg = fp.get_double("KRORG"); // Krog(So=1-Sgcr-Swl) = Krog(So=0.898996)
const auto kro = fp.get_double("KRO"); // Krow(So=Somax) = Krog(So=Somax)
BOOST_CHECK_CLOSE(krorw[0], 1.0, 1.0e-10); // TOLCRIT currently not supported.
BOOST_CHECK_CLOSE(krorw[0], 0.882459, 1.0e-10);
BOOST_CHECK_CLOSE(krorg[0], 0.896942, 1.0e-10);
BOOST_CHECK_CLOSE(kro [0], 1.0, 1.0e-10);
BOOST_CHECK_CLOSE(kro [0], 1.0 , 1.0e-10);
}
// Gas end-points
{
const auto sgl = fp.get_double("SGL");
const auto sgcr = fp.get_double("SGCR");
const auto sgcr = fp.get_double("SGCR"); // Max Sg for which Krg(Sg) <= TOLCRIT
const auto sgu = fp.get_double("SGU");
BOOST_CHECK_CLOSE(sgl [0], 0.0, 1.0e-10);
BOOST_CHECK_CLOSE(sgcr[0], 0.03, 1.0e-10);
@@ -995,61 +1428,122 @@ BOOST_AUTO_TEST_CASE(SatFunc_EndPts_Family_II) {
}
}
BOOST_AUTO_TEST_CASE(OPERATE) {
std::string deck_string = R"(
GRID
BOOST_AUTO_TEST_CASE(SatFunc_EndPts_Family_I_TolCrit_Large) {
// TOLCRIT = 0.01
const auto es = ::Opm::EclipseState {
::Opm::Parser{}.parseString(satfunc_model_setup() + tolCrit(0.01) + satfunc_family_I() + end())
};
PORO
6*1.0 /
auto fp = es.fieldProps();
OPERATE
PORO 1 3 1 1 1 1 'MAXLIM' PORO 0.50 /
PORO 1 3 2 2 1 1 'MAXLIM' PORO 0.25 /
/
// Water end-points
{
const auto swl = fp.get_double("SWL");
const auto swcr = fp.get_double("SWCR");
const auto swu = fp.get_double("SWU");
BOOST_CHECK_CLOSE(swl [0], 0.071004, 1.0e-10);
BOOST_CHECK_CLOSE(swcr[0], 0.231004, 1.0e-10); // Max Sw for which Krw(Sw) <= TOLCRIT
BOOST_CHECK_CLOSE(swu [0], 1.0 , 1.0e-10);
const auto krwr = fp.get_double("KRWR"); // Krw(Sw=1-Sowcr-Sgl) = Krw(Sw=0.551004)
const auto krw = fp.get_double("KRW"); // Krw(Swmax) = Krw(Sw=1)
BOOST_CHECK_CLOSE(krwr[0], 0.261115, 1.0e-10);
BOOST_CHECK_CLOSE(krw [0], 1.0 , 1.0e-10);
PERMX
6*1/
const auto pcw = fp.get_double("PCW");
BOOST_CHECK_CLOSE(pcw[0], 7.847999*unit::barsa, 1.0e-10);
}
PERMY
6*1000/
// Oil end-points
{
const auto sowcr = fp.get_double("SOWCR");
const auto sogcr = fp.get_double("SOGCR");
BOOST_CHECK_CLOSE(sowcr[0], 0.448996, 1.0e-10); // Max So for which Krow(So) <= TOLCRIT
BOOST_CHECK_CLOSE(sogcr[0], 0.238996, 1.0e-10); // Max So for which Krog(So) <= TOLCRIT
OPERATE
PERMX 1 3 1 1 1 1 'MINLIM' PERMX 2 /
PERMX 1 3 2 2 1 1 'MINLIM' PERMX 4 /
PERMY 1 3 1 1 1 1 'MAXLIM' PERMY 100 /
PERMY 1 3 2 2 1 1 'MAXLIM' PERMY 200 /
PERMZ 1 3 1 1 1 1 'MULTA' PERMY 2 1000 /
PERMZ 1 3 2 2 1 1 'MULTA' PERMX 3 300 /
/
const auto krorw = fp.get_double("KRORW"); // Krow(So=1-Swcr-Sgl) = Krow(So=0.768996)
const auto krorg = fp.get_double("KRORG"); // Krog(So=1-Sgcr-Swl) = Krog(So=0.838996)
const auto kro = fp.get_double("KRO"); // Krow(So=Somax) = Krog(So=Somax)
BOOST_CHECK_CLOSE(krorw[0], 0.328347, 1.0e-10);
BOOST_CHECK_CLOSE(krorg[0], 0.712749, 1.0e-10);
BOOST_CHECK_CLOSE(kro [0], 1.0, 1.0e-10);
}
// Gas end-points
{
const auto sgl = fp.get_double("SGL");
const auto sgcr = fp.get_double("SGCR");
const auto sgu = fp.get_double("SGU");
BOOST_CHECK_CLOSE(sgl [0], 0.0, 1.0e-10);
BOOST_CHECK_CLOSE(sgcr[0], 0.090000, 1.0e-10); // Max Sg for which Krg(Sg) <= TOLCRIT
BOOST_CHECK_CLOSE(sgu [0], 0.928996, 1.0e-10);
const auto krgr = fp.get_double("KRGR"); // Krg(Sg=1-Sogcr-Swl) = Krg(Sg=0.690000)
const auto krg = fp.get_double("KRG"); // Krg(Sgmax) = Krg(Sg=0.928996)
BOOST_CHECK_CLOSE(krgr[0], 0.578171, 1.0e-10);
BOOST_CHECK_CLOSE(krg [0], 1.0 , 1.0e-10);
)";
UnitSystem unit_system(UnitSystem::UnitType::UNIT_TYPE_METRIC);
auto to_si = [&unit_system](double raw_value) { return unit_system.to_si(UnitSystem::measure::permeability, raw_value); };
EclipseGrid grid(3,2,1);
Deck deck = Parser{}.parseString(deck_string);
FieldPropsManager fpm(deck, Phases{true, true, true}, grid, TableManager());
const auto& poro = fpm.get_double("PORO");
BOOST_CHECK_EQUAL(poro[0], 0.50);
BOOST_CHECK_EQUAL(poro[3], 0.25);
const auto& permx = fpm.get_double("PERMX");
BOOST_CHECK_EQUAL(permx[0], to_si(2));
BOOST_CHECK_EQUAL(permx[3], to_si(4));
const auto& permy = fpm.get_double("PERMY");
BOOST_CHECK_EQUAL(permy[0], to_si(100));
BOOST_CHECK_EQUAL(permy[3], to_si(200));
const auto& permz = fpm.get_double("PERMZ");
for (std::size_t i = 0; i < 3; i++) {
BOOST_CHECK_EQUAL(permz[i] , 2*permy[i] + to_si(1000));
BOOST_CHECK_EQUAL(permz[i+3], 3*permx[i+3] + to_si(300));
const auto pcg = fp.get_double("PCG");
BOOST_CHECK_CLOSE(pcg[0], 0.0, 1.0e-10);
}
}
BOOST_AUTO_TEST_CASE(SatFunc_EndPts_Family_II_TolCrit_Large) {
// TOLCRIT = 0.01
const auto es = ::Opm::EclipseState {
::Opm::Parser{}.parseString(satfunc_model_setup() + tolCrit(0.01) + satfunc_family_II() + end())
};
auto fp = es.fieldProps();
// Water end-points
{
const auto swl = fp.get_double("SWL");
const auto swcr = fp.get_double("SWCR");
const auto swu = fp.get_double("SWU");
BOOST_CHECK_CLOSE(swl [0], 0.071004, 1.0e-10);
BOOST_CHECK_CLOSE(swcr[0], 0.231004, 1.0e-10); // Max Sw for which Krw(Sw) <= TOLCRIT
BOOST_CHECK_CLOSE(swu [0], 1.0 , 1.0e-10);
const auto krwr = fp.get_double("KRWR"); // Krw(Sw=1-Sowcr-Sgl) = Krw(Sw=0.551004)
const auto krw = fp.get_double("KRW"); // Krw(Swmax) = Krw(Sw=1)
BOOST_CHECK_CLOSE(krwr[0], 0.261115, 1.0e-10);
BOOST_CHECK_CLOSE(krw [0], 1.0 , 1.0e-10);
const auto pcw = fp.get_double("PCW");
BOOST_CHECK_CLOSE(pcw[0], 7.847999*unit::barsa, 1.0e-10);
}
// Oil end-points
{
const auto sowcr = fp.get_double("SOWCR");
const auto sogcr = fp.get_double("SOGCR");
BOOST_CHECK_CLOSE(sowcr[0], 0.448996, 1.0e-10); // Max So for which Krow(So) <= TOLCRIT
BOOST_CHECK_CLOSE(sogcr[0], 0.248996, 1.0e-10); // Max So for which Krog(So) <= TOLCRIT
const auto krorw = fp.get_double("KRORW"); // Krow(So=1-Swcr-Sgl) = Krow(So=0.768996)
const auto krorg = fp.get_double("KRORG"); // Krog(So=1-Sgcr-Swl) = Krog(So=0.838996)
const auto kro = fp.get_double("KRO"); // Krow(So=Somax) = Krog(So=Somax)
BOOST_CHECK_CLOSE(krorw[0], 0.328347, 1.0e-10);
BOOST_CHECK_CLOSE(krorg[0], 0.712749, 1.0e-10);
BOOST_CHECK_CLOSE(kro [0], 1.0, 1.0e-10);
}
// Gas end-points
{
const auto sgl = fp.get_double("SGL");
const auto sgcr = fp.get_double("SGCR");
const auto sgu = fp.get_double("SGU");
BOOST_CHECK_CLOSE(sgl [0], 0.0, 1.0e-10);
BOOST_CHECK_CLOSE(sgcr[0], 0.090000, 1.0e-10); // Max Sg for which Krg(Sg) <= TOLCRIT
BOOST_CHECK_CLOSE(sgu [0], 0.928996, 1.0e-10);
const auto krgr = fp.get_double("KRGR"); // Krg(Sg=1-Sogcr-Swl) = Krg(Sg=0.680000)
const auto krg = fp.get_double("KRG"); // Krg(Sgmax) = Krg(Sg=0.928996)
BOOST_CHECK_CLOSE(krgr[0], 0.562914, 1.0e-10);
BOOST_CHECK_CLOSE(krg [0], 1.0 , 1.0e-10);
const auto pcg = fp.get_double("PCG");
BOOST_CHECK_CLOSE(pcg[0], 0.0, 1.0e-10);
}
}