diff --git a/opm/parser/eclipse/EclipseState/Grid/FieldPropsManager.hpp b/opm/parser/eclipse/EclipseState/Grid/FieldPropsManager.hpp index b82e1d7eb..38704b132 100644 --- a/opm/parser/eclipse/EclipseState/Grid/FieldPropsManager.hpp +++ b/opm/parser/eclipse/EclipseState/Grid/FieldPropsManager.hpp @@ -186,6 +186,7 @@ public: virtual bool has_int(const std::string& keyword) const { return this->has(keyword); } virtual bool has_double(const std::string& keyword) const { return this->has(keyword); } + virtual void apply_tran(const std::string& keyword, std::vector& tranx) const; private: /* diff --git a/src/opm/parser/eclipse/EclipseState/Grid/FieldProps.cpp b/src/opm/parser/eclipse/EclipseState/Grid/FieldProps.cpp index 0d658d630..02c47207c 100644 --- a/src/opm/parser/eclipse/EclipseState/Grid/FieldProps.cpp +++ b/src/opm/parser/eclipse/EclipseState/Grid/FieldProps.cpp @@ -265,21 +265,21 @@ std::string make_region_name(const std::string& deck_value) { throw std::invalid_argument("The input string: " + deck_value + " was invalid. Expected: O/F/M"); } -FieldProps::ScalarOperation fromString(const std::string& keyword) { +ScalarOperation fromString(const std::string& keyword) { if (keyword == ParserKeywords::ADD::keywordName || keyword == ParserKeywords::ADDREG::keywordName) - return FieldProps::ScalarOperation::ADD; + return ScalarOperation::ADD; if (keyword == ParserKeywords::EQUALS::keywordName || keyword == ParserKeywords::EQUALREG::keywordName) - return FieldProps::ScalarOperation::EQUAL; + return ScalarOperation::EQUAL; if (keyword == ParserKeywords::MULTIPLY::keywordName || keyword == ParserKeywords::MULTIREG::keywordName) - return FieldProps::ScalarOperation::MUL; + return ScalarOperation::MUL; if (keyword == ParserKeywords::MINVALUE::keywordName) - return FieldProps::ScalarOperation::MIN; + return ScalarOperation::MIN; if (keyword == ParserKeywords::MAXVALUE::keywordName) - return FieldProps::ScalarOperation::MAX; + return ScalarOperation::MAX; throw std::invalid_argument("Keyword operation not recognized"); } @@ -324,6 +324,10 @@ FieldProps::FieldProps(const Deck& deck, const Phases& phases, const EclipseGrid grid_ptr(&grid), tables(tables_arg) { + this->tran.emplace( "TRANX", TranCalculator("TRANX") ); + this->tran.emplace( "TRANY", TranCalculator("TRANY") ); + this->tran.emplace( "TRANZ", TranCalculator("TRANZ") ); + if (deck.hasKeyword()) { const DeckKeyword& multregpKeyword = deck.getKeyword("MULTREGP"); for (const auto& record : multregpKeyword) { @@ -470,12 +474,11 @@ bool FieldProps::supported(const std::string& keyword) { template <> -FieldProps::FieldData& FieldProps::init_get(const std::string& keyword) { +FieldProps::FieldData& FieldProps::init_get(const std::string& keyword, const keywords::keyword_info& kw_info) { auto iter = this->double_data.find(keyword); if (iter != this->double_data.end()) return iter->second; - const keywords::keyword_info& kw_info = keywords::global_kw_info(keyword); this->double_data[keyword] = FieldData(kw_info, this->active_size, kw_info.global ? this->global_size : 0); if (keyword == ParserKeywords::PORV::keywordName) @@ -490,26 +493,33 @@ FieldProps::FieldData& FieldProps::init_get(const std::string& keyword) return this->double_data[keyword]; } - +template <> +FieldProps::FieldData& FieldProps::init_get(const std::string& keyword) { + keywords::keyword_info kw_info = keywords::global_kw_info(keyword); + return this->init_get(keyword, kw_info); +} template <> -FieldProps::FieldData& FieldProps::init_get(const std::string& keyword) { +FieldProps::FieldData& FieldProps::init_get(const std::string& keyword, const keywords::keyword_info& kw_info) { auto iter = this->int_data.find(keyword); if (iter != this->int_data.end()) return iter->second; + this->int_data[keyword] = FieldData(kw_info, this->active_size, kw_info.global ? this->global_size : 0); + return this->int_data[keyword]; +} + +template <> +FieldProps::FieldData& FieldProps::init_get(const std::string& keyword) { if (keywords::isFipxxx(keyword)) { auto kw_info = keywords::keyword_info{}; kw_info.init(1); - this->int_data[keyword] = FieldData(kw_info, this->active_size, 0); + return this->init_get(keyword, kw_info); } else { const keywords::keyword_info& kw_info = keywords::global_kw_info(keyword); - this->int_data[keyword] = FieldData(kw_info, this->active_size, kw_info.global ? this->global_size : 0); + return this->init_get(keyword, kw_info); } - - - return this->int_data[keyword]; } @@ -629,8 +639,8 @@ void FieldProps::handle_int_keyword(const keywords::keyword_info& kw_info, } -void FieldProps::handle_double_keyword(Section section, const keywords::keyword_info& kw_info, const DeckKeyword& keyword, const Box& box) { - auto& field_data = this->init_get(keyword.name()); +void FieldProps::handle_double_keyword(Section section, const keywords::keyword_info& kw_info, const DeckKeyword& keyword, const std::string& keyword_name, const Box& box) { + auto& field_data = this->init_get(keyword_name, kw_info); const auto& deck_data = keyword.getSIDoubleData(); const auto& deck_status = keyword.getValueStatus(); @@ -649,6 +659,9 @@ void FieldProps::handle_double_keyword(Section section, const keywords::keyword_ } } +void FieldProps::handle_double_keyword(Section section, const keywords::keyword_info& kw_info, const DeckKeyword& keyword, const Box& box) { + this->handle_double_keyword(section, kw_info, keyword, keyword.name(), box ); +} @@ -762,20 +775,37 @@ void FieldProps::handle_operation(const DeckKeyword& keyword, Box box) { const std::string& target_kw = record.getItem(0).get(0); box.update(record); - if (FieldProps::supported(target_kw)) { - auto& field_data = this->init_get(target_kw); - + if (FieldProps::supported(target_kw) || this->tran.count(target_kw) > 0) { if (keyword.name() == ParserKeywords::OPERATE::keywordName) { + auto& field_data = this->init_get(target_kw); const std::string& src_kw = record.getItem("ARRAY").get(0); const auto& src_data = this->init_get(src_kw); FieldProps::operate(record, field_data, src_data, box.index_list()); } else { + std::string unique_name = target_kw; + auto operation = fromString(keyword.name()); double scalar_value = record.getItem(1).get(0); + keywords::keyword_info kw_info; + auto tran_iter = this->tran.find(target_kw); + if (tran_iter != this->tran.end()) { + kw_info = tran_iter->second.make_kw_info(operation); + unique_name = tran_iter->second.next_name(); + tran_iter->second.add_action(operation, unique_name); + } else + kw_info = keywords::global_kw_info(target_kw); + + auto& field_data = this->init_get(unique_name, kw_info); + if (keyword.name() != ParserKeywords::MULTIPLY::keywordName) scalar_value = this->getSIValue(target_kw, scalar_value); - FieldProps::apply(fromString(keyword.name()), field_data.data, field_data.value_status, scalar_value, box.index_list()); - if (field_data.global_data) - FieldProps::apply(fromString(keyword.name()), *field_data.global_data, *field_data.global_value_status, scalar_value, box.global_index_list()); + + if (tran_iter != this->tran.end()) { + assign_scalar(field_data.data, field_data.value_status, scalar_value, box.index_list()); + } else { + FieldProps::apply(operation, field_data.data, field_data.value_status, scalar_value, box.index_list()); + if (field_data.global_data) + FieldProps::apply(operation, *field_data.global_data, *field_data.global_value_status, scalar_value, box.global_index_list()); + } } continue; @@ -974,6 +1004,17 @@ void FieldProps::scanEDITSection(const EDITSection& edit_section) { Box box(*this->grid_ptr); for (const auto& keyword : edit_section) { const std::string& name = keyword.name(); + + auto tran_iter = this->tran.find(name); + if (tran_iter!= this->tran.end()) { + auto& tran_calc = tran_iter->second; + auto unique_name = tran_calc.next_name(); + keywords::keyword_info kw_info; + this->handle_double_keyword(Section::EDIT, kw_info, keyword, unique_name, box); + tran_calc.add_action( ScalarOperation::EQUAL, unique_name ); + continue; + } + if (keywords::EDIT::double_keywords.count(name) == 1) { this->handle_double_keyword(Section::EDIT, keywords::EDIT::double_keywords.at(name), keyword, box); continue; @@ -1087,6 +1128,37 @@ const std::string& FieldProps::default_region() const { } +void FieldProps::apply_tran(const std::string& keyword, std::vector& data) { + const auto& calculator = this->tran.at(keyword); + for (const auto& action : calculator) { + const auto& action_data = this->double_data.at(action.field); + + for (std::size_t index = 0; index < this->active_size; index++) { + + if (action_data.value_status[index] != value::status::deck_value) + continue; + + switch (action.op) { + case ScalarOperation::EQUAL: + data[index] = action_data.data[index]; + break; + + case ScalarOperation::MUL: + data[index] *= action_data.data[index]; + break; + + case ScalarOperation::MAX: + data[index] = std::max(action_data.data[index], data[index]); + break; + + default: + throw std::logic_error("Unhandled value in switch"); + } + } + } +} + + template std::vector FieldProps::defaulted(const std::string& keyword); template std::vector FieldProps::defaulted(const std::string& keyword); } diff --git a/src/opm/parser/eclipse/EclipseState/Grid/FieldProps.hpp b/src/opm/parser/eclipse/EclipseState/Grid/FieldProps.hpp index 17e9b4985..f5b31d8d3 100644 --- a/src/opm/parser/eclipse/EclipseState/Grid/FieldProps.hpp +++ b/src/opm/parser/eclipse/EclipseState/Grid/FieldProps.hpp @@ -20,6 +20,7 @@ #define FIELDPROPS_HPP #include +#include #include #include #include @@ -162,9 +163,9 @@ static const std::unordered_map> int_keywords = { namespace EDIT { static const std::unordered_map> double_keywords = {{"MULTPV", keyword_info{}.init(1.0)}, {"PORV", keyword_info{}.unit_string("ReservoirVolume")}, - {"TRANX", keyword_info{}.unit_string("Transmissibility").init(1.0)}, - {"TRANY", keyword_info{}.unit_string("Transmissibility").init(1.0)}, - {"TRANZ", keyword_info{}.unit_string("Transmissibility").init(1.0)}, + {"TRANX", keyword_info{}.unit_string("Transmissibility")}, + {"TRANY", keyword_info{}.unit_string("Transmissibility")}, + {"TRANZ", keyword_info{}.unit_string("Transmissibility")}, {"MULTX", keyword_info{}.init(1.0).mult(true)}, {"MULTX-", keyword_info{}.init(1.0).mult(true)}, {"MULTY", keyword_info{}.init(1.0).mult(true)}, @@ -275,10 +276,87 @@ keyword_info global_kw_info(const std::string& name); } +enum class ScalarOperation { + ADD = 1, + EQUAL = 2, + MUL = 3, + MIN = 4, + MAX = 5 +}; + + +class TranCalculator { +public: + + struct TranAction { + ScalarOperation op; + std::string field; + }; + + + TranCalculator(const std::string& name_arg) : + m_name(name_arg) + {} + + std::string next_name() const { + return this->m_name + std::to_string( this->actions.size() ); + } + + std::vector::const_iterator begin() const { + return this->actions.begin(); + } + + std::vector::const_iterator end() const { + return this->actions.end(); + } + + void add_action(ScalarOperation op, const std::string& field) { + this->actions.push_back(TranAction{op, field}); + } + + std::size_t size() const { + return this->actions.size(); + } + + const std::string& name() const { + return this->m_name; + } + + + keywords::keyword_info make_kw_info(ScalarOperation op) { + keywords::keyword_info kw_info; + switch (op) { + case ScalarOperation::MUL: + kw_info.init(1); + break; + case ScalarOperation::ADD: + kw_info.init(0); + break; + case ScalarOperation::MAX: + kw_info.init(std::numeric_limits::min()); + break; + case ScalarOperation::MIN: + kw_info.init(std::numeric_limits::max()); + break; + default: + break; + } + return kw_info; + } + +private: + std::string m_name; + std::vector actions; +}; + + + class FieldProps { public: + + struct MultregpRecord { int region_value; double multiplier; @@ -293,13 +371,6 @@ public: }; - enum class ScalarOperation { - ADD = 1, - EQUAL = 2, - MUL = 3, - MIN = 4, - MAX = 5 - }; template static void compress(std::vector& data, const std::vector& active_map) { @@ -584,6 +655,8 @@ public: return this->double_data.size(); } + void apply_tran(const std::string& keyword, std::vector& data); + private: void scanGRIDSection(const GRIDSection& grid_section); void scanEDITSection(const EDITSection& edit_section); @@ -607,6 +680,9 @@ private: template FieldData& init_get(const std::string& keyword); + template + FieldData& init_get(const std::string& keyword, const keywords::keyword_info& kw_info); + std::string region_name(const DeckItem& region_item); std::vector region_index( const std::string& region_name, int region_value ); void handle_operation(const DeckKeyword& keyword, Box box); @@ -617,6 +693,7 @@ private: double get_alpha(const std::string& func_name, const std::string& target_array, double raw_alpha); void handle_keyword(const DeckKeyword& keyword, Box& box); + void handle_double_keyword(Section section, const keywords::keyword_info& kw_info, const DeckKeyword& keyword, const std::string& keyword_name, const Box& box); void handle_double_keyword(Section section, const keywords::keyword_info& kw_info, const DeckKeyword& keyword, const Box& box); void handle_int_keyword(const keywords::keyword_info& kw_info, const DeckKeyword& keyword, const Box& box); void init_satfunc(const std::string& keyword, FieldData& satfunc); @@ -637,6 +714,8 @@ private: std::vector multregp; std::unordered_map> int_data; std::unordered_map> double_data; + + std::unordered_map tran; }; } diff --git a/src/opm/parser/eclipse/EclipseState/Grid/FieldPropsManager.cpp b/src/opm/parser/eclipse/EclipseState/Grid/FieldPropsManager.cpp index 0089987c0..bbde7ea30 100644 --- a/src/opm/parser/eclipse/EclipseState/Grid/FieldPropsManager.cpp +++ b/src/opm/parser/eclipse/EclipseState/Grid/FieldPropsManager.cpp @@ -109,6 +109,11 @@ std::size_t FieldPropsManager::active_size() const { return this->fp->active_size; } +void FieldPropsManager::apply_tran(const std::string& keyword, std::vector& data) const { + this->fp->apply_tran(keyword, data); +} + + template bool FieldPropsManager::supported(const std::string&); template bool FieldPropsManager::supported(const std::string&); diff --git a/tests/parser/FieldPropsTests.cpp b/tests/parser/FieldPropsTests.cpp index 52c5e2f13..81b5f6df7 100644 --- a/tests/parser/FieldPropsTests.cpp +++ b/tests/parser/FieldPropsTests.cpp @@ -1778,3 +1778,67 @@ OPERATE BOOST_CHECK_THROW(make_fp(invalid_region), std::logic_error); BOOST_CHECK_THROW(make_fp(invalid_operate), std::logic_error); } + + + +BOOST_AUTO_TEST_CASE(TRAN_Calculator) { + std::string deck_string = R"( +GRID + +PORO + 1000*0.10 / + +EDIT + +TRANX + 1000*0.10 / + +BOX + 1 10 1 10 1 1 / + +TRANX + 100*0.20 / + +ENDBOX + +MULTIPLY + TRANY 2.0 / +/ + +MAXVALUE + TRANZ 0 1 10 1 10 1 1 / +/ + + +)"; + UnitSystem unit_system(UnitSystem::UnitType::UNIT_TYPE_METRIC); + auto to_si = [&unit_system](double raw_value) { return unit_system.to_si(UnitSystem::measure::transmissibility, raw_value); }; + std::vector actnum(1000, 1); + for (std::size_t i=0; i< 1000; i += 2) + actnum[i] = 0; + EclipseGrid grid(EclipseGrid(10,10,10), actnum); + Deck deck = Parser{}.parseString(deck_string); + FieldPropsManager fpm(deck, Phases{true, true, true}, grid, TableManager()); + std::vector tranx( grid.getNumActive() ); + std::vector trany( grid.getNumActive(), to_si(1.0) ); + std::vector tranz( grid.getNumActive(), to_si(1.0) ); + BOOST_CHECK(!fpm.has_double("TRANX")); + + BOOST_CHECK_THROW(fpm.apply_tran("TRANA", tranx), std::out_of_range); + + fpm.apply_tran("TRANX", tranx); + fpm.apply_tran("TRANY", trany); + fpm.apply_tran("TRANZ", tranz); + + + for (std::size_t i=0; i < 50; i++) + BOOST_CHECK_EQUAL(tranx[i], to_si(0.20)); + + for (std::size_t i=100; i < tranx.size(); i++) + BOOST_CHECK_EQUAL(tranx[i], to_si(0.10)); + + BOOST_CHECK_EQUAL(trany[0], to_si(2.0)); + for (std::size_t i=0; i < trany.size(); i++) + BOOST_CHECK_EQUAL(trany[i], to_si(2.0)); + +}