Support for keywords ENPTVD and ENKRVD.

This commit is contained in:
osae 2012-12-17 14:02:30 +01:00
parent cf5eab8794
commit d1aada8d81
5 changed files with 292 additions and 63 deletions

View File

@ -118,6 +118,7 @@ namespace EclipseKeywords
string("PLYMAX"), string("TLMIXPAR"), string("WPOLYMER"),
string("GRUPTREE"), string("GCONINJE"), string("GCONPROD"),
string("WGRUPCON"), string("ENDSCALE"), string("SCALECRS"),
string("ENPTVD"), string("ENKRVD"),
// The following fields only have a dummy implementation
// that allows us to ignore them.
string("SWFN"),

View File

@ -194,6 +194,8 @@ public:
SPECIAL_FIELD(WGRUPCON)
SPECIAL_FIELD(ENDSCALE)
SPECIAL_FIELD(SCALECRS)
SPECIAL_FIELD(ENPTVD)
SPECIAL_FIELD(ENKRVD)
// The following fields only have a dummy implementation
// that allows us to ignore them.

View File

@ -378,7 +378,10 @@ namespace
}
// Replace default values -1 by linear interpolation
inline void insertDefaultValues(std::vector<std::vector<double> >& table, int ncol)
inline void insertDefaultValues(std::vector<std::vector<double> >& table,
int ncol,
double defaultOut = 0.0,
bool defaultInterpolation = true)
{
const int sz = table[0].size();
for (int k=1; k<ncol; ++k) {
@ -400,14 +403,19 @@ namespace
}
if (xv.empty()) {
// Nothing specified, the entire column is defaulted.
// We insert zeros.
// We insert default value.
for (int i=0; i<int(indx.size()); ++i) {
table[k][indx[i]] = 0.0;
table[k][indx[i]] = defaultOut;
}
} else if (defaultInterpolation) {
// Interpolate
for (int i=0; i<int(indx.size()); ++i) {
table[k][indx[i]] = linearInterpolationExtrap(xv, yv, x[i]);
}
} else {
// Interpolate
for (int i=0; i<int(indx.size()); ++i) {
table[k][indx[i]] = linearInterpolationExtrap(xv, yv, x[i]);
table[k][indx[i]] = linearInterpolation(xv, yv, x[i]);
}
}
}

View File

@ -2266,7 +2266,6 @@ struct ENDSCALE : public SpecialBase {
}
};
struct SCALECRS : public SpecialBase {
std::string scalecrs_;
@ -2303,6 +2302,181 @@ struct SCALECRS : public SpecialBase {
}
};
struct ENPTVD : public SpecialBase {
std::vector<std::vector<std::vector<double> > > table_;
std::vector<bool> mask_;
virtual std::string name() const {
return std::string("ENPTVD");
}
virtual void read(std::istream & is) {
table_.resize(5);
std::vector<std::vector<double> > sub_table(5);
while (!is.eof()) {
if(is.peek() == int('/')) {
if (sub_table[0].empty() && !(table_[0].empty())) {
is >> ignoreLine;
break;
} else {
THROW("Error reading ENPTVD data - none or incomplete table.");
}
}
std::vector<double> data(9,-1.0);
int nread = readDefaultedVectorData(is, data, 9);
if (nread != 9) {
THROW("Error reading ENPTVD data - depth and 8 saturations pr line.");
}
if (data[0] == -1.0) {
THROW("Error reading ENPTVD data - depth can not be defaulted.");
}
if ((data[4] != -1.0) || (data[5] != -1.0) || (data[6] != -1.0) || (data[8] != -1.0)) {
THROW("Error reading ENPTVD data - non-default values in column 5-7,9 not supported.");
}
sub_table[0].push_back(data[0]); //depth
sub_table[1].push_back(data[1]); //swl
sub_table[2].push_back(data[2]); //swcr
sub_table[3].push_back(data[3]); //swu
sub_table[4].push_back(data[7]); //sowcr
is >> ignoreWhitespace;
if(is.peek() == int('/')) {
is >> ignoreLine;
if (sub_table[0].size() >= 2) {
insertDefaultValues(sub_table, 5, -1.0, false);
std::vector<std::vector<double> >::iterator it_sub = sub_table.begin();
for(std::vector<std::vector<std::vector<double> > >::size_type i=0; i<table_.size(); ++i) {
table_[i].push_back(*it_sub);
(*it_sub).clear();
++it_sub;
}
} else {
THROW("Error reading ENPTVD data - minimum 2 lines pr sub-table.");
}
}
}
for (std::vector<std::vector<std::vector<double> > >::size_type i=1; i<table_.size(); ++i) {
mask_.push_back(false);
for (std::vector<std::vector<double> >::size_type j=0; j<table_[0].size(); ++j) {
if (table_[i][j][0] != -1.0) {
mask_.back() = true;
break;
}
}
}
}
virtual void write(std::ostream & os) const {
os << name() << '\n';
std::cout << "-----depth-------swl------swcr-------swu-----sowcr" << std::endl;
for (std::vector<std::vector<double> >::size_type j=0; j<table_[0].size(); ++j) {
std::cout << "--------------------------------------------------" << std::endl;
for (std::vector<double>::size_type k=0; k<table_[0][j].size(); ++k) {
for (std::vector<std::vector<std::vector<double> > >::size_type i=0; i<table_.size(); ++i) {
std::cout << std::setw(10) << table_[i][j][k];
}
std::cout << std::endl;
}
}
}
virtual void convertToSI(const EclipseUnits& units) {
//depths
for (std::vector<std::vector<double> >::size_type j=0; j<table_[0].size(); ++j) {
for (std::vector<double>::size_type k=0; k<table_[0][j].size(); ++k) {
table_[0][j][k] *= units.length;
}
}
}
};
struct ENKRVD : public SpecialBase {
std::vector<std::vector<std::vector<double> > > table_;
std::vector<bool> mask_;
virtual std::string name() const {
return std::string("ENKRVD");
}
virtual void read(std::istream & is) {
table_.resize(5);
std::vector<std::vector<double> > sub_table(5);
while (!is.eof()) {
if(is.peek() == int('/')) {
if (sub_table[0].empty() && !(table_[0].empty())) {
is >> ignoreLine;
break;
} else {
THROW("Error reading ENKRVD data - none or incomplete table.");
}
}
std::vector<double> data(8,-1.0);
int nread = readDefaultedVectorData(is, data, 8);
if (nread != 8) {
THROW("Error reading ENKRVD data - depth and 7 relperms pr line.");
}
if (data[0] == -1.0) {
THROW("Error reading ENKRVD data - depth can not be defaulted.");
}
if ((data[2] != -1.0) || (data[5] != -1.0) || (data[6] != -1.0)) {
THROW("Error reading ENKRVD data - non-default values in column 3,6-7 not supported.");
}
sub_table[0].push_back(data[0]); //depth
sub_table[1].push_back(data[1]); //krw
sub_table[2].push_back(data[3]); //kro
sub_table[3].push_back(data[4]); //krw(sowcr)
sub_table[4].push_back(data[7]); //kro(swcr)
is >> ignoreWhitespace;
if(is.peek() == int('/')) {
is >> ignoreLine;
if (sub_table[0].size() >= 2) {
insertDefaultValues(sub_table, 5, -1.0, false);
std::vector<std::vector<double> >::iterator it_sub = sub_table.begin();
for(std::vector<std::vector<std::vector<double> > >::size_type i=0; i<table_.size(); ++i) {
table_[i].push_back(*it_sub);
(*it_sub).clear();
++it_sub;
}
} else {
THROW("Error reading ENKRVD data - minimum 2 lines pr sub-table.");
}
}
}
for (std::vector<std::vector<std::vector<double> > >::size_type i=1; i<table_.size(); ++i) {
mask_.push_back(false);
for (std::vector<std::vector<double> >::size_type j=0; j<table_[0].size(); ++j) {
if (table_[i][j][0] != -1.0) {
mask_.back() = true;
break;
}
}
}
}
virtual void write(std::ostream & os) const {
os << name() << '\n';
std::cout << "-----depth-------krw------krow------krwr-----krorw" << std::endl;
for (std::vector<std::vector<double> >::size_type j=0; j<table_[0].size(); ++j) {
std::cout << "--------------------------------------------------" << std::endl;
for (std::vector<double>::size_type k=0; k<table_[0][j].size(); ++k) {
for (std::vector<std::vector<std::vector<double> > >::size_type i=0; i<table_.size(); ++i) {
std::cout << std::setw(10) << table_[i][j][k];
}
std::cout << std::endl;
}
}
}
virtual void convertToSI(const EclipseUnits& units) {
//depths
for (std::vector<std::vector<double> >::size_type j=0; j<table_[0].size(); ++j) {
for (std::vector<double>::size_type k=0; k<table_[0][j].size(); ++k) {
table_[0][j][k] *= units.length;
}
}
}
};
// The following fields only have a dummy implementation
// that allows us to ignore them.

View File

@ -124,32 +124,6 @@ namespace Opm
initEPS(deck, grid, std::string("KRWR"), eps_.krwr_);
initEPS(deck, grid, std::string("KRO"), eps_.kro_);
initEPS(deck, grid, std::string("KRORW"), eps_.krorw_);
/*
double ss[PhaseUsage::MaxNumPhases], kr[PhaseUsage::MaxNumPhases];
int oldP = std::cout.precision();
std::cout.precision(4);
for (unsigned int i=0; i<=100; ++i) {
ss[phase_usage_.phase_pos[Aqua]] = i*0.01;
ss[phase_usage_.phase_pos[Liquid]] = 1.0 - ss[phase_usage_.phase_pos[Aqua]];
endScaling(ss, 15, kr);
std::cout << std::showpoint
<< std::setw(10) << ss[phase_usage_.phase_pos[Aqua]]
<< std::setw(12) << kr[phase_usage_.phase_pos[Aqua]]
<< std::setw(10) << kr[phase_usage_.phase_pos[Liquid]];
endScaling(ss, 45, kr);
std::cout << std::setw(12) << kr[phase_usage_.phase_pos[Aqua]]
<< std::setw(10) << kr[phase_usage_.phase_pos[Liquid]];
endScaling(ss, 75, kr);
std::cout << std::setw(12) << kr[phase_usage_.phase_pos[Aqua]]
<< std::setw(10) << kr[phase_usage_.phase_pos[Liquid]];
endScaling(ss, 105, kr);
std::cout << std::setw(12) << kr[phase_usage_.phase_pos[Aqua]]
<< std::setw(10) << kr[phase_usage_.phase_pos[Liquid]]
<< std::noshowpoint << std::endl;
}
std::cout.precision(oldP);
*/
}
}
@ -280,44 +254,94 @@ namespace Opm
// Initialize saturation scaling parameter
template <class SatFuncSet>
void SaturationPropsFromDeck<SatFuncSet>::initEPS(const EclipseGridParser& deck,
const UnstructuredGrid& grid,
const std::string& keyword,
std::vector<double>& scaleparam)
const UnstructuredGrid& grid,
const std::string& keyword,
std::vector<double>& scaleparam)
{
if (deck.hasField(keyword)) {
bool useKeyword = deck.hasField(keyword);
bool hasENPTVD = deck.hasField("ENPTVD");
bool hasENKRVD = deck.hasField("ENKRVD");
int itab = 0;
std::vector<std::vector<std::vector<double> > > table_dummy;
std::vector<std::vector<std::vector<double> > >& table = table_dummy;
// Active keyword assigned default values for each cell (in case of possible box-wise assignment)
scaleparam.resize(grid.number_of_cells);
int phase_pos_aqua = phase_usage_.phase_pos[BlackoilPhases::Aqua];
if (keyword == std::string("SWCR")) {
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).swcr_;
} else if (keyword == std::string("SWL")) {
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).smin_[phase_pos_aqua];
// Active keyword assigned default values for each cell (in case of possible box-wise assignment)
int phase_pos_aqua = phase_usage_.phase_pos[BlackoilPhases::Aqua];
if (keyword[0] == 'S' && (useKeyword || hasENPTVD)) {
if (keyword == std::string("SWL")) {
if (useKeyword || deck.getENPTVD().mask_[0]) {
itab = 1;
scaleparam.resize(grid.number_of_cells);
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).smin_[phase_pos_aqua];
}
} else if (keyword == std::string("SWCR")) {
if (useKeyword || deck.getENPTVD().mask_[1]) {
itab = 2;
scaleparam.resize(grid.number_of_cells);
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).swcr_;
}
} else if (keyword == std::string("SWU")) {
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).smax_[phase_pos_aqua];
if (useKeyword || deck.getENPTVD().mask_[2]) {
itab = 3;
scaleparam.resize(grid.number_of_cells);
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).smax_[phase_pos_aqua];
}
} else if (keyword == std::string("SOWCR")) {
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).sowcr_;
} else if (keyword == std::string("KRW")) {
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).krwmax_;
} else if (keyword == std::string("KRWR")) {
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).krwr_;
} else if (keyword == std::string("KRO")) {
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).kromax_;
} else if (keyword == std::string("KRORW")) {
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).krorw_;
} else {
THROW("SaturationPropsFromDeck::initEndscale() -- unknown keyword: '" << keyword << "'");
if (useKeyword || deck.getENPTVD().mask_[3]) {
itab = 4;
scaleparam.resize(grid.number_of_cells);
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).sowcr_;
}
}else {
THROW(" -- unknown keyword: '" << keyword << "'");
}
if (!useKeyword && itab > 0) {
table = deck.getENPTVD().table_;
}
} else if (keyword[0] == 'K' && (useKeyword || hasENKRVD)) {
if (keyword == std::string("KRW")) {
if (useKeyword || deck.getENKRVD().mask_[0]) {
itab = 1;
scaleparam.resize(grid.number_of_cells);
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).krwmax_;
}
} else if (keyword == std::string("KRO")) {
if (useKeyword || deck.getENKRVD().mask_[1]) {
itab = 2;
scaleparam.resize(grid.number_of_cells);
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).kromax_;
}
} else if (keyword == std::string("KRWR")) {
if (useKeyword || deck.getENKRVD().mask_[2]) {
itab = 3;
scaleparam.resize(grid.number_of_cells);
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).krwr_;
}
} else if (keyword == std::string("KRORW")) {
if (useKeyword || deck.getENKRVD().mask_[3]) {
itab = 4;
scaleparam.resize(grid.number_of_cells);
for (int i=0; i<grid.number_of_cells; ++i)
scaleparam[i] = funcForCell(i).krorw_;
}
} else {
THROW(" -- unknown keyword: '" << keyword << "'");
}
if (!useKeyword && itab > 0) {
table = deck.getENKRVD().table_;
}
}
if (scaleparam.empty()) {
return;
} else if (useKeyword) {
// Keyword values from deck
std::cout << "--- Scaling parameter '" << keyword << "' assigned." << std::endl;
const int* gc = grid.global_cell;
@ -326,8 +350,28 @@ namespace Opm
const int deck_pos = (gc == NULL) ? c : gc[c];
scaleparam[c] = val[deck_pos];
}
} else {
std::cout << "--- Scaling parameter '" << keyword << "' assigned via ";
if (keyword[0] == 'S')
deck.getENPTVD().write(std::cout);
else
deck.getENKRVD().write(std::cout);
const double* cc = grid.cell_centroids;
const int dim = grid.dimensions;
for (int cell = 0; cell < grid.number_of_cells; ++cell) {
int jtab = cell_to_func_.empty() ? 0 : cell_to_func_[cell];
if (table[itab][jtab][0] != -1.0) {
std::vector<double>& depth = table[0][jtab];
std::vector<double>& val = table[itab][jtab];
double zc = cc[dim*cell+dim-1];
if (zc >= depth.front() && zc <= depth.back()) { //don't want extrap outside depth interval
scaleparam[cell] = linearInterpolation(depth, val, zc);
}
}
}
}
}
// Saturation scaling
template <class SatFuncSet>