diff --git a/opm/input/eclipse/Schedule/UDQ/UDQAssign.cpp b/opm/input/eclipse/Schedule/UDQ/UDQAssign.cpp index 88cb1bd1824..dc245c06696 100644 --- a/opm/input/eclipse/Schedule/UDQ/UDQAssign.cpp +++ b/opm/input/eclipse/Schedule/UDQ/UDQAssign.cpp @@ -19,6 +19,8 @@ #include +#include + #include #include @@ -36,7 +38,6 @@ namespace Opm { void UDQAssign::AssignRecord::eval(UDQSet& values) const { if (this->input_selector.empty() && - this->rst_selector.empty() && this->numbered_selector.empty()) { values.assign(this->value); @@ -44,11 +45,6 @@ void UDQAssign::AssignRecord::eval(UDQSet& values) const else if (! this->input_selector.empty()) { values.assign(this->input_selector[0], this->value); } - else if (! this->rst_selector.empty()) { - for (const auto& wgname : this->rst_selector) { - values.assign(wgname, this->value); - } - } else { for (const auto& item : this->numbered_selector) { for (const auto& number : item.numbers) { @@ -61,12 +57,13 @@ void UDQAssign::AssignRecord::eval(UDQSet& values) const bool UDQAssign::AssignRecord::operator==(const AssignRecord& data) const { return (this->input_selector == data.input_selector) - && (this->rst_selector == data.rst_selector) && (this->numbered_selector == data.numbered_selector) && (this->report_step == data.report_step) && (this->value == data.value); } +// --------------------------------------------------------------------------- + UDQAssign::UDQAssign(const std::string& keyword, const std::vector& input_selector, const double value, @@ -77,16 +74,6 @@ UDQAssign::UDQAssign(const std::string& keyword, this->add_record(input_selector, value, report_step); } -UDQAssign::UDQAssign(const std::string& keyword, - const std::unordered_set& rst_selector, - const double value, - const std::size_t report_step) - : m_keyword (keyword) - , m_var_type(UDQ::varType(keyword)) -{ - this->add_record(rst_selector, value, report_step); -} - UDQAssign::UDQAssign(const std::string& keyword, const std::vector& selector, double value, @@ -107,6 +94,15 @@ UDQAssign::UDQAssign(const std::string& keyword, this->add_record(std::move(selector), value, report_step); } +UDQAssign::UDQAssign(const std::string& keyword, + const RestartIO::RstUDQ& assignRst, + const std::size_t report_step) + : m_keyword { keyword } + , m_var_type { assignRst.category } +{ + this->add_record(assignRst, report_step); +} + UDQAssign UDQAssign::serializationTestObject() { UDQAssign result; @@ -114,8 +110,6 @@ UDQAssign UDQAssign::serializationTestObject() result.m_var_type = UDQVarType::CONNECTION_VAR; result.records.emplace_back(std::vector{"test1"}, 1.0, 0); - result.records.emplace_back(std::unordered_set{ "I-45" }, 3.1415, 3); - // Class-template argument deduction for the vector element type. result.records.emplace_back(std::vector { UDQSet::EnumeratedItems::serializationTestObject() }, 2.71828, 42); @@ -129,13 +123,6 @@ void UDQAssign::add_record(const std::vector& input_selector, this->records.emplace_back(input_selector, value, report_step); } -void UDQAssign::add_record(const std::unordered_set& rst_selector, - const double value, - const std::size_t report_step) -{ - this->records.emplace_back(rst_selector, value, report_step); -} - void UDQAssign::add_record(const std::vector& selector, const double value, const std::size_t report_step) @@ -150,6 +137,34 @@ void UDQAssign::add_record(std::vector&& selector, this->records.emplace_back(std::move(selector), value, report_step); } +void UDQAssign::add_record(const RestartIO::RstUDQ& assignRst, + const std::size_t report_step) +{ + if (assignRst.name != this->m_keyword) { + throw std::invalid_argument { + fmt::format("ASSIGN UDQ '{}' must not attempt to include " + "information for unrelated UDQ '{}' from restart file.", + this->m_keyword, assignRst.name) + }; + } + + switch (assignRst.category) { + case UDQVarType::SCALAR: + case UDQVarType::FIELD_VAR: + this->add_record(assignRst.entityNames(), + assignRst.scalarValue(), report_step); + break; + + case UDQVarType::WELL_VAR: + case UDQVarType::GROUP_VAR: + this->add_well_or_group_records(assignRst, report_step); + break; + + default: + break; + } +} + const std::string& UDQAssign::keyword() const { return this->m_keyword; @@ -225,4 +240,30 @@ bool UDQAssign::operator==(const UDQAssign& data) const && (this->records == data.records); } +// --------------------------------------------------------------------------- +// Private member functions below separator +// --------------------------------------------------------------------------- + +void UDQAssign::add_well_or_group_records(const RestartIO::RstUDQ& assignRst, + const std::size_t report_step) +{ + const auto& wgnames = assignRst.entityNames(); + const auto& nameIdx = assignRst.nameIndex(); + const auto n = assignRst.numEntities(); + + // Note: We intentionally allocate a single selector string and reuse + // that for every add_record() call. The loop here guarantees that we + // handle the case of different values for well or group, albeit at the + // cost of the 'records' data member being larger than necessary if all + // entities do have the same value. + auto selector = std::vector(1); + + for (auto i = 0*n; i < n; ++i) { + for (const auto& subValuePair : assignRst[i]) { + selector.front() = wgnames[nameIdx[i]]; + this->add_record(selector, subValuePair.second, report_step); + } + } +} + } // namespace Opm diff --git a/opm/input/eclipse/Schedule/UDQ/UDQAssign.hpp b/opm/input/eclipse/Schedule/UDQ/UDQAssign.hpp index 65e2703efce..ecbe7229bef 100644 --- a/opm/input/eclipse/Schedule/UDQ/UDQAssign.hpp +++ b/opm/input/eclipse/Schedule/UDQ/UDQAssign.hpp @@ -29,25 +29,283 @@ #include #include +namespace Opm::RestartIO { + class RstUDQ; +} + namespace Opm { +/// Representation of a UDQ ASSIGN statement class UDQAssign { public: + /// Default constructor. + UDQAssign() = default; + + /// Constructor. + /// + /// \param[in] keyword UDQ name. + /// + /// \param[in] selector Collection of entity names to which this + /// assignment applies. Might, for instance, be a selection of well or + /// group names for a well/group level UDQ, or an empty vector for a + /// scalar/field level UDQ. + /// + /// \param[in] value Numeric UDQ value for the entities named in the \p + /// selector. + /// + /// \param[in] report_step Time at which this assignment happens. + /// Assignments should be performed exactly once and the time value + /// ensures this behaviour. + UDQAssign(const std::string& keyword, + const std::vector& selector, + double value, + std::size_t report_step); + + /// Constructor. + /// + /// \param[in] keyword UDQ name. + /// + /// \param[in] selector Collection of named and numbered entities to + /// which this assignment applies. Might, for instance, be a selection + /// of well segments for a segment level UDQ. + /// + /// \param[in] value Numeric UDQ value for the entities identified in + /// the \p selector. + /// + /// \param[in] report_step Time at which this assignment happens. + /// Assignments should be performed exactly once and the time value + /// ensures this behaviour. + UDQAssign(const std::string& keyword, + const std::vector& selector, + double value, + std::size_t report_step); + + /// Constructor. + /// + /// Assumes ownership over the selector. + /// + /// \param[in] keyword UDQ name. + /// + /// \param[in] selector Collection of named and numbered entities to + /// which this assignment applies. Might, for instance, be a selection + /// of well segments for a segment level UDQ. + /// + /// \param[in] value Numeric UDQ value for the entities identified in + /// the \p selector. + /// + /// \param[in] report_step Time at which this assignment happens. + /// Assignments should be performed exactly once and the time value + /// ensures this behaviour. + UDQAssign(const std::string& keyword, + std::vector&& selector, + double value, + std::size_t report_step); + + /// Constructor. + /// + /// Reconstitutes an assignment from restart file information + /// + /// \param[in] keyword UDQ name. + /// + /// \param[in] assignRst Aggregate UDQ assignment information restored + /// from restart file information. + /// + /// \param[in] report_step Time at which this assignment happens. + /// Assignments should be performed exactly once and the time value + /// ensures this behaviour. + UDQAssign(const std::string& keyword, + const RestartIO::RstUDQ& assignRst, + const std::size_t report_step); + + /// Create a serialisation test object. + static UDQAssign serializationTestObject(); + + /// Name of UDQ to which this assignment applies. + const std::string& keyword() const; + + /// Kind of UDQ to which this assignment applies. + UDQVarType var_type() const; + + /// Add new record to existing UDQ assignment. + /// + /// \param[in] selector Collection of entity names to which this + /// assignment applies. Might, for instance, be a selection of well or + /// group names for a well/group level UDQ, or an empty vector for a + /// scalar/field level UDQ. + /// + /// \param[in] value Numeric UDQ value for the entities named in the \p + /// selector. + /// + /// \param[in] report_step Time at which this assignment happens. + /// Assignments should be performed exactly once and the time value + /// ensures this behaviour. + void add_record(const std::vector& selector, + double value, + std::size_t report_step); + + /// Add new record to existing UDQ assignment. + /// + /// \param[in] selector Collection of named and numbered entities to + /// which this assignment applies. Might, for instance, be a selection + /// of well segments for a segment level UDQ. + /// + /// \param[in] value Numeric UDQ value for the entities identified in + /// the \p selector. + /// + /// \param[in] report_step Time at which this assignment happens. + /// Assignments should be performed exactly once and the time value + /// ensures this behaviour. + void add_record(const std::vector& selector, + double value, + std::size_t report_step); + + /// Add new record to existing UDQ assignment. + /// + /// Assumes ownership over the selector. + /// + /// \param[in] selector Collection of named and numbered entities to + /// which this assignment applies. Might, for instance, be a selection + /// of well segments for a segment level UDQ. + /// + /// \param[in] value Numeric UDQ value for the entities identified in + /// the \p selector. + /// + /// \param[in] report_step Time at which this assignment happens. + /// Assignments should be performed exactly once and the time value + /// ensures this behaviour. + void add_record(std::vector&& selector, + double value, + std::size_t report_step); + + /// Add new record to existing UDQ assignment. + /// + /// Reconstitutes an assignment from restart file information. Mostly + /// needed for interface compatibility in generic code. + /// + /// \param[in] assignRst Aggregate UDQ assignment information restored + /// from restart file information. + /// + /// \param[in] report_step Time at which this assignment happens. + /// Assignments should be performed exactly once and the time value + /// ensures this behaviour. + void add_record(const RestartIO::RstUDQ& assignRst, + const std::size_t report_step); + + /// Apply current assignment to a selection of enumerated items. + /// + /// \param[in] items Collection of named and numbered entities to which + /// apply this assignment. Might, for instance, be a selection of well + /// segments for a segment level UDQ. + /// + /// \return UDQ set of the current assignment applied to \p items. + /// Items known at construction time, or defined in subsequent calls to + /// member function add_record(), will have a defined value in the + /// resulting UDQ set while unrecognised items will have an undefined + /// value. + UDQSet eval(const std::vector& items) const; + + /// Apply current assignment to a selection of named items. + /// + /// \param[in] wells Collection of named entities to which apply this + /// assignment. Might, for instance, be a selection of well names for a + /// well level UDQ. + /// + /// \return UDQ set of the current assignment applied to \p wells. + /// Named items known at construction time, or defined in subsequent + /// calls to member function add_record(), will have a defined value in + /// the resulting UDQ set while unrecognised items will have an + /// undefined value. + UDQSet eval(const std::vector& wells) const; + + /// Construct scalar UDQ set for a scalar UDQ assignment + /// + /// Throws an exception of type \code std::invalid_argument \endcode if + /// this assignment statement does not pertain to a scalar or field + /// level UDQ. + /// + /// \return UDQ set of the current assignment of a scalar UDQ. + UDQSet eval() const; + + /// Time at which this assignment happens. + /// + /// Assignments should be performed exactly once and the time value + /// ensures this behaviour. + /// + /// \return Constructor argument of the same name. + std::size_t report_step() const; + + /// Equality predicate. + /// + /// \param[in] data Object against which \code *this \endcode will be + /// tested for equality. + /// + /// \return Whether or not \code *this \endcode is the same as \p data. + bool operator==(const UDQAssign& data) const; + + /// Convert between byte array and object representation. + /// + /// \tparam Serializer Byte array conversion protocol. + /// + /// \param[in,out] serializer Byte array conversion object. + template + void serializeOp(Serializer& serializer) + { + serializer(m_keyword); + serializer(m_var_type); + serializer(records); + } + +private: // If the same keyword is assigned several times the different // assignment records are assembled in one UDQAssign instance. This is // an attempt to support restart in a situation where a full UDQ ASSIGN // statement can be swapped with a UDQ DEFINE statement. struct AssignRecord { + /// Collection of entity names to which this assignment applies. + /// + /// Might, for instance, be a selection of well or group names for a + /// well/group level UDQ, or an empty vector for a scalar/field + /// level UDQ. + /// + /// Empty for an enumerated assignment record. std::vector input_selector{}; - std::unordered_set rst_selector{}; + + /// Collection of named and numbered entities to which this + /// assignment applies. + /// + /// Might, for instance, be a selection of well segments for a + /// segment level UDQ. + /// + /// Empty for a named assignment record. std::vector numbered_selector{}; + + /// Numeric UDQ value for the entities identified in the selector. double value{}; + + /// Time at which this assignment happens. + /// + /// Assignments should be performed exactly once and the time value + /// ensures this behaviour. std::size_t report_step{}; + /// Default constructor. AssignRecord() = default; + /// Constructor. + /// + /// \param[in] selector Collection of entity names to which this + /// assignment applies. Might, for instance, be a selection of well + /// or group names for a well/group level UDQ, or an empty vector + /// for a scalar/field level UDQ. + /// + /// \param[in] value_arg Numeric UDQ value for the entities named in + /// the \p selector. + /// + /// \param[in] report_step_arg Time at which this assignment + /// happens. Assignments should be performed exactly once and the + /// time value ensures this behaviour. AssignRecord(const std::vector& selector, const double value_arg, const std::size_t report_step_arg) @@ -56,14 +314,18 @@ class UDQAssign , report_step (report_step_arg) {} - AssignRecord(const std::unordered_set& selector, - const double value_arg, - const std::size_t report_step_arg) - : rst_selector(selector) - , value (value_arg) - , report_step (report_step_arg) - {} - + /// Constructor. + /// + /// \param[in] selector Collection of named and numbered entities to + /// which this assignment applies. Might, for instance, be a + /// selection of well segments for a segment level UDQ. + /// + /// \param[in] value_arg Numeric UDQ value for the entities + /// identified in the \p selector. + /// + /// \param[in] report_step_arg Time at which this assignment + /// happens. Assignments should be performed exactly once and the + /// time value ensures this behaviour. AssignRecord(const std::vector& selector, const double value_arg, const std::size_t report_step_arg) @@ -72,6 +334,20 @@ class UDQAssign , report_step (report_step_arg) {} + /// Constructor. + /// + /// Assumes ownership over the selector. + /// + /// \param[in] selector Collection of named and numbered entities to + /// which this assignment applies. Might, for instance, be a + /// selection of well segments for a segment level UDQ. + /// + /// \param[in] value_arg Numeric UDQ value for the entities + /// identified in the \p selector. + /// + /// \param[in] report_step_arg Time at which this assignment happens. + /// Assignments should be performed exactly once and the time value + /// ensures this behaviour. AssignRecord(std::vector&& selector, const double value_arg, const std::size_t report_step_arg) @@ -80,82 +356,57 @@ class UDQAssign , report_step (report_step_arg) {} + /// Apply assignment record to existing UDQ set + /// + /// Populates members of UDQ set that are known to the current + /// assignment record. + /// + /// \param[in,out] values UDQ set for which to assign elements. void eval(UDQSet& values) const; + /// Equality predicate + /// + /// \param[in] data Record against which the current object will be + /// compared for equality. + /// + /// \return \code *this == data \endcode. bool operator==(const AssignRecord& data) const; + /// Convert between byte array and object representation. + /// + /// \tparam Serializer Byte array conversion protocol. + /// + /// \param[in,out] serializer Byte array conversion object. template void serializeOp(Serializer& serializer) { serializer(this->input_selector); - serializer(this->rst_selector); serializer(this->numbered_selector); serializer(this->value); serializer(this->report_step); } }; - UDQAssign() = default; - UDQAssign(const std::string& keyword, - const std::vector& selector, - double value, - std::size_t report_step); - - UDQAssign(const std::string& keyword, - const std::unordered_set& selector, - double value, - std::size_t report_step); - - UDQAssign(const std::string& keyword, - const std::vector& selector, - double value, - std::size_t report_step); - - UDQAssign(const std::string& keyword, - std::vector&& selector, - double value, - std::size_t report_step); - - static UDQAssign serializationTestObject(); - - const std::string& keyword() const; - UDQVarType var_type() const; - - void add_record(const std::vector& selector, - double value, - std::size_t report_step); - - void add_record(const std::unordered_set& rst_selector, - double value, - std::size_t report_step); - - void add_record(const std::vector& selector, - double value, - std::size_t report_step); - - void add_record(std::vector&& selector, - double value, - std::size_t report_step); - - UDQSet eval(const std::vector& items) const; - UDQSet eval(const std::vector& wells) const; - UDQSet eval() const; - std::size_t report_step() const; - - bool operator==(const UDQAssign& data) const; - - template - void serializeOp(Serializer& serializer) - { - serializer(m_keyword); - serializer(m_var_type); - serializer(records); - } - -private: + /// Name of UDQ to which this assignment applies. std::string m_keyword{}; + + /// Kind of UDQ to which this assignment applies. UDQVarType m_var_type{UDQVarType::NONE}; + + /// Assignment records for this UDQ assignment. std::vector records{}; + + /// Reconstitute well or group level assignment from restart file + /// information + /// + /// \param[in] assignRst Aggregate UDQ assignment information restored + /// from restart file information. + /// + /// \param[in] report_step Time at which this assignment happens. + /// Assignments should be performed exactly once and the time value + /// ensures this behaviour. + void add_well_or_group_records(const RestartIO::RstUDQ& assignRst, + const std::size_t report_step); }; } // namespace Opm diff --git a/opm/input/eclipse/Schedule/UDQ/UDQConfig.cpp b/opm/input/eclipse/Schedule/UDQ/UDQConfig.cpp index 3f0a9ffcfad..d390113707f 100644 --- a/opm/input/eclipse/Schedule/UDQ/UDQConfig.cpp +++ b/opm/input/eclipse/Schedule/UDQ/UDQConfig.cpp @@ -20,6 +20,7 @@ #include #include +#include #include #include @@ -142,35 +143,23 @@ namespace { namespace Opm { UDQConfig::UDQConfig(const UDQParams& params) - : udq_params(params) - , udqft(this->udq_params) + : udq_params { params } + , udqft { udq_params } {} - UDQConfig::UDQConfig(const UDQParams& params, const RestartIO::RstState& rst_state) - : UDQConfig(params) + UDQConfig::UDQConfig(const UDQParams& params, + const RestartIO::RstState& rst_state) + : UDQConfig { params } { - for (const auto& rst_udq : rst_state.udqs) { - if (rst_udq.is_define()) { - KeywordLocation location("UDQ", "Restart file", 0); - this->add_define(rst_udq.name, - location, - { rst_udq.expression() }, - rst_state.header.report_step); - - auto pos = this->m_definitions.find(rst_udq.name); - assert (pos != this->m_definitions.end()); - - pos->second.update_status(rst_udq.updateStatus(), - rst_state.header.report_step); + for (const auto& udq : rst_state.udqs) { + if (udq.isDefine()) { + this->add_define(udq, rst_state.header.report_step); } else { - this->add_assign(rst_udq.name, - rst_udq.assign_selector(), - rst_udq.assign_value(), - rst_state.header.report_step); + this->add_assign(udq, rst_state.header.report_step); } - this->add_unit(rst_udq.name, rst_udq.unit); + this->add_unit(udq.name, udq.unit); } } @@ -190,89 +179,62 @@ namespace Opm { return result; } - const UDQParams& UDQConfig::params() const + const std::string& UDQConfig::unit(const std::string& key) const { - return this->udq_params; + const auto pair_ptr = this->units.find(key); + if (pair_ptr == this->units.end()) { + throw std::invalid_argument("No such UDQ quantity: " + key); + } + + return pair_ptr->second; } - void UDQConfig::add_node(const std::string& quantity, const UDQAction action) + bool UDQConfig::has_unit(const std::string& keyword) const { - auto index_iter = this->input_index.find(quantity); - if (this->input_index.find(quantity) == this->input_index.end()) { - auto var_type = UDQ::varType(quantity); - auto insert_index = this->input_index.size(); + return this->units.find(keyword) != this->units.end(); + } - this->type_count[var_type] += 1; - this->input_index[quantity] = UDQIndex(insert_index, this->type_count[var_type], action, var_type); - } - else { - index_iter->second.action = action; - } + bool UDQConfig::has_keyword(const std::string& keyword) const + { + return (this->m_assignments.find(keyword) != this->m_assignments.end()) + || (this->m_definitions.find(keyword) != this->m_definitions.end()); } - void UDQConfig::add_assign(const std::string& quantity, - SegmentMatcherFactory create_segment_matcher, - const std::vector& selector, - const double value, - const std::size_t report_step) + void UDQConfig::add_record(SegmentMatcherFactory create_segment_matcher, + const DeckRecord& record, + const KeywordLocation& location, + const std::size_t report_step) { - this->add_node(quantity, UDQAction::ASSIGN); + using KW = ParserKeywords::UDQ; - switch (UDQ::varType(quantity)) { - case UDQVarType::SEGMENT_VAR: - this->add_enumerated_assign(quantity, - std::move(create_segment_matcher), - selector, value, report_step); - break; + const auto action = UDQ::actionType(record.getItem().get(0)); + const auto& quantity = record.getItem().get(0); + const auto data = RawString::strings(record.getItem().getData()); - default: - this->add_named_assign(quantity, selector, value, report_step); - break; + if (action == UDQAction::UPDATE) { + this->add_update(quantity, report_step, location, data); } - - if (this->m_assignments.find(quantity) != this->m_assignments.end()) { - this->pending_assignments_.push_back(quantity); + else if (action == UDQAction::UNITS) { + this->add_unit(quantity, data.front()); } - } - - void UDQConfig::add_assign(const std::string& quantity, - const std::unordered_set& selector, - const double value, - const std::size_t report_step) - { - this->add_node(quantity, UDQAction::ASSIGN); - - auto assignment = this->m_assignments.find(quantity); - if (assignment == this->m_assignments.end()) { - this->m_assignments.emplace(std::piecewise_construct, - std::forward_as_tuple(quantity), - std::forward_as_tuple(quantity, selector, - value, report_step)); + else if (action == UDQAction::ASSIGN) { + auto selector = std::vector(data.begin(), data.end() - 1); + std::transform(selector.cbegin(), selector.cend(), selector.begin(), strip_quotes); + const auto value = std::stod(data.back()); + this->add_assign(quantity, + std::move(create_segment_matcher), + selector, value, report_step); + } + else if (action == UDQAction::DEFINE) { + this->add_define(quantity, location, data, report_step); } else { - assignment->second.add_record(selector, value, report_step); + throw std::runtime_error { + "Unknown UDQ Operation " + std::to_string(static_cast(action)) + }; } } - void UDQConfig::add_define(const std::string& quantity, - const KeywordLocation& location, - const std::vector& expression, - const std::size_t report_step) - { - this->add_node(quantity, UDQAction::DEFINE); - - this->m_definitions.insert_or_assign(quantity, - UDQDefine { - this->udq_params, - quantity, - report_step, - location, - expression - }); - - this->define_order.insert(quantity); - } - void UDQConfig::add_unit(const std::string& keyword, const std::string& quoted_unit) { const auto Unit = strip_quotes(quoted_unit); @@ -318,41 +280,50 @@ namespace Opm { def.update_status(update_status, report_step); } - void UDQConfig::add_record(SegmentMatcherFactory create_segment_matcher, - const DeckRecord& record, - const KeywordLocation& location, - const std::size_t report_step) + void UDQConfig::add_assign(const std::string& quantity, + SegmentMatcherFactory create_segment_matcher, + const std::vector& selector, + const double value, + const std::size_t report_step) { - using KW = ParserKeywords::UDQ; + this->add_node(quantity, UDQAction::ASSIGN); - const auto action = UDQ::actionType(record.getItem().get(0)); - const auto& quantity = record.getItem().get(0); - const auto data = RawString::strings(record.getItem().getData()); + switch (UDQ::varType(quantity)) { + case UDQVarType::SEGMENT_VAR: + this->add_enumerated_assign(quantity, + std::move(create_segment_matcher), + selector, value, report_step); + break; - if (action == UDQAction::UPDATE) { - this->add_update(quantity, report_step, location, data); - } - else if (action == UDQAction::UNITS) { - this->add_unit(quantity, data.front()); - } - else if (action == UDQAction::ASSIGN) { - auto selector = std::vector(data.begin(), data.end() - 1); - std::transform(selector.cbegin(), selector.cend(), selector.begin(), strip_quotes); - const auto value = std::stod(data.back()); - this->add_assign(quantity, - std::move(create_segment_matcher), - selector, value, report_step); - } - else if (action == UDQAction::DEFINE) { - this->add_define(quantity, location, data, report_step); + default: + this->add_assign_impl(quantity, selector, value, report_step); + break; } - else { - throw std::runtime_error { - "Unknown UDQ Operation " + std::to_string(static_cast(action)) - }; + + if (this->m_assignments.find(quantity) != this->m_assignments.end()) { + this->pending_assignments_.push_back(quantity); } } + void UDQConfig::add_define(const std::string& quantity, + const KeywordLocation& location, + const std::vector& expression, + const std::size_t report_step) + { + this->add_node(quantity, UDQAction::DEFINE); + + this->m_definitions.insert_or_assign(quantity, + UDQDefine { + this->udq_params, + quantity, + report_step, + location, + expression + }); + + this->define_order.insert(quantity); + } + void UDQConfig::add_table(const std::string& name, UDT udt) { m_tables.emplace(name, std::move(udt)); @@ -366,9 +337,43 @@ namespace Opm { return update; } - const UDQAssign& UDQConfig::assign(const std::string& key) const + void UDQConfig::eval_assign(const std::size_t report_step, + const Schedule& sched, + const WellMatcher& wm, + SegmentMatcherFactory create_segment_matcher, + SummaryState& st, + UDQState& udq_state) const { - return this->m_assignments.at(key); + auto factories = UDQContext::MatcherFactories{}; + factories.segments = std::move(create_segment_matcher); + + auto context = UDQContext { + this->function_table(), wm, this->m_tables, + std::move(factories), st, udq_state + }; + + this->eval_assign(report_step, sched, context); + } + + void UDQConfig::eval(const std::size_t report_step, + const Schedule& sched, + const WellMatcher& wm, + SegmentMatcherFactory create_segment_matcher, + RegionSetMatcherFactory create_region_matcher, + SummaryState& st, + UDQState& udq_state) const + { + auto factories = UDQContext::MatcherFactories {}; + factories.segments = std::move(create_segment_matcher); + factories.regions = std::move(create_region_matcher); + + auto context = UDQContext { + this->function_table(), wm, this->m_tables, + std::move(factories), st, udq_state + }; + + this->eval_assign(report_step, sched, context); + this->eval_define(report_step, udq_state, context); } const UDQDefine& UDQConfig::define(const std::string& key) const @@ -376,10 +381,9 @@ namespace Opm { return this->m_definitions.at(key); } - UDQAction UDQConfig::action_type(const std::string& udq_key) const + const UDQAssign& UDQConfig::assign(const std::string& key) const { - auto action_iter = this->input_index.find(udq_key); - return action_iter->second.action; + return this->m_assignments.at(key); } std::vector UDQConfig::definitions() const @@ -430,18 +434,6 @@ namespace Opm { return res; } - std::size_t UDQConfig::size() const - { - return std::count_if(this->input_index.begin(), this->input_index.end(), - [](const auto& index_pair) - { - const auto action = index_pair.second.action; - - return (action == UDQAction::DEFINE) - || (action == UDQAction::ASSIGN); - }); - } - void UDQConfig::exportTypeCount(std::array(UDQVarType::NumTypes)>& count) const { count.fill(0); @@ -451,54 +443,16 @@ namespace Opm { } } - std::vector UDQConfig::assignments() const - { - std::vector ret; - - for (const auto& [key, Input] : this->input_index) { - if (Input.action == UDQAction::ASSIGN) { - ret.push_back(this->m_assignments.at(key)); - } - } - - return ret; - } - - std::vector UDQConfig::assignments(const UDQVarType var_type) const - { - std::vector filtered_assigns; - - for (const auto& index_pair : this->input_index) { - auto assign_iter = this->m_assignments.find(index_pair.first); - if ((assign_iter != this->m_assignments.end()) && - (assign_iter->second.var_type() == var_type)) - { - filtered_assigns.push_back(assign_iter->second); - } - } - - return filtered_assigns; - } - - const std::string& UDQConfig::unit(const std::string& key) const - { - const auto pair_ptr = this->units.find(key); - if (pair_ptr == this->units.end()) { - throw std::invalid_argument("No such UDQ quantity: " + key); - } - - return pair_ptr->second; - } - - bool UDQConfig::has_unit(const std::string& keyword) const + std::size_t UDQConfig::size() const { - return this->units.find(keyword) != this->units.end(); - } + return std::count_if(this->input_index.begin(), this->input_index.end(), + [](const auto& index_pair) + { + const auto action = index_pair.second.action; - bool UDQConfig::has_keyword(const std::string& keyword) const - { - return (this->m_assignments.find(keyword) != this->m_assignments.end()) - || (this->m_definitions.find(keyword) != this->m_definitions.end()); + return (action == UDQAction::DEFINE) + || (action == UDQAction::ASSIGN); + }); } UDQInput UDQConfig::operator[](const std::string& keyword) const @@ -550,6 +504,40 @@ namespace Opm { throw std::logic_error("Internal error - should not be here"); } + std::vector UDQConfig::assignments() const + { + std::vector ret; + + for (const auto& [key, Input] : this->input_index) { + if (Input.action == UDQAction::ASSIGN) { + ret.push_back(this->m_assignments.at(key)); + } + } + + return ret; + } + + std::vector UDQConfig::assignments(const UDQVarType var_type) const + { + std::vector filtered_assigns; + + for (const auto& index_pair : this->input_index) { + auto assign_iter = this->m_assignments.find(index_pair.first); + if ((assign_iter != this->m_assignments.end()) && + (assign_iter->second.var_type() == var_type)) + { + filtered_assigns.push_back(assign_iter->second); + } + } + + return filtered_assigns; + } + + const UDQParams& UDQConfig::params() const + { + return this->udq_params; + } + const UDQFunctionTable& UDQConfig::function_table() const { return this->udqft; @@ -574,6 +562,59 @@ namespace Opm { ; } + void UDQConfig::required_summary(std::unordered_set& summary_keys) const + { + for (const auto& def_pair : this->m_definitions) { + def_pair.second.required_summary(summary_keys); + } + } + + // =========================================================================== + // Private member functions below separator + // =========================================================================== + + void UDQConfig::add_node(const std::string& quantity, const UDQAction action) + { + auto index_iter = this->input_index.find(quantity); + if (this->input_index.find(quantity) == this->input_index.end()) { + auto var_type = UDQ::varType(quantity); + auto insert_index = this->input_index.size(); + + this->type_count[var_type] += 1; + this->input_index[quantity] = UDQIndex(insert_index, this->type_count[var_type], action, var_type); + } + else { + index_iter->second.action = action; + } + } + + UDQAction UDQConfig::action_type(const std::string& udq_key) const + { + auto action_iter = this->input_index.find(udq_key); + return action_iter->second.action; + } + + void UDQConfig::add_assign(const RestartIO::RstUDQ& udq, + const std::size_t report_step) + { + this->add_node(udq.name, UDQAction::ASSIGN); + + this->add_assign_impl(udq.name, udq, report_step); + } + + void UDQConfig::add_define(const RestartIO::RstUDQ& udq, + const std::size_t report_step) + { + this->add_define(udq.name, KeywordLocation { "UDQ", "Restart file", 0 }, + std::vector { std::string { udq.definingExpression() } }, + report_step); + + auto pos = this->m_definitions.find(udq.name); + assert (pos != this->m_definitions.end()); + + pos->second.update_status(udq.currentUpdateStatus(), report_step); + } + void UDQConfig::eval_assign(const std::size_t report_step, const Schedule& sched, UDQContext& context) const @@ -655,67 +696,6 @@ namespace Opm { } } - void UDQConfig::eval(const std::size_t report_step, - const Schedule& sched, - const WellMatcher& wm, - SegmentMatcherFactory create_segment_matcher, - RegionSetMatcherFactory create_region_matcher, - SummaryState& st, - UDQState& udq_state) const - { - auto factories = UDQContext::MatcherFactories {}; - factories.segments = std::move(create_segment_matcher); - factories.regions = std::move(create_region_matcher); - - auto context = UDQContext { - this->function_table(), wm, this->m_tables, - std::move(factories), st, udq_state - }; - - this->eval_assign(report_step, sched, context); - this->eval_define(report_step, udq_state, context); - } - - void UDQConfig::eval_assign(const std::size_t report_step, - const Schedule& sched, - const WellMatcher& wm, - SegmentMatcherFactory create_segment_matcher, - SummaryState& st, - UDQState& udq_state) const - { - auto factories = UDQContext::MatcherFactories{}; - factories.segments = std::move(create_segment_matcher); - - auto context = UDQContext { - this->function_table(), wm, this->m_tables, - std::move(factories), st, udq_state - }; - - this->eval_assign(report_step, sched, context); - } - - void UDQConfig::required_summary(std::unordered_set& summary_keys) const - { - for (const auto& def_pair : this->m_definitions) { - def_pair.second.required_summary(summary_keys); - } - } - - void UDQConfig::add_named_assign(const std::string& quantity, - const std::vector& selector, - const double value, - const std::size_t report_step) - { - auto [asgnPos, inserted] = this->m_assignments - .emplace(std::piecewise_construct, - std::forward_as_tuple(quantity), - std::forward_as_tuple(quantity, selector, value, report_step)); - - if (! inserted) { - asgnPos->second.add_record(selector, value, report_step); - } - } - void UDQConfig::add_enumerated_assign(const std::string& quantity, SegmentMatcherFactory create_segment_matcher, const std::vector& selector, @@ -735,14 +715,7 @@ namespace Opm { auto items = UDQSet:: enumerateItems(segmentMatcher->findSegments(setDescriptor)); - auto [asgnPos, inserted] = this->m_assignments - .emplace(std::piecewise_construct, - std::forward_as_tuple(quantity), - std::forward_as_tuple(quantity, items, value, report_step)); - - if (! inserted) { - asgnPos->second.add_record(std::move(items), value, report_step); - } + this->add_assign_impl(quantity, std::move(items), value, report_step); } } // namespace Opm diff --git a/opm/input/eclipse/Schedule/UDQ/UDQConfig.hpp b/opm/input/eclipse/Schedule/UDQ/UDQConfig.hpp index fc522b64d22..b2a4a4c443c 100644 --- a/opm/input/eclipse/Schedule/UDQ/UDQConfig.hpp +++ b/opm/input/eclipse/Schedule/UDQ/UDQConfig.hpp @@ -39,6 +39,7 @@ #include #include #include +#include #include namespace Opm { @@ -54,61 +55,203 @@ namespace Opm { } // namespace Opm -namespace Opm { namespace RestartIO { +namespace Opm::RestartIO { struct RstState; -}} // namespace Opm::RestartIO + class RstUDQ; +} // namespace Opm::RestartIO namespace Opm { + /// Collection of all user-defined quantities in the current simulation run class UDQConfig { public: + /// Factory function for constructing region set matchers. using RegionSetMatcherFactory = std::function()>; + + /// Factory function for constructing segment set matchers. using SegmentMatcherFactory = std::function()>; + /// Default constructor UDQConfig() = default; + + /// Constructor + /// + /// Main constructor for a base run. + /// + /// \param[in] params UDQ parameters from UDQPARAM keyword. explicit UDQConfig(const UDQParams& params); + + /// Constructor + /// + /// Main constructor for a restarted simulation run. + /// + /// \param[in] params UDQ parameters from UDQPARAM keyword. + /// + /// \param[in] rst_state Object state from restart file information. UDQConfig(const UDQParams& params, const RestartIO::RstState& rst_state); + /// Create a serialisation test object. static UDQConfig serializationTestObject(); + /// Retrieve unit string for a particular UDQ + /// + /// \param[in] key Unqualified UDQ name such as FUNNY, GUITAR, + /// WURST, or SUSHI. + /// + /// \return Input unit string associated to \p key. Throws an + /// exception of type \code std::invalid_argument \endcode if no + /// unit string exists for \p key. const std::string& unit(const std::string& key) const; + + /// Query whether or not a particular UDQ has an associated unit + /// string. + /// + /// \param[in] keyword Unqualified UDQ name such as FUNNY, GUITAR, + /// WURST, or SUSHI. + /// + /// \return Whether or not there is a unit string associated to \p + /// keyword. bool has_unit(const std::string& keyword) const; + + /// Query whether or not a particular UDQ exists in the collection. + /// + /// \param[in] keyword Unqualified UDQ name such as FUNNY, GUITAR, + /// WURST, or SUSHI. + /// + /// \return whether or not \p keyword exists in the collection. bool has_keyword(const std::string& keyword) const; + /// Incorporate a single UDQ record into the known collection + /// + /// \param[in] create_segment_matcher Factory function for + /// constructing segment set matchers. + /// + /// \param[in] record UDQ keyword record, such as a DEFINE, ASSIGN, + /// UPDATE, or UNIT statement. + /// + /// \param[in] location Input file/line information for the UDQ + /// record. Mostly for diagnostic purposes. + /// + /// \param[in] report_step Time at which this record is encountered. void add_record(SegmentMatcherFactory create_segment_matcher, const DeckRecord& record, const KeywordLocation& location, std::size_t report_step); + /// Incorporate a unit string for a UDQ + /// + /// Implements the UNIT statement. + /// + /// \param[in] keyword Unqualified UDQ name such as FUNNY, GUITAR, + /// WURST, or SUSHI. + /// + /// \param[in] unit Unit string for UDQ \p keyword. void add_unit(const std::string& keyword, const std::string& unit); + /// Incorporate update status change for a UDQ + /// + /// Implements the UPDATE statement + /// + /// \param[in] keyword Unqualified UDQ name such as FUNNY, GUITAR, + /// WURST, or SUSHI. + /// + /// \param[in] report_step Time at which this record is encountered. + /// + /// \param[in] location Input file/line information for the UDQ + /// record. Mostly for diagnostic purposes. + /// + /// \param[in] data Update status. Should be a single element + /// vector containing one of the status strings ON, OFF, or NEXT. void add_update(const std::string& keyword, std::size_t report_step, const KeywordLocation& location, const std::vector& data); + /// Incorporate a UDQ assignment. + /// + /// Implements the ASSIGN statement. + /// + /// \param[in] quantity Unqualified UDQ name such as FUNNY, GUITAR, + /// WURST, or SUSHI. + /// + /// \param[in] create_segment_matcher Factory function for + /// constructing segment set matchers. + /// + /// \param[in] selector UDQ set element selection to which this + /// assignment applies. Might be a well name pattern if \p quantity + /// is a well level UDQ. + /// + /// \param[in] report_step Time at which this assignment statement + /// is encountered. void add_assign(const std::string& quantity, SegmentMatcherFactory create_segment_matcher, const std::vector& selector, double value, std::size_t report_step); - void add_assign(const std::string& quantity, - const std::unordered_set& selector, - double value, - std::size_t report_step); - + /// Incorporate a UDQ defining expressions. + /// + /// Implements the DEFINE statement. + /// + /// \param[in] quantity Unqualified UDQ name such as FUNNY, GUITAR, + /// WURST, or SUSHI. + /// + /// \param[in] location Input file/line information for the UDQ + /// record. Mostly for diagnostic purposes. + /// + /// \param[in] expression Defining expression for the UDQ. Function + /// add_define() parses this expression into an abstract syntax tree + /// which is used in subsequent evaluation contexts. + /// + /// \param[in] report_step Time at which this assignment statement + /// is encountered. void add_define(const std::string& quantity, const KeywordLocation& location, const std::vector& expression, std::size_t report_step); + /// Incorporate a user defined table. + /// + /// Implements the UDT keyword. + /// + /// \param[in] name Name of user defined table. + /// + /// \param[in] udt Tabulated values. void add_table(const std::string& name, UDT udt); + /// Clear all pending assignments + /// + /// Clears all internal data structures of any assignment records. + /// Typically called at the end of a report step in order to signify + /// that all assignments have been applied. + /// + /// \return Whether or not there were any active assignments in the + /// internal representation. Allows client code to take action, if + /// needed, based on the knowledge that all assignments have been + /// applied and to prepare for the next report step. bool clear_pending_assignments(); + /// Apply all pending assignments. + /// + /// Assigns new UDQ values to both the summary and UDQ state objects. + /// + /// \param[in] report_step Current report step. + /// + /// \param[in] sched Full dynamic input schedule. + /// + /// \param[in] wm Well name pattern matcher. + /// + /// \param[in] create_segment_matcher Factory function for + /// constructing segment set matchers. + /// + /// \param[in,out] st Summary vectors. For output and evaulating + /// ACTION condition purposes. Values pertaining to UDQs being + /// assigned here will be updated. + /// + /// \param[in,out] udq_state Dynamic values for all known UDQs. + /// Values pertaining to UDQs being assigned here will be updated. void eval_assign(std::size_t report_step, const Schedule& sched, const WellMatcher& wm, @@ -116,6 +259,30 @@ namespace Opm { SummaryState& st, UDQState& udq_state) const; + /// Compute new values for all UDQs + /// + /// Uses both assignment and defining expressions as applicable. + /// Assigns new UDQ values to both the summary and UDQ state + /// objects. + /// + /// \param[in] report_step Current report step. + /// + /// \param[in] sched Full dynamic input schedule. + /// + /// \param[in] wm Well name pattern matcher. + /// + /// \param[in] create_segment_matcher Factory function for + /// constructing segment set matchers. + /// + /// \param[in] create_region_matcher Factory function for + /// constructing region set matchers. + /// + /// \param[in,out] st Summary vectors. For output and evaulating + /// ACTION condition purposes. Values pertaining to UDQs being + /// assigned here will be updated. + /// + /// \param[in,out] udq_state Dynamic values for all known UDQs. + /// Values pertaining to UDQs being assigned here will be updated. void eval(std::size_t report_step, const Schedule& sched, const WellMatcher& wm, @@ -124,32 +291,119 @@ namespace Opm { SummaryState& st, UDQState& udq_state) const; + /// Retrieve defining expression and evaluation object for a single + /// UDQ + /// + /// \param[in] key Unqualified UDQ name such as FUNNY, GUITAR, + /// WURST, or SUSHI. + /// + /// \return Defining expression and evaluation object for \p key. + /// Throws an exception of type \code std::out_of_range \endcode if + /// no such object exists for the UDQ \p key. const UDQDefine& define(const std::string& key) const; + + /// Retrieve any pending assignment object for a single UDQ + /// + /// \param[in] key Unqualified UDQ name such as FUNNY, GUITAR, + /// WURST, or SUSHI. + /// + /// \return Pending assingment object for \p key. Throws an + /// exception of type \code std::out_of_range \endcode if no such + /// object exists for the UDQ \p key. const UDQAssign& assign(const std::string& key) const; + /// Retrieve defining expressions and evaluation objects for all + /// known UDQs. std::vector definitions() const; + + /// Retrieve defining expressions and evaluation objects for all + /// known UDQs of a particular category. + /// + /// \param[in] var_type UDQ category. + /// + /// \return Defining expressions and evaluation objects for all + /// known UDQs of category \p var_type. std::vector definitions(UDQVarType var_type) const; + + /// Retrieve unprocessed input objects for all UDQs + /// + /// Needed for restart file output purposes. std::vector input() const; + /// Export count of all known UDQ categories in the current run. + /// + /// \param[out] count Count of all active UDQs of all categories in + /// the current run. void exportTypeCount(std::array(UDQVarType::NumTypes)>& count) const; - // The size() method will return the number of active DEFINE and ASSIGN - // statements; this will correspond to the length of the vector returned - // from input(). + /// Total number of active DEFINE and ASSIGN statements. + /// + /// Corresponds to the length of the vector returned from input(). std::size_t size() const; + /// Unprocessed input object for named quantity. + /// + /// \param[in] keyword Unqualified UDQ name such as FUNNY, GUITAR, + /// WURST, or SUSHI. + /// + /// \return Unprocessed input object for \p keyword. Throws an + /// exception of type \code std::invalid_argument \endcode if no + /// such named UDQ exists. UDQInput operator[](const std::string& keyword) const; + + /// Unprocessed input object for enumerated quantity. + /// + /// \param[in] insert_index Linear index in order of appearance for + /// an individual UDQ. + /// + /// \return Unprocessed input object for \p keyword. Throws an + /// exception of type \code std::invalid_argument \endcode if no + /// such numbered UDQ exists. UDQInput operator[](std::size_t insert_index) const; + /// Retrieve pending assignment objects for all known UDQs. std::vector assignments() const; + + /// Retrieve pending assignment objects for all known UDQs of a + /// particular category. + /// + /// \param[in] var_type UDQ category. + /// + /// \return Pending assignment objects for all known UDQs of + /// category \p var_type. std::vector assignments(UDQVarType var_type) const; + + /// Retrieve run's active UDQ parameters const UDQParams& params() const; + + /// Retrieve run's active UDQ function table. const UDQFunctionTable& function_table() const; + + /// Retrieve run's active user defined tables. const std::unordered_map& tables() const; + /// Equality predicate. + /// + /// \param[in] config Object against which \code *this \endcode will + /// be tested for equality. + /// + /// \return Whether or not \code *this \endcode is the same as \p config. bool operator==(const UDQConfig& config) const; + + /// Export all summary vectors needed to compute values for the + /// current collection of user defined quantities. + /// + /// \param[in,out] summary_keys Named summary vectors. Upon + /// completion, any additional summary vectors needed to evaluate + /// the current set of user defined quantities will be included in + /// this set. void required_summary(std::unordered_set& summary_keys) const; + /// Convert between byte array and object representation. + /// + /// \tparam Serializer Byte array conversion protocol. + /// + /// \param[in,out] serializer Byte array conversion object. template void serializeOp(Serializer& serializer) { @@ -170,48 +424,177 @@ namespace Opm { } private: + /// Run's active UDQ parameters. + /// + /// Initialised from the run's UDQPARAM keyword. UDQParams udq_params; - UDQFunctionTable udqft; - // The choices of data structures are strongly motivated by the - // constraints imposed by the ECLIPSE formatted restart files; for - // writing restart files it is essential to keep meticulous control - // over the ordering of the keywords. In this class the ordering is - // mainly maintained by the input_index map which keeps track of the - // insert order of each keyword, and whether the keyword is - // currently DEFINE'ed or ASSIGN'ed. - std::unordered_map m_definitions; - std::unordered_map m_assignments; - std::unordered_map m_tables; - std::unordered_map units; + /// Run's active function table table. + UDQFunctionTable udqft; + // The following data structures are constrained by compatibility + // requirements in our simulation restart files. In particuar we + // need to control the keyword ordering. In this class the ordering + // is maintained mainly by the input_index map which tracks the + // insertion order of each keyword and whether the keyword (UDQ + // name) is currently DEFINE'ed or ASSIGN'ed. + + /// Defining expressions and evaluation objects for all pertinent + /// quantities. + /// + /// Keyed by UDQ name. + std::unordered_map m_definitions{}; + + /// Run's UDQ assignment statements. + /// + /// Keyed by UDQ name. + std::unordered_map m_assignments{}; + + /// Run's user defined input tables. + /// + /// Keyed by table name (i.e., TU* strings). + std::unordered_map m_tables{}; + + /// Unit strings for some or all input UDQs. + /// + /// Defined only for those UDQs which have an associate UNIT + /// statement. + /// + /// Keyed by UDQ name. + std::unordered_map units{}; + + /// Ordered set of DEFINE statements. + /// + /// Mostly unused and should probably be removed. IOrderSet define_order; + + /// Ordered set of UDQ inputs. OrderedMap input_index; - std::map type_count; + /// Number of UDQs of each category currently active. + std::map type_count{}; + + /// List of pending assignment statements. + /// + /// Marked mutable because this will be modified in member function + /// + /// UDQConfig::eval_assign(step, sched, context) const mutable std::vector pending_assignments_{}; + /// Incorporate operation for new or existing UDQ + /// + /// Preserves order of operations in input_index. + /// + /// \param[in] quantity Unqualified UDQ name such as FUNNY, GUITAR, + /// WURST, or SUSHI. + /// + /// \param[in] action Kind of operation. void add_node(const std::string& quantity, UDQAction action); + + /// Retrieve operation currently associated a UDQ + /// + /// \param[in] udq_key Unqualified UDQ name such as FUNNY, GUITAR, + /// WURST, or SUSHI. Shall be the name of an existing UDQ since + /// this function otherwise will dereference \code input_index.end() + /// \endcode. + /// + /// \return Operation currently associated to \p udq_key. UDQAction action_type(const std::string& udq_key) const; + /// Reconstitute an assignment statement from restart file information. + /// + /// \param[in] udq Restart file representation of an assignment statement. + /// + /// \param[in] report_step Time at which this assignment is + /// encountered. Typically corresponds to the restart time. + void add_assign(const RestartIO::RstUDQ& udq, const std::size_t report_step); + + /// Reconstitute a definition statement from restart file information. + /// + /// \param[in] udq Restart file representation of a UDQ definition + /// statement. + /// + /// \param[in] report_step Time at which this definition is + /// encountered. Typically corresponds to the restart time. + void add_define(const RestartIO::RstUDQ& udq, const std::size_t report_step); + + /// Apply all pending assignments. + /// + /// Assigns new UDQ values to both the summary and UDQ state objects. + /// + /// \param[in] report_step Current report step. + /// + /// \param[in] schedule Full dynamic input schedule. + /// + /// \param[in,out] context Pattern matchers and state objects. + /// Values pertaining to UDQs being assigned here will be updated. void eval_assign(std::size_t report_step, const Schedule& sched, UDQContext& context) const; - void eval_define(std::size_t report_step, + /// Compute new values for all UDQs + /// + /// Evaluates all applicable defining expressions. Assigns new UDQ + /// values to both the summary and UDQ state objects. + /// + /// \param[in] report_step Current report step. + /// + /// \param[in] schedule Full dynamic input schedule. + /// + /// \param[in,out] context Pattern matchers and state objects. + /// Values pertaining to UDQs being evaluated here will be updated. + void eval_define(std::size_t report_step, const UDQState& udq_state, - UDQContext& context) const; - - void add_named_assign(const std::string& quantity, - const std::vector& selector, - double value, - std::size_t report_step); + UDQContext& context) const; + /// Incorporate an enumerated assignment statement into known UDQ + /// collection. + /// + /// Typically assigns a segment level UDQ for one or more segments + /// in a single multi-segmented well. + /// + /// \param[in] quantity Unqualified UDQ name. Typically names a + /// segment level UDQ such as SUSHI. + /// + /// \param[in] create_segment_matcher Factory function for + /// constructing segment set matchers. + /// + /// \param[in] selector UDQ set element selection to which this + /// assignment applies. Might be a well name and a segment number. + /// + /// \param[in] report_step Time at which this assignment statement + /// is encountered. void add_enumerated_assign(const std::string& quantity, SegmentMatcherFactory create_segment_matcher, const std::vector& selector, double value, std::size_t report_step); + + /// Common implementation of all add*assign() member functions. + /// + /// \tparam Args Parameter pack type for UDQAssign::add_record() + /// + /// \param[in] quantity Unqualified UDQ name such as FUNNY, GUITAR, + /// WURST, or SUSHI. + /// + /// \param[in] args Constructor or add_record() arguments for + /// UDQAssign object. + template + void add_assign_impl(const std::string& quantity, + Args&&... args) + { + // Note: Duplicate 'quantity' arguments is intentional here. + // The first is the key in 'm_assignments', while the second is + // the first argument to the UDQAssign constructor. + const auto& [asgnPos, inserted] = + this->m_assignments.try_emplace(quantity, quantity, args...); + + if (! inserted) { + // We already have an assignment object for 'quantity'. Add + // a new assignment record to that object. + asgnPos->second.add_record(std::forward(args)...); + } + } }; } // namespace Opm diff --git a/opm/input/eclipse/Schedule/UDQ/UDQState.cpp b/opm/input/eclipse/Schedule/UDQ/UDQState.cpp index c21f4918aeb..8d1cebebc64 100644 --- a/opm/input/eclipse/Schedule/UDQ/UDQState.cpp +++ b/opm/input/eclipse/Schedule/UDQ/UDQState.cpp @@ -19,13 +19,15 @@ #include +#include + #include #include #include #include -#include #include +#include #include @@ -129,6 +131,22 @@ void add_results(const std::string& udq_key, } } +// Load restart values for UDQs defined at the group or well levels. +void load_restart_values(const Opm::RestartIO::RstUDQ& udq, + SMap& values) +{ + const auto& wgnames = udq.entityNames(); + const auto& nameIdx = udq.nameIndex(); + const auto n = udq.numEntities(); + + for (auto i = 0*n; i < n; ++i) { + for (const auto& subValuePair : udq[i]) { + values.insert_or_assign(wgnames[nameIdx[i]], + subValuePair.second); + } + } +} + double get_scalar(const SMap& values, const std::string& udq_key, const double undef_value) @@ -166,47 +184,30 @@ namespace Opm { void UDQState::load_rst(const RestartIO::RstState& rst_state) { for (const auto& udq : rst_state.udqs) { - if (udq.is_define()) { - if (udq.var_type == UDQVarType::WELL_VAR) { - auto& values = this->well_values[udq.name]; - for (const auto& [wname, value] : udq.values()) { - values.insert_or_assign(wname, value); - } + // Note: Cases listed in order of increasing enumerator values from + // the UDQEnums.hpp header file (Opm::UDQVarType). + switch (udq.category) { + case UDQVarType::SCALAR: + case UDQVarType::FIELD_VAR: + if (udq.isScalar()) { + // There is a well defined scalar value in the 'udq' object + // for this scalar or field-level UDQ. + this->scalar_values.insert_or_assign(udq.name, udq.scalarValue()); } + break; - if (udq.var_type == UDQVarType::GROUP_VAR) { - auto& values = this->group_values[udq.name]; - for (const auto& [gname, value] : udq.values()) { - values.insert_or_assign(gname, value); - } - } - if (const auto& field_value = udq.field_value(); field_value.has_value()) { - this->scalar_values[udq.name] = field_value.value(); - } - } - else { - const auto value = udq.assign_value(); - - if ((udq.var_type == UDQVarType::WELL_VAR) && - ! udq.assign_selector().empty()) - { - auto& values = this->well_values[udq.name]; - for (const auto& wname : udq.assign_selector()) { - values.insert_or_assign(wname, value); - } - } + case UDQVarType::WELL_VAR: + load_restart_values(udq, this->well_values[udq.name]); + break; - if (udq.var_type == UDQVarType::GROUP_VAR) { - auto& values = this->group_values[udq.name]; - for (const auto& gname : udq.assign_selector()) { - values.insert_or_assign(gname, value); - } - } + case UDQVarType::GROUP_VAR: + load_restart_values(udq, this->group_values[udq.name]); + break; - if (udq.var_type == UDQVarType::FIELD_VAR) { - this->scalar_values[udq.name] = value; - } + default: + // Not currently supported + break; } } } diff --git a/opm/io/eclipse/rst/state.cpp b/opm/io/eclipse/rst/state.cpp index 67cb13f1870..f476d632e73 100644 --- a/opm/io/eclipse/rst/state.cpp +++ b/opm/io/eclipse/rst/state.cpp @@ -34,24 +34,40 @@ #include #include -#include -#include - -#include -#include - #include + #include #include #include #include + #include +#include +#include + +#include + +#include + #include +#include +#include +#include #include #include #include #include +#include +#include +#include +#include + +#include + +#include + +namespace VI = ::Opm::RestartIO::Helpers::VectorItems; namespace { std::string @@ -60,12 +76,12 @@ namespace { { auto zudl_begin = zudl.begin(); auto zudl_end = zudl_begin; - std::advance( zudl_begin, (udq_index + 0) * Opm::UDQDims::entriesPerZUDL() ); - std::advance( zudl_end , (udq_index + 1) * Opm::UDQDims::entriesPerZUDL() ); + std::advance(zudl_begin, (udq_index + 0) * Opm::UDQDims::entriesPerZUDL()); + std::advance(zudl_end , (udq_index + 1) * Opm::UDQDims::entriesPerZUDL()); - auto define = std::accumulate(zudl_begin, zudl_end, std::string{}, std::plus<>()); - if (!define.empty() && (define[0] == '~')) { - define[0] = '-'; + auto define = fmt::format("{}", fmt::join(zudl_begin, zudl_end, "")); + if (!define.empty() && (define.front() == '~')) { + define.front() = '-'; } return define; @@ -76,11 +92,165 @@ namespace { { return Opm::UDQ::updateType(iudq[udq_index * Opm::UDQDims::entriesPerIUDQ()]); } -} -namespace VI = ::Opm::RestartIO::Helpers::VectorItems; + template + boost::iterator_range::const_iterator> + getDataWindow(const std::vector& arr, + const std::size_t windowSize, + const std::size_t entity, + const std::size_t subEntity = 0, + const std::size_t maxSubEntitiesPerEntity = 1) + { + const auto off = + windowSize * (subEntity + maxSubEntitiesPerEntity*entity); + + auto begin = arr.begin() + off; + auto end = begin + windowSize; + + return { begin, end }; + } + + class UDQVectors + { + public: + UDQVectors(std::shared_ptr rst_view) + : rstView_{ std::move(rst_view) } + { + const auto& intehead = this->rstView_->getKeyword("INTEHEAD"); + + this->maxNumMsWells_ = intehead[VI::intehead::NSWLMX]; + this->numGroups_ = intehead[VI::intehead::NGMAXZ]; // Including FIELD + this->numWells_ = intehead[VI::intehead::NWMAXZ]; + } + + enum Type : std::size_t { + Field, Group, Segment, Well, + // ---- Implementation Helper ---- + NumTypes + }; + + void prepareNext(const Type t) { ++this->varIx_[t]; } + + const std::vector& zudn() const + { + return this->rstView_->getKeyword("ZUDN"); + } + + bool hasGroup() const { return this->rstView_->hasKeyword("DUDG"); } + bool hasWell() const { return this->rstView_->hasKeyword("DUDW"); } + + double currentFieldUDQValue() const + { + return this->rstView_->getKeyword("DUDF")[ this->varIx_[Type::Field] ]; + } + + auto currentGroupUDQValue() const + { + return getDataWindow(this->rstView_->getKeyword("DUDG"), + this->numGroups_, this->varIx_[Type::Group]); + } + + auto currentWellUDQValue() const + { + return getDataWindow(this->rstView_->getKeyword("DUDW"), + this->numWells_, this->varIx_[Type::Well]); + } + + private: + std::shared_ptr rstView_; + + std::size_t maxNumMsWells_{0}; + std::size_t numGroups_{0}; + std::size_t numWells_{0}; + + std::array varIx_{}; + }; + + bool isDefaultedUDQ(const double x) + { + return x == Opm::UDQ::restart_default; + } + + void restoreFieldUDQValue(const UDQVectors& udqs, + Opm::RestartIO::RstUDQ& udq) + { + const auto dudf = udqs.currentFieldUDQValue(); + + if (! isDefaultedUDQ(dudf)) { + udq.assignScalarValue(dudf); + } + } -namespace Opm { namespace RestartIO { + void restoreGroupUDQValue(const UDQVectors& udqs, + const std::vector& groups, + Opm::RestartIO::RstUDQ& udq) + { + const auto nGrp = groups.size(); + const auto dudg = udqs.currentGroupUDQValue(); + + const auto subEntity = 0; + auto entity = 0; + for (auto iGrp = 0*nGrp; iGrp < nGrp; ++iGrp) { + if (isDefaultedUDQ(dudg[iGrp])) { + continue; + } + + udq.addValue(entity++, subEntity, dudg[iGrp]); + udq.addEntityName(groups[iGrp].name); + } + } + + void restoreWellUDQValue(const UDQVectors& udqs, + const std::vector& wells, + Opm::RestartIO::RstUDQ& udq) + { + const auto nWell = wells.size(); + const auto dudw = udqs.currentWellUDQValue(); + + const auto subEntity = 0; + auto entity = 0; + for (auto iWell = 0*nWell; iWell < nWell; ++iWell) { + if (isDefaultedUDQ(dudw[iWell])) { + continue; + } + + udq.addValue(entity++, subEntity, dudw[iWell]); + udq.addEntityName(wells[iWell].name); + } + } + + void restoreSingleUDQ(const Opm::RestartIO::RstState& rst, + UDQVectors& udqValues, + Opm::RestartIO::RstUDQ& udq) + { + udq.prepareValues(); + + // Note: Categories ordered by enumerator values in UDQEnums.hpp. + switch (udq.category) { + case Opm::UDQVarType::FIELD_VAR: + restoreFieldUDQValue(udqValues, udq); + udqValues.prepareNext(UDQVectors::Type::Field); + break; + + case Opm::UDQVarType::WELL_VAR: + restoreWellUDQValue(udqValues, rst.wells, udq); + udqValues.prepareNext(UDQVectors::Type::Well); + break; + + case Opm::UDQVarType::GROUP_VAR: + restoreGroupUDQValue(udqValues, rst.groups, udq); + udqValues.prepareNext(UDQVectors::Type::Group); + break; + + default: + break; + } + + udq.commitValues(); + } +} + +namespace Opm::RestartIO { RstState::RstState(std::shared_ptr rstView, const Runspec& runspec, @@ -186,7 +356,8 @@ void RstState::add_wells(const std::vector& zwel, const std::vector& xwel, const std::vector& icon, const std::vector& scon, - const std::vector& xcon) { + const std::vector& xcon) +{ for (int iw = 0; iw < this->header.num_wells; iw++) { std::size_t zwel_offset = iw * this->header.nzwelz; @@ -223,7 +394,8 @@ void RstState::add_msw(const std::vector& zwel, const std::vector& scon, const std::vector& xcon, const std::vector& iseg, - const std::vector& rseg) { + const std::vector& rseg) +{ for (int iw = 0; iw < this->header.num_wells; iw++) { std::size_t zwel_offset = iw * this->header.nzwelz; @@ -251,59 +423,35 @@ void RstState::add_msw(const std::vector& zwel, } } -void RstState::add_udqs(const std::vector& iudq, - const std::vector& zudn, - const std::vector& zudl, - const std::vector& dudw, - const std::vector& dudg, - const std::vector& dudf) +void RstState::add_udqs(std::shared_ptr rstView) { - auto well_var = 0*this->header.num_wells; - auto group_var = 0*this->header.ngroup; - auto field_var = 0*this->header.num_udq(); + const auto& iudq = rstView->getKeyword("IUDQ"); + const auto& zudn = rstView->getKeyword("ZUDN"); + const auto& zudl = rstView->getKeyword("ZUDL"); + + if (rstView->hasKeyword("IUAD")) { + const auto& iuad = rstView->getKeyword("IUAD"); + const auto& iuap = rstView->getKeyword("IUAP"); + const auto& igph = rstView->getKeyword("IGPH"); - for (auto udq_index = 0*this->header.num_udq(); udq_index < this->header.num_udq(); ++udq_index) { + this->udq_active = RstUDQActive(iuad, iuap, igph); + } + + auto udqValues = UDQVectors { std::move(rstView) }; + + for (auto udq_index = 0*this->header.num_udq(); + udq_index < this->header.num_udq(); ++udq_index) + { const auto& name = zudn[udq_index*UDQDims::entriesPerZUDN() + 0]; const auto& unit = zudn[udq_index*UDQDims::entriesPerZUDN() + 1]; const auto define = udq_define(zudl, udq_index); + auto& udq = define.empty() ? this->udqs.emplace_back(name, unit) : this->udqs.emplace_back(name, unit, define, udq_update(iudq, udq_index)); - if (udq.var_type == UDQVarType::WELL_VAR) { - for (std::size_t well_index = 0; well_index < this->wells.size(); well_index++) { - auto well_value = dudw[ well_var * this->header.max_wells_in_field + well_index ]; - if (well_value == UDQ::restart_default) - continue; - - const auto& well_name = this->wells[well_index].name; - udq.add_value(well_name, well_value); - } - - ++well_var; - } - - if (udq.var_type == UDQVarType::GROUP_VAR) { - for (std::size_t group_index = 0; group_index < this->groups.size(); group_index++) { - auto group_value = dudg[ group_var * this->header.max_groups_in_field + group_index ]; - if (group_value == UDQ::restart_default) - continue; - - const auto& group_name = this->groups[group_index].name; - udq.add_value(group_name, group_value); - } - - ++group_var; - } - - if (udq.var_type == UDQVarType::FIELD_VAR) { - auto field_value = dudf[ field_var ]; - if (field_value != UDQ::restart_default) - udq.add_value(field_value); - - ++field_var; - } + restoreSingleUDQ(*this, udqValues, udq); } } @@ -402,7 +550,8 @@ void RstState::add_wlist(const std::vector& zwls, } } -const RstWell& RstState::get_well(const std::string& wname) const { +const RstWell& RstState::get_well(const std::string& wname) const +{ const auto well_iter = std::find_if(this->wells.begin(), this->wells.end(), [&wname] (const auto& well) { @@ -428,6 +577,7 @@ RstState RstState::load(std::shared_ptr rstView, const auto& igrp = rstView->getKeyword("IGRP"); const auto& sgrp = rstView->getKeyword("SGRP"); const auto& xgrp = rstView->getKeyword("XGRP"); + state.add_groups(zgrp, igrp, sgrp, xgrp); } @@ -442,15 +592,19 @@ RstState RstState::load(std::shared_ptr rstView, const auto& xcon = rstView->getKeyword("XCON"); if (rstView->hasKeyword("ISEG")) { + // Multi-segmented wells in restart file. const auto& iseg = rstView->getKeyword("ISEG"); const auto& rseg = rstView->getKeyword("RSEG"); state.add_msw(zwel, iwel, swel, xwel, icon, scon, xcon, iseg, rseg); - } else + } + else { + // Standard wells only. state.add_wells(zwel, iwel, swel, xwel, icon, scon, xcon); + } if (rstView->hasKeyword("IWLS")) { const auto& iwls = rstView->getKeyword("IWLS"); @@ -461,22 +615,7 @@ RstState RstState::load(std::shared_ptr rstView, } if (state.header.num_udq() > 0) { - const auto& iudq = rstView->getKeyword("IUDQ"); - const auto& zudn = rstView->getKeyword("ZUDN"); - const auto& zudl = rstView->getKeyword("ZUDL"); - - const auto& dudw = state.header.nwell_udq > 0 ? rstView->getKeyword("DUDW") : std::vector{}; - const auto& dudg = state.header.ngroup_udq > 0 ? rstView->getKeyword("DUDG") : std::vector{}; - const auto& dudf = state.header.nfield_udq > 0 ? rstView->getKeyword("DUDF") : std::vector{}; - - state.add_udqs(iudq, zudn, zudl, dudw, dudg, dudf); - - if (rstView->hasKeyword("IUAD")) { - const auto& iuad = rstView->getKeyword("IUAD"); - const auto& iuap = rstView->getKeyword("IUAP"); - const auto& igph = rstView->getKeyword("IGPH"); - state.udq_active = RstUDQActive(iuad, iuap, igph); - } + state.add_udqs(rstView); } if (state.header.num_action > 0) { @@ -490,8 +629,7 @@ RstState RstState::load(std::shared_ptr rstView, state.add_actions(parser, runspec, state.header.sim_time(), zact, iact, sact, zacn, iacn, sacn, zlact); } - return state; } -}} // namespace Opm::RestartIO +} // namespace Opm::RestartIO diff --git a/opm/io/eclipse/rst/state.hpp b/opm/io/eclipse/rst/state.hpp index 0bbd7b1658f..c8ad5423940 100644 --- a/opm/io/eclipse/rst/state.hpp +++ b/opm/io/eclipse/rst/state.hpp @@ -34,6 +34,7 @@ #include +#include #include #include #include @@ -108,12 +109,7 @@ struct RstState const std::vector& iseg, const std::vector& rseg); - void add_udqs(const std::vector& iudq, - const std::vector& zudn, - const std::vector& zudl, - const std::vector& dudw, - const std::vector& dudg, - const std::vector& dudf); + void add_udqs(std::shared_ptr rstView); void add_actions(const Parser& parser, const Runspec& runspec, diff --git a/opm/io/eclipse/rst/udq.cpp b/opm/io/eclipse/rst/udq.cpp index 46f5eec0105..3e78c970e8a 100644 --- a/opm/io/eclipse/rst/udq.cpp +++ b/opm/io/eclipse/rst/udq.cpp @@ -3,7 +3,8 @@ This file is part of the Open Porous Media project (OPM). - OPM is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. @@ -16,137 +17,279 @@ along with OPM. If not, see . */ -#include - #include + #include -namespace Opm { -namespace RestartIO { +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include -RstUDQ::RstDefine::RstDefine(const std::string& expression_arg, UDQUpdate status_arg) : - expression(expression_arg), - status(status_arg) -{} +#include +namespace detail { + template + struct Overloaded : public Ts... { using Ts::operator()...; }; -void RstUDQ::RstAssign::update_value(const std::string& name_arg, double new_value) { - auto current_value = this->value.value_or(new_value); - if (current_value != new_value) - throw std::logic_error(fmt::format("Internal error: the UDQ {} changes value {} -> {} during restart load", name_arg, current_value, new_value)); + template + Overloaded(Ts...) -> Overloaded; +} + +void Opm::RestartIO::RstUDQ::ValueRange::Iterator::setValue() +{ + this->deref_value_.first = this->i_[this->ix_]; - this->value = new_value; + this->deref_value_.second = + std::visit(detail::Overloaded { + [] (const double v) { return v; }, + [ix = this->ix_](const double* v) { return v[ix]; } + }, this->value_); } +// --------------------------------------------------------------------------- -RstUDQ::RstUDQ(const std::string& name_arg, const std::string& unit_arg, const std::string& define_arg, UDQUpdate update_arg) - : name(name_arg) - , unit(unit_arg) - , var_type(UDQ::varType(name_arg)) - , data(RstDefine{define_arg,update_arg}) +Opm::RestartIO::RstUDQ::RstUDQ(const std::string& name_arg, + const std::string& unit_arg, + const std::string& define_arg, + const UDQUpdate update_arg) + : name { name_arg } + , unit { unit_arg } + , category { UDQ::varType(name_arg) } + , definition_{ std::in_place, define_arg, update_arg } +{} + +Opm::RestartIO::RstUDQ::RstUDQ(const std::string& name_arg, + const std::string& unit_arg) + : name { name_arg } + , unit { unit_arg } + , category { UDQ::varType(name_arg) } +{} + +void Opm::RestartIO::RstUDQ::prepareValues() { + this->entityMap_.clear(); } -RstUDQ::RstUDQ(const std::string& name_arg, const std::string& unit_arg) - : name(name_arg) - , unit(unit_arg) - , var_type(UDQ::varType(name_arg)) - , data(RstAssign{}) +void Opm::RestartIO::RstUDQ:: +addValue(const int entity, const int subEntity, const double value) { -} + if (this->isScalar()) { + throw std::logic_error { + fmt::format("UDQ {} cannot be defined as a scalar " + "and then used as UDQ set at restart time", + this->name) + }; + } + + if (std::holds_alternative(this->sa_)) { + this->sa_.emplace>(); + } + + this->entityMap_.addConnection(entity, subEntity); + std::get>(this->sa_).push_back(value); -bool RstUDQ::is_define() const { - return std::holds_alternative(this->data); + this->maxEntityIdx_ = ! this->maxEntityIdx_.has_value() + ? entity + : std::max(entity, *this->maxEntityIdx_); } -UDQUpdate RstUDQ::updateStatus() const +void Opm::RestartIO::RstUDQ::commitValues() { - return this->is_define() - ? std::get(this->data).status - : UDQUpdate::OFF; -} + if (! this->isUDQSet()) { + // Scalar or monostate. Nothing to do. + return; + } + + // If we get here this is a UDQ set. Form compressed UDQ mapping and + // reorder the values to match this compressed mapping. -void RstUDQ::add_value(const std::string& wgname, double value) { - if (this->is_define()) { - auto& def = std::get(this->data); - def.values.emplace_back(wgname, value); - } else { - auto& assign = std::get(this->data); - assign.update_value(this->name, value); - assign.selector.insert(wgname); + this->entityMap_.compress(this->maxEntityIdx_.value() + 1); + + auto newSA = std::vector(this->entityMap_.numEdges(), 0.0); + { + const auto& currSA = std::get>(this->sa_); + + const auto& ix = this->entityMap_.compressedIndexMap(); + + for (auto i = 0*currSA.size(); i < currSA.size(); ++i) { + newSA[ix[i]] += currSA[i]; + } } + + std::get>(this->sa_).swap(newSA); } -void RstUDQ::add_value(double value) { - if (this->is_define()) { - auto& def = std::get(this->data); - def.field_value = value; - } else { - auto& assign = std::get(this->data); - assign.update_value(this->name, value); +void Opm::RestartIO::RstUDQ::assignScalarValue(const double value) +{ + if (this->isUDQSet()) { + throw std::logic_error { + fmt::format("UDQ {} cannot be defined as a UDQ set " + "and then used as a scalar at restart time", + this->name) + }; } -} -double RstUDQ::assign_value() const { - const auto& assign = std::get(this->data); - return assign.value.value(); + this->sa_ = value; } -const std::unordered_set& RstUDQ::assign_selector() const { - const auto& assign = std::get(this->data); - return assign.selector; +void Opm::RestartIO::RstUDQ::addEntityName(std::string_view wgname) +{ + this->wgnames_.emplace_back(wgname); } -const std::string& RstUDQ::expression() const { - const auto& define = std::get(this->data); - return define.expression; +double Opm::RestartIO::RstUDQ::scalarValue() const +{ + if (std::holds_alternative(this->sa_)) { + throw std::logic_error { + fmt::format("Cannot request scalar value from UDQ " + "{} when no values have been assigned", + this->name) + }; + } + + if (! this->isScalar()) { + throw std::logic_error { + fmt::format("Cannot request scalar value " + "from non-scalar UDQ set {}", + this->name) + }; + } + + return std::get(this->sa_); } -const std::vector>& RstUDQ::values() const { - const auto& define = std::get(this->data); - return define.values; +std::size_t Opm::RestartIO::RstUDQ::numEntities() const +{ + return this->maxEntityIdx_.has_value() + ? static_cast(1 + *this->maxEntityIdx_) + : this->wgnames_.size(); } -std::optional RstUDQ::field_value() const { - const auto& define = std::get(this->data); - return define.field_value; +Opm::RestartIO::RstUDQ::ValueRange +Opm::RestartIO::RstUDQ::operator[](const std::size_t i) const +{ + if (std::holds_alternative(this->sa_)) { + throw std::logic_error { + fmt::format("Cannot request values for " + "entity {} from UDQ {} when " + "no values have been assigned", + i, this->name) + }; + } + + if (this->isScalar()) { + return this->scalarRange(i); + } + + return this->udqSetRange(i); } -RstUDQActive::RstRecord::RstRecord(UDAControl c, std::size_t i, std::size_t u1, std::size_t u2) - : control(c) - , input_index(i) - , use_count(u1) - , wg_offset(u2) -{} +Opm::UDQUpdate Opm::RestartIO::RstUDQ::currentUpdateStatus() const +{ + return this->isDefine() + ? this->definition_->status + : UDQUpdate::OFF; +} -RstUDQActive::RstUDQActive(const std::vector& iuad_arg, const std::vector& iuap, const std::vector& igph) +const std::vector& Opm::RestartIO::RstUDQ::nameIndex() const { - auto uda_size = UDQDims::entriesPerIUAD(); - for (std::size_t iuad_index = 0; iuad_index < iuad_arg.size() / uda_size; iuad_index++) { - auto offset = iuad_index * uda_size; - this->iuad.emplace_back( UDQ::udaControl(iuad_arg[offset + 0]), - iuad_arg[offset + 1] - 1, - iuad_arg[offset + 3], - iuad_arg[offset + 4] - 1); - } + this->ensureValidNameIndex(); - std::transform(iuap.begin(), iuap.end(), std::back_inserter(this->wg_index), [](const int& value) { return value - 1;}); + return *this->wgNameIdx_; +} - for (const auto& int_phase : igph) { - Phase phase{Phase::OIL}; +std::string_view Opm::RestartIO::RstUDQ::definingExpression() const +{ + return this->isDefine() + ? std::string_view { this->definition_->expression } + : std::string_view {}; +} - if (int_phase == 1) - phase = Phase::OIL; +Opm::RestartIO::RstUDQ::ValueRange +Opm::RestartIO::RstUDQ::scalarRange(const std::size_t i) const +{ + this->ensureValidNameIndex(); - if (int_phase == 2) - phase = Phase::WATER; + return ValueRange { + i, i + 1, this->wgNameIdx_->data(), + std::get(this->sa_) + }; +} - if (int_phase == 3) - phase = Phase::GAS; +Opm::RestartIO::RstUDQ::ValueRange +Opm::RestartIO::RstUDQ::udqSetRange(const std::size_t i) const +{ + const auto& start = this->entityMap_.startPointers(); + const auto& cols = this->entityMap_.columnIndices(); + const auto& sa = std::get>(this->sa_); - this->ig_phase.push_back(phase); + return ValueRange { + start[i + 0], start[i + 1], cols.data(), sa.data() }; } +void Opm::RestartIO::RstUDQ::ensureValidNameIndex() const +{ + if (! this->wgNameIdx_.has_value()) { + this->wgNameIdx_.emplace(this->wgnames_.size()); + + std::iota(this->wgNameIdx_->begin(), this->wgNameIdx_->end(), + std::vector::size_type{0}); + } } + +// --------------------------------------------------------------------------- + +Opm::RestartIO::RstUDQActive:: +RstRecord::RstRecord(const UDAControl c, + const std::size_t i, + const std::size_t u1, + const std::size_t u2) + : control (c) + , input_index(i) + , use_count (u1) + , wg_offset (u2) +{} + +Opm::RestartIO::RstUDQActive:: +RstUDQActive(const std::vector& iuad_arg, + const std::vector& iuap, + const std::vector& igph) + : wg_index { iuap } + , ig_phase (igph.size(), Phase::OIL) +{ + for (auto offset = 0*iuad_arg.size(); + offset < iuad_arg.size(); + offset += UDQDims::entriesPerIUAD()) + { + const auto* uda = &iuad_arg[offset]; + + this->iuad.emplace_back(UDQ::udaControl(uda[0]), + uda[1] - 1, + uda[3], + uda[4] - 1); + } + + std::transform(this->wg_index.begin(), + this->wg_index.end(), + this->wg_index.begin(), + [](const int wgObjectID) { return wgObjectID - 1; }); + + std::transform(igph.begin(), igph.end(), this->ig_phase.begin(), + [](const int phase) + { + if (phase == 1) { return Phase::OIL; } + if (phase == 2) { return Phase::WATER; } + if (phase == 3) { return Phase::GAS; } + + return Phase::OIL; + }); } diff --git a/opm/io/eclipse/rst/udq.hpp b/opm/io/eclipse/rst/udq.hpp index e443382241a..d4fda482d68 100644 --- a/opm/io/eclipse/rst/udq.hpp +++ b/opm/io/eclipse/rst/udq.hpp @@ -3,7 +3,8 @@ This file is part of the Open Porous Media project (OPM). - OPM is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. @@ -19,96 +20,508 @@ #ifndef RST_UDQ #define RST_UDQ +#include + +#include + +#include + #include +#include #include #include -#include +#include #include #include #include -#include - -#include -#include - -namespace Opm { - -namespace RestartIO { - -struct RstHeader; -class RstUDQ { +namespace Opm::RestartIO { + +/// Container for a single user defined quantity (UDQ) reconstituted from +/// restart file information. +/// +/// The producing side is expected to construct an RstUDQ object, and to +/// signal if the UDQ represents a quantity with a defining expression (four +/// parameter constructor) or without a defining expression (two parameter +/// constructor), the latter typically representing an assigned quantity +/// only. Moreover, the producing side should call prepareValues() prior to +/// incorporating numeric values. Then to assign the numeric values loaded +/// from a restart file using one of the functions assignScalarValue() or +/// addValue(). The former is intended for scalar UDQs, e.g., those at +/// field level, while the latter is intended for UDQ sets such as those +/// pertaining to wells, groups, or well segments. Mixing +/// assignScalarValue() and addValue() for a single UDQ object will make +/// class RstUDQ throw an exception. Once all values have been added, the +/// producing side is expected to call commitValues(). +/// +/// The consuming side, typically client code which forms \c UDQConfig and +/// \c UDQState objects, is expected to query the object for how to +/// interpret the values and then to call \code operator[]() \endcode to +/// retrieve a range of numeric values for the UDQ pertaining to a sing +/// entity. +class RstUDQ +{ public: - struct RstDefine { - RstDefine(const std::string& expression_arg, UDQUpdate status_arg); - - std::string expression; - UDQUpdate status; - std::vector> values; - std::optional field_value; - }; - - struct RstAssign { - void update_value(const std::string& name_arg, double new_value); - - std::optional value; - std::unordered_set selector; + /// Sequence of sub-entity IDs and values pertaining to a single entity. + class ValueRange + { + public: + /// Simple forward iterator over a UDQ set's values pertaining to a + /// single entity (e.g., a well or a group). + class Iterator + { + public: + /// Iterator's category (forward iterator) + using iterator_category = std::forward_iterator_tag; + + /// Iterator's value type + /// + /// Pair's \code .first \endcode element is the zero-based + /// sub-entity index, while \code .second \endcode is the + /// sub-entity's value. + using value_type = std::pair; + + /// Iterator's difference type. + using difference_type = int; + + /// Iterator's pointer type (return type from operator->()) + using pointer = const value_type*; + + /// Iterator's reference type (return type from operator*()) + using reference = const value_type&; + + /// Pre-increment operator. + /// + /// \return *this. + Iterator& operator++() + { + ++this->ix_; + + return *this; + } + + /// Post-increment operator. + /// + /// \return Iterator pointing to element prior to increment result. + Iterator operator++(int) + { + auto iter = *this; + + ++(*this); + + return iter; + } + + /// Dereference operator + /// + /// \return Element at current position. + reference operator*() + { + this->setValue(); + return this->deref_value_; + } + + /// Indirection operator + /// + /// \return Pointer to element at current position. + pointer operator->() + { + this->setValue(); + return &this->deref_value_; + } + + /// Equality predicate + /// + /// \param[in] that Object to which \c *this will be compared for equality. + /// + /// \return Whether or not \c *this equals \p that. + bool operator==(const Iterator& that) const + { + return (this->ix_ == that.ix_) + && (this->i_ == that.i_) + && (this->value_ == that.value_); + } + + /// Inequality predicate + /// + /// \param[in] that Object to which \c *this will be compared for inequality. + /// + /// \return \code ! (*this == that) + bool operator!=(const Iterator& that) const + { + return ! (*this == that); + } + + friend class ValueRange; + + private: + /// Constructor + /// + /// Accessible to RegionIndexRange only. + /// + /// \param[in] index range element value. + Iterator(std::size_t ix, + const int* i, + const std::variant& value) + : ix_{ix}, i_{i}, value_{value} + {} + + /// Index range element value + std::size_t ix_; + const int* i_; + std::variant value_; + + value_type deref_value_{}; + + void setValue(); + }; + + /// Start of value range + Iterator begin() const + { return this->makeIterator(this->begin_); } + + /// End of value range + Iterator end() const + { return this->makeIterator(this->end_); } + + friend class RstUDQ; + + private: + /// Constructor + /// + /// Scalar version. + /// + /// \param[in] begin_arg Index of start of range. + /// + /// \param[in] end_arg Index of element one past the end of the + /// range. + /// + /// \param[in] i Range of sub-entity indices. + /// + /// \param[in] value Constant value applicable to every element in + /// the range. + ValueRange(std::size_t begin_arg, + std::size_t end_arg, + const int* i, + const double value) + : begin_{begin_arg}, end_{end_arg}, i_{i}, value_{value} + {} + + /// Constructor + /// + /// Array version. + /// + /// \param[in] begin_arg Index of start of range. + /// + /// \param[in] end_arg Index of element one past the end of the + /// range. + /// + /// \param[in] i Range of sub-entity indices. + /// + /// \param[in] value Range of values, one scalar value for each + /// sub-entity in the value range. + ValueRange(std::size_t begin_arg, + std::size_t end_arg, + const int* i, + const double* value) + : begin_{begin_arg}, end_{end_arg}, i_{i}, value_{value} + {} + + /// Index of start of the range. + std::size_t begin_{}; + + /// Index of element one past the end of the range. + std::size_t end_{}; + + /// Sub-entities. + const int* i_{}; + + /// Values pertaining to each sub-entity of the range. + /// + /// Holds a \c double in the scalar version, and a \code const + /// double* \end in the array version. + std::variant value_{}; + + /// Form an iterator from a sequence index + /// + /// \param[in] ix Sequence index + /// + /// \return Range iterator pointing to position \p ix. + Iterator makeIterator(const std::size_t ix) const + { + return { ix, this->i_, this->value_ }; + } }; - + /// UDQ name. + /// + /// First argument to a DEFINE or an ASSIGN statement. + std::string name{}; + + /// UDQ's unit string. + std::string unit{}; + + /// UDQ's category, i.e., which level this UDQ applies to. + /// + /// Exmples include, the FIELD (FU*), group (GU*), well (WU*), or well + /// segments (SU*). + UDQVarType category{UDQVarType::NONE}; + + /// Constructor + /// + /// Quantity given by DEFINE statement. + /// + /// \param{in] name_arg UDQ name. + /// + /// \param{in] unit_arg UDQ's unit string. + /// + /// \param{in] define_arg UDQ's defining expression, previously entered + /// in a DEFINE statement. + /// + /// \param{in] status_arg UDQ's update status, i.e., when to calculate + /// new values for the quantity. RstUDQ(const std::string& name_arg, const std::string& unit_arg, const std::string& define_arg, - UDQUpdate status_arg); - + UDQUpdate status_arg); + + /// Constructor + /// + /// Quantity given by ASSIGN statement. + /// + /// \param{in] name_arg UDQ name. + /// + /// \param{in] unit_arg UDQ's unit string. RstUDQ(const std::string& name_arg, const std::string& unit_arg); - void add_value(double value); - void add_value(const std::string& wgname, double value); - - bool is_define() const; - UDQUpdate updateStatus() const; - double assign_value() const; - const std::unordered_set& assign_selector() const; - const std::string& expression() const; - const std::vector>& values() const; - std::optional field_value() const; - - std::string name; - std::string unit; - UDQVarType var_type; + /// Prepare inclusion of entities and values. + /// + /// Used by producing code, i.e., the layer which reads UDQ information + /// directly from file. + void prepareValues(); + + /// Assign numeric UDQ value for a entity/sub-entity pair. + /// + /// Conflicts with assignScalarValue(). If assignScalarValue() has been + /// called previously, then this function will throw an exception of + /// type \code std::logic_error \endcode. + /// + /// \param[in] entity Numeric entity index. Assumed to be zero-based + /// linear index of an entity, such as a well or a group. + /// + /// \param[in] subEntity Numeric sub-entity index. Assumed to be a + /// zero-based linear index of a sub-entity such as a well segment. For + /// well or group level UDQs, pass zero for the sub-entity. + /// + /// \param[in] value Numeric UDQ value for this entity/sub-entity pair. + void addValue(const int entity, const int subEntity, const double value); + + /// End value accumulation. + /// + /// Builds entity to sub-entity mapping, and sorts the corresponding + /// numeric UDQ values. Should normally be the final member function + /// call on the producing side. + void commitValues(); + + /// Assign a scalar value for the UDQ. + /// + /// Conflicts with addValue(). If addValue() has been called + /// previously, then this function will throw an exception of type \code + /// std::logic_error \endcode. + /// + /// \param[in] value Single numeric value pertaining to all entities in + /// this UDQ set. + void assignScalarValue(const double value); + + /// Add a name for the last (or next) entity used in a call to addValue() + /// + /// \param[in] wgname Entity, typically well or group, name. + void addEntityName(std::string_view wgname); + + /// Retrieve UDQ's scalar value. + /// + /// Throws a \code std::logic_error \endcode unless isScalar() returns + /// \c true. + /// + /// Used by consuming code, i.e., the layer which forms \c UDQConfig and + /// \c UDQState objects from restart file information. + double scalarValue() const; + + /// Retrieve number of of entities known to this UDQ. + std::size_t numEntities() const; + + /// Get read-only access to sub-entities and values associated to a + /// single top-level entity. + /// + /// \param[in] i Entity index. Must be in the range [0..numEntities()). + /// + /// \return Sequence of values pertaining to entity \p i. + ValueRange operator[](const std::size_t i) const; + + /// Get sequence of UDQ's entity names. + /// + /// \return Sequence of entity names entered in earlier addEntityName() + /// calls. Order not guaranteed. + const std::vector& entityNames() const + { + return this->wgnames_; + } + + /// Index map for entity names. + /// + /// The entity name of entity 'i' is + /// \code + /// entityNames()[ nameIndex()[i] ] + /// \endcode + /// provided named entities are meaningful for this UDQ--i.e., if it + /// pertains to the well, group, connection, or segment levels. + const std::vector& nameIndex() const; + + /// UDQ's defining expression + /// + /// Empty character sequence if this UDQ does not have a defining + /// expression--i.e., if it was entered using ASSIGN statements. + std::string_view definingExpression() const; + + /// UDQ's current update status. + UDQUpdate currentUpdateStatus() const; + + /// Whether or not this UDQ has a defining expression + /// + /// Needed to properly reconstitute the UDQConfig object from restart + /// file information. + bool isDefine() const + { + return this->definition_.has_value(); + } + + /// Predicate for whether or not this UDQ is a scalar + /// quantity--typically a FIELD level value. + bool isScalar() const + { + return std::holds_alternative(this->sa_); + } private: - std::variant data; -}; - - + /// Entity mapping type. + /// + /// VertexID = int + /// TrackCompressedIdx = true (need SA mapping) + /// PermitSelfConnections = true (MS well 5 may have segment number 5). + using Graph = utility::CSRGraphFromCoordinates; + + /// Wrapper for a DEFINE expression + struct Definition + { + /// Constructor + /// + /// Mostly to enable using optional<>::emplace(). + /// + /// \param[in] expression_arg UDQ's defining expression + /// + /// \param[in] status_arg UDQ's update status. + Definition(const std::string& expression_arg, + const UDQUpdate status_arg) + : expression { expression_arg } + , status { status_arg } + {} + + /// UDQ's defining expression. + std::string expression{}; + + /// UDQ's update status. + UDQUpdate status{}; + }; -class RstUDQActive { + /// Map entities to range of sub-entities. + /// + /// Examples include + /// + /// -# Well and its segments (SU*) + /// -# Well and its connections (CU*) + /// + /// For field, group, and well level UDQs, the sub-entities are usually + /// trivial (i.e., one sub-entity of numeric index zero), although one + /// *might* choose to represent aquifer and block level UDQs as + /// sub-entities of the FIELD. + Graph entityMap_{}; + + /// UDQ values. + /// + /// Single scalar if all entities in this UDQ share the same constant + /// value--typically only for field level UDQs. Range of scalars (\code + /// vector \endcode) if each entity/sub-entity potentially has a + /// different value. Monostate unless any values have been assigned. + std::variant> sa_; + + /// Entity names. + /// + /// Typically well or group names. + std::vector wgnames_{}; + + /// Largest entity index seen in all addValue() calls so far. + /// + /// Nullopt before first call to addValue(). + std::optional maxEntityIdx_{}; + + /// Map entity indices to entity names. + mutable std::optional> wgNameIdx_{}; + + /// UDQ's definition. + /// + /// Nullopt unless this UDQ has a defining expression (four parameter + /// constructor used to form the object). + std::optional definition_{}; + + /// Whether or not this UDQ represents a non-scalar UDQ set. + bool isUDQSet() const + { + return std::holds_alternative>(this->sa_); + } + + /// Construct a value range for a set of sub-entities that share a + /// common, constant numeric value (\c sa_ is \c double). + /// + /// \param[in] i Entity index. + /// + /// \return Value range for entity \p i. + ValueRange scalarRange(const std::size_t i) const; + + /// Construct a value range for a set of entities whose numeric values + /// may differ between sub-entities (\c sa_ is \code vector + /// \endcode). + /// + /// \param[in] i Entity index. + /// + /// \return Value range for entity \p i. + ValueRange udqSetRange(const std::size_t i) const; + + void ensureValidNameIndex() const; +}; - struct RstRecord { +struct RstUDQActive +{ + struct RstRecord + { RstRecord(UDAControl c, std::size_t i, std::size_t u1, std::size_t u2); - - UDAControl control; + UDAControl control; std::size_t input_index; std::size_t use_count; std::size_t wg_offset; }; - -public: RstUDQActive() = default; - RstUDQActive(const std::vector& iuad, const std::vector& iuap, const std::vector& igph); + RstUDQActive(const std::vector& iuad, + const std::vector& iuap, + const std::vector& igph); - std::vector iuad; - std::vector wg_index; - std::vector ig_phase; + std::vector iuad{}; + std::vector wg_index{}; + std::vector ig_phase{}; }; -} -} - +} // namespace Opm::RestartIO -#endif +#endif // RST_UDQ diff --git a/tests/parser/UDQTests.cpp b/tests/parser/UDQTests.cpp index c3572fd455f..9a085591e34 100644 --- a/tests/parser/UDQTests.cpp +++ b/tests/parser/UDQTests.cpp @@ -24,6 +24,8 @@ #include #include +#include + #include #include #include @@ -50,14 +52,14 @@ #include +#include +#include + #include #include #include #include -#include -#include - #include #include #include @@ -2857,14 +2859,32 @@ UDQ BOOST_CHECK_NO_THROW( make_schedule(valid) ); } -BOOST_AUTO_TEST_CASE(UDQ_ASSIGN_RST) { - std::unordered_set selector{"W1", "W2"}; - UDQAssign assign("WUBHP", selector, 100, 2); - auto res = assign.eval( {"W1", "W2", "W3"}); - BOOST_CHECK_EQUAL(res.size(), 3); - BOOST_CHECK_EQUAL(res["W1"].get(), 100); - BOOST_CHECK_EQUAL(res["W2"].get(), 100); - BOOST_CHECK_EQUAL(res["W3"].defined(), false); +BOOST_AUTO_TEST_CASE(UDQ_ASSIGN_RST) +{ + using namespace std::string_literals; + + auto assignRst = RestartIO::RstUDQ { "WUBHP", "BARSA" }; + + assignRst.prepareValues(); + assignRst.addValue(0, 0, 100.0); + assignRst.addEntityName("W1"); + + assignRst.addValue(1, 0, 123.4); + assignRst.addEntityName("W2"); + assignRst.commitValues(); + + const auto report_step = std::size_t{2}; + + const auto assign = UDQAssign { "WUBHP", assignRst, report_step }; + const auto res = assign.eval(std::vector { "W1"s, "W2"s, "W3"s, }); + + BOOST_REQUIRE_EQUAL(res.size(), 3); + + BOOST_CHECK_CLOSE(res["W1"].get(), 100.0, 1.0e-8); + BOOST_CHECK_CLOSE(res["W2"].get(), 123.4, 1.0e-8); + + BOOST_CHECK_MESSAGE(! res["W3"].defined(), + R"(Assignment UDQ set must NOT have a defined value for well "W3")"); } BOOST_AUTO_TEST_CASE(UDQ_ASSIGN_SEGMENT)