Merge pull request #1129 from totto82/test_swatinit_fix

Make it optinal to apply SWATINIT
This commit is contained in:
Atgeirr Flø Rasmussen 2017-01-06 10:23:39 +01:00 committed by GitHub
commit 3e3a76d37d
5 changed files with 201 additions and 19 deletions

View File

@ -168,6 +168,7 @@ list (APPEND TEST_DATA_FILES
tests/liveoil.DATA
tests/capillary.DATA
tests/capillary_overlap.DATA
tests/capillarySwatinit.DATA
tests/compressed_gridproperty.data
tests/deadfluids.DATA
tests/equil_livegas.DATA

View File

@ -68,6 +68,8 @@ namespace Opm
* \param[in] props Property object, pvt and capillary properties are used.
* \param[in] deck Simulation deck, used to obtain EQUIL and related data.
* \param[in] gravity Acceleration of gravity, assumed to be in Z direction.
* \param[in] applySwatInit Make it possible to not apply SWATINIT even if it
* is present in the deck
*/
template<class Grid>
void initStateEquil(const Grid& grid,
@ -75,7 +77,8 @@ namespace Opm
const Opm::Deck& deck,
const Opm::EclipseState& eclipseState,
const double gravity,
BlackoilState& state);
BlackoilState& state,
bool applySwatInit = true);
/**
@ -246,13 +249,17 @@ namespace Opm
const Opm::Deck& deck,
const Opm::EclipseState& eclipseState,
const Grid& G ,
const double grav = unit::gravity)
const double grav = unit::gravity,
const std::vector<double>& swat_init = {}
)
: pp_(props.numPhases(),
std::vector<double>(UgGridHelpers::numCells(G))),
sat_(props.numPhases(),
std::vector<double>(UgGridHelpers::numCells(G))),
rs_(UgGridHelpers::numCells(G)),
rv_(UgGridHelpers::numCells(G))
rv_(UgGridHelpers::numCells(G)),
swat_init_(swat_init)
{
// Get the equilibration records.
const std::vector<EquilRecord> rec = getEquil(eclipseState);
@ -340,20 +347,6 @@ namespace Opm
}
}
// Check for presence of kw SWATINIT
if (deck.hasKeyword("SWATINIT")) {
const std::vector<double>& swat_init = eclipseState.
get3DProperties().getDoubleGridProperty("SWATINIT").getData();
const int nc = UgGridHelpers::numCells(G);
swat_init_.resize(nc);
const int* gc = UgGridHelpers::globalCell(G);
for (int c = 0; c < nc; ++c) {
const int deck_pos = (gc == NULL) ? c : gc[c];
swat_init_[c] = swat_init[deck_pos];
}
}
// Compute pressures, saturations, rs and rv factors.
calcPressSatRsRv(eqlmap, rec, props, G, grav);

View File

@ -882,6 +882,8 @@ namespace Opm
* \param[in] props Property object, pvt and capillary properties are used.
* \param[in] deck Simulation deck, used to obtain EQUIL and related data.
* \param[in] gravity Acceleration of gravity, assumed to be in Z direction.
* \param[in] applySwatInit Make it possible to not apply SWATINIT even if it
* is present in the deck
*/
template<class Grid>
void initStateEquil(const Grid& grid,
@ -889,10 +891,26 @@ namespace Opm
const Opm::Deck& deck,
const Opm::EclipseState& eclipseState,
const double gravity,
BlackoilState& state)
BlackoilState& state,
bool applySwatinit = true)
{
typedef EQUIL::DeckDependent::InitialStateComputer ISC;
ISC isc(props, deck, eclipseState, grid, gravity);
//Check for presence of kw SWATINIT
std::vector<double> swat_init = {};
if (eclipseState.get3DProperties().hasDeckDoubleGridProperty("SWATINIT") && applySwatinit) {
const std::vector<double>& swat_init_ecl = eclipseState.
get3DProperties().getDoubleGridProperty("SWATINIT").getData();
const int nc = UgGridHelpers::numCells(grid);
swat_init.resize(nc);
const int* gc = UgGridHelpers::globalCell(grid);
for (int c = 0; c < nc; ++c) {
const int deck_pos = (gc == NULL) ? c : gc[c];
swat_init[c] = swat_init_ecl[deck_pos];
}
}
ISC isc(props, deck, eclipseState, grid, gravity, swat_init);
const auto pu = props.phaseUsage();
const int ref_phase = pu.phase_used[BlackoilPhases::Liquid]
? pu.phase_pos[BlackoilPhases::Liquid]

View File

@ -0,0 +1,86 @@
-- Most of the following sections are not actually needed by the test,
-- but it is required by the Eclipse reference manual that they are
-- present. Also, the higher level opm-parser classes
-- (i.e. Opm::EclipseState et al.) assume that they are present.
-------------------------------------
RUNSPEC
WATER
OIL
GAS
DIMENS
40 40 40 /
TABDIMS
1 1 40 20 1 20 /
EQLDIMS
-- NTEQUL
1 /
-------------------------------------
GRID
-- Opm::EclipseState assumes that _some_ grid gets defined, so let's
-- specify a fake one...
DXV
40*1 /
DYV
40*1 /
DZV
40*1 /
DEPTHZ
1681*123.456 /
-------------------------------------
PROPS
PVDO
100 1.0 1.0
200 0.9 1.0
/
PVDG
100 0.010 0.1
200 0.005 0.2
/
PVTW
1.0 1.0 4.0E-5 0.96 0.0
/
SWOF
0.2 0 1 0.4
1 1 0 0.1
/
SGOF
0 0 1 0.2
0.8 1 0 0.5
/
DENSITY
700 1000 1
/
SWATINIT
5*0
10*0.5
5*1 /
-------------------------------------
SOLUTION
EQUIL
50 150 50 0.25 20 0.35 1* 1* 0
/
-------------------------------------
SCHEDULE
-- empty section

View File

@ -809,4 +809,88 @@ BOOST_AUTO_TEST_CASE (DeckWithRSVDAndRVVD)
}
}
BOOST_AUTO_TEST_CASE (DeckWithSwatinit)
{
Opm::GridManager gm(1, 1, 20, 1.0, 1.0, 5.0);
const UnstructuredGrid& grid = *(gm.c_grid());
Opm::Parser parser;
Opm::ParseContext parseContext;
Opm::Deck deck = parser.parseFile("capillarySwatinit.DATA" , parseContext);
Opm::EclipseState eclipseState(deck , parseContext);
Opm::BlackoilPropertiesFromDeck props(deck, eclipseState, grid, false);
Opm::BlackoilState state( Opm::UgGridHelpers::numCells( grid ) , Opm::UgGridHelpers::numFaces( grid ) , 3);
// reference saturations
const std::vector<double> s[3]{
{ 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.42192000000000002, 0.77802666666666664, 1, 1, 1, 1, 1, 1, 1, 1, 1 },
{ 0, 0, 0, 0.00736, 0.792746666666, 0.8, 0.8, 0.8, 0.8, 0.57807999999999993, 0.22197333333333336, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
{ 0.8, 0.8, 0.8, 0.79264, 0.007253333333, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 }
};
// sw in cell 13 and 14 is forced to be swu=1 since P_oil - P_wat < 0.
const std::vector<double> swatinit[3]{
{ 0.2, 0.2, 0.2, 0.2, 0.2, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 1, 1, 1, 1, 1, 1, 1, 1 },
{ 0, 0, 0, 0.00736, 0.792746666666, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0, 0, 0, 0, 0, 0, 0, 0 },
{ 0.8, 0.8, 0.8, 0.79264, 0.007253333333, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 }
};
std::vector<double> sats = state.saturation();
for (int phase = 0; phase < 3; ++phase) {
for (size_t i = 0; i < 20; ++i) {
sats[3*i + phase] = s[phase][i];
}
}
// reference pcs
const int numCells = Opm::UgGridHelpers::numCells(grid);
std::vector<int> cells(numCells);
for (int c = 0; c < numCells; ++c) { cells[c] = c; }
std::vector<double> pc_original = state.saturation();
props.capPress(numCells, sats.data(), cells.data(), pc_original.data(), nullptr);
std::vector<double> pc_scaled_truth = pc_original;
// only modify pcow
// sw = 0.2
for (size_t i = 0; i < 5; ++i) {
pc_scaled_truth[3*i + 0] = 40000;
}
// sw = 0.5
for (size_t i = 5; i < 12; ++i) {
pc_scaled_truth[3*i + 0] = 28750;
}
// sw = 1
for (size_t i = 12; i < 20; ++i) {
pc_scaled_truth[3*i + 0] = 10000;
}
// compute the initial state
// apply swatinit
Opm::BlackoilState state_scaled = state;
initStateEquil(grid, props, deck, eclipseState, 10.0, state_scaled, true);
// don't apply swatinit
Opm::BlackoilState state_unscaled = state;
initStateEquil(grid, props, deck, eclipseState, 10.0, state_unscaled, false);
// compute pc
std::vector<double> pc_scaled= state.saturation();
props.capPress(numCells, state_scaled.saturation().data(), cells.data(), pc_scaled.data(), nullptr);
std::vector<double> pc_unscaled= state.saturation();
props.capPress(numCells, state_unscaled.saturation().data(), cells.data(), pc_unscaled.data(), nullptr);
// test
const double reltol = 1.0e-6;
for (int phase = 0; phase < 3; ++phase) {
for (size_t i = 0; i < 20; ++i) {
CHECK( pc_original[3*i + phase ], pc_unscaled[3*i + phase ], reltol);
CHECK( pc_scaled_truth[3*i + phase], pc_scaled[3*i + phase ], reltol);
}
}
for (int phase = 0; phase < 3; ++phase) {
for (size_t i = 0; i < 20; ++i) {
CHECK(state_unscaled.saturation()[3*i + phase], s[phase][i], reltol);
CHECK(state_scaled.saturation()[3*i + phase], swatinit[phase][i], reltol);
}
}
}
BOOST_AUTO_TEST_SUITE_END()