PORV is treated more as a regular keyword

This commit is contained in:
Joakim Hove
2019-12-10 14:06:14 +01:00
parent 0f609083c3
commit 2fb0cd599b
5 changed files with 61 additions and 65 deletions
@@ -419,8 +419,6 @@ void FieldProps::reset_actnum(const std::vector<int>& new_actnum) {
FieldProps::compress(this->cell_volume, active_map);
FieldProps::compress(this->cell_depth, active_map);
if (this->porv_ptr)
this->porv_ptr.reset( nullptr );
this->m_actnum = std::move(new_actnum);
this->active_size = new_active_size;
@@ -507,6 +505,9 @@ FieldProps::FieldData<double>& FieldProps::get(const std::string& keyword) {
if (keywords::PROPS::satfunc.count(keyword) == 1)
this->init_satfunc(keyword, this->double_data[keyword]);
if (keyword == ParserKeywords::PORV::keywordName)
this->init_porv(this->double_data[keyword]);
return this->double_data[keyword];
} else
throw std::out_of_range("Double keyword: " + keyword + " is not supported");
@@ -531,8 +532,7 @@ FieldProps::FieldData<int>& FieldProps::get(const std::string& keyword) {
throw std::out_of_range("Integer keyword " + keyword + " is not supported");
}
std::vector<Box::cell_index> FieldProps::region_index( const DeckItem& region_item, int region_value ) {
std::string region_name = region_item.defaultApplied(0) ? this->m_default_region : make_region_name(region_item.get<std::string>(0));
std::vector<Box::cell_index> FieldProps::region_index( const std::string& region_name, int region_value ) {
const auto& region = this->get<int>(region_name);
if (!region.valid())
throw std::invalid_argument("Trying to work with invalid region: " + region_name);
@@ -550,6 +550,11 @@ std::vector<Box::cell_index> FieldProps::region_index( const DeckItem& region_it
return index_list;
}
std::vector<Box::cell_index> FieldProps::region_index( const DeckItem& region_item, int region_value ) {
std::string region_name = region_item.defaultApplied(0) ? this->m_default_region : make_region_name(region_item.get<std::string>(0));
return this->region_index(region_name, region_value);
}
template <>
bool FieldProps::has<double>(const std::string& keyword) const {
@@ -678,8 +683,10 @@ void FieldProps::apply(const DeckRecord& record, FieldData<T>& target_data, cons
if ((check_target == false) || (value::has_value(target_data.value_status[cell_index.active_index]))) {
target_data.data[cell_index.active_index] = func(target_data.data[cell_index.active_index], src_data.data[cell_index.active_index]);
target_data.value_status[cell_index.active_index] = src_data.value_status[cell_index.active_index];
}
}
} else
throw std::invalid_argument("Tried to use unset property value in OPERATE/OPERATER keyword");
} else
throw std::invalid_argument("Tried to use unset property value in OPERATE/OPERATER keyword");
}
}
@@ -687,17 +694,20 @@ void FieldProps::handle_region_operation(const DeckKeyword& keyword) {
for (const auto& record : keyword) {
const std::string& target_kw = record.getItem(0).get<std::string>(0);
int region_value = record.getItem("REGION_NUMBER").get<int>(0);
const auto& index_list = this->region_index(record.getItem("REGION_NAME"), region_value);
if (FieldProps::supported<double>(target_kw)) {
auto& field_data = this->get<double>(target_kw);
if (keyword.name() == ParserKeywords::OPERATER::keywordName) {
// For the OPERATER keyword we fetch the region name from the deck record
// with no extra hoops.
const auto& index_list = this->region_index(record.getItem("REGION_NAME").get<std::string>(0), region_value);
const std::string& src_kw = record.getItem("ARRAY_PARAMETER").get<std::string>(0);
const auto& src_data = this->get<double>(src_kw);
FieldProps::apply(record, field_data, src_data, index_list);
} else {
double value = record.getItem(1).get<double>(0);
const auto& index_list = this->region_index(record.getItem("REGION_NAME"), region_value);
if (keyword.name() != ParserKeywords::MULTIPLY::keywordName)
value = this->getSIValue(target_kw, value);
FieldProps::apply(fromString(keyword.name()), field_data, value, index_list);
@@ -811,56 +821,34 @@ void FieldProps::handle_keyword(const DeckKeyword& keyword, Box& box) {
/**********************************************************************/
std::vector<double> FieldProps::porv(bool global) {
if (!this->porv_ptr) {
FieldProps::FieldData<double> porv(this->active_size);
if (this->has<double>("PORV"))
porv = this->get<double>("PORV");
auto& porv_data = porv.data;
auto& porv_status = porv.value_status;
void FieldProps::init_porv(FieldData<double>& porv) {
auto& porv_data = porv.data;
auto& porv_status = porv.value_status;
if (!porv.valid()) {
const auto& poro = this->get<double>("PORO");
const auto& poro_status = poro.value_status;
const auto& poro_data = poro.data;
const auto& poro = this->get<double>("PORO");
const auto& poro_status = poro.value_status;
const auto& poro_data = poro.data;
for (std::size_t active_index = 0; active_index < this->active_size; active_index++) {
if (value::has_value(porv_status[active_index]))
continue;
if (value::has_value(poro_status[active_index])) {
porv_data[active_index] = this->cell_volume[active_index] * poro_data[active_index];
porv_status[active_index] = value::status::valid_default;
}
}
for (std::size_t active_index = 0; active_index < this->active_size; active_index++) {
if (value::has_value(poro_status[active_index])) {
porv_data[active_index] = this->cell_volume[active_index] * poro_data[active_index];
porv_status[active_index] = value::status::valid_default;
}
if (!porv.valid())
throw std::invalid_argument("Do not have enough information to create PORV");
// The NTG multiplication is only done one the cells which have PORV caclulated as PORO * V
if (this->has<double>("NTG")) {
const auto& ntg = this->get_valid_data<double>("NTG");
for (std::size_t active_index = 0; active_index < this->active_size; active_index++) {
if (porv_status[active_index] == value::status::valid_default)
porv_data[active_index] *= ntg[active_index];
}
}
// The MULTPV multiplication is done on all cells
if (this->has<double>("MULTPV")) {
const auto& multpv = this->get_valid_data<double>("MULTPV");
std::transform(porv_data.begin(), porv_data.end(), multpv.begin(), porv_data.begin(), std::multiplies<double>());
}
this->porv_ptr = std::make_unique<std::vector<double>>(porv_data);
}
if (global)
return this->global_copy(*this->porv_ptr);
return *this->porv_ptr;
if (this->has<double>("NTG")) {
const auto& ntg = this->get_valid_data<double>("NTG");
for (std::size_t active_index = 0; active_index < this->active_size; active_index++)
porv_data[active_index] *= ntg[active_index];
}
if (this->has<double>("MULTPV")) {
const auto& multpv = this->get_valid_data<double>("MULTPV");
printf("Doing MULTPV multiplication");
std::transform(porv_data.begin(), porv_data.end(), multpv.begin(), porv_data.begin(), std::multiplies<double>());
}
}
@@ -885,7 +873,6 @@ std::vector<double> FieldProps::porv(bool global) {
std::vector<int> FieldProps::actnum() {
auto actnum = this->m_actnum;
const auto& deck_actnum = this->get<int>("ACTNUM");
const auto& porv_data = this->porv(false);
std::vector<int> global_map(this->active_size);
{
@@ -899,6 +886,8 @@ std::vector<int> FieldProps::actnum() {
}
const auto& porv = this->get<double>("PORV");
const auto& porv_data = porv.data;
for (std::size_t active_index = 0; active_index < this->active_size; active_index++) {
auto global_index = global_map[active_index];
actnum[global_index] = deck_actnum.data[active_index];