From 33a27e4e88ccfd3c2c78af5d91f21c8198679e51 Mon Sep 17 00:00:00 2001 From: Vegard Kippe Date: Fri, 15 Sep 2023 15:18:52 +0200 Subject: [PATCH] Facilitate UDA for WSEGVALV item 4 --- opm/input/eclipse/Schedule/MSW/Segment.hpp | 1 + opm/input/eclipse/Schedule/MSW/Valve.hpp | 23 +++++++++- .../eclipse/Schedule/MSW/WellSegments.hpp | 1 + opm/input/eclipse/Schedule/Well/Well.hpp | 1 + .../input/eclipse/Schedule/MSW/Segment.cpp | 4 ++ src/opm/input/eclipse/Schedule/MSW/Valve.cpp | 44 ++++++++++++++++--- .../eclipse/Schedule/MSW/WellSegments.cpp | 4 ++ src/opm/input/eclipse/Schedule/Well/Well.cpp | 8 ++++ .../share/keywords/000_Eclipse100/W/WSEGVALV | 2 +- src/opm/output/eclipse/AggregateMSWData.cpp | 2 +- src/opm/output/eclipse/WriteRFT.cpp | 2 +- 11 files changed, 81 insertions(+), 11 deletions(-) diff --git a/opm/input/eclipse/Schedule/MSW/Segment.hpp b/opm/input/eclipse/Schedule/MSW/Segment.hpp index 0177582ad..1e4c09b24 100644 --- a/opm/input/eclipse/Schedule/MSW/Segment.hpp +++ b/opm/input/eclipse/Schedule/MSW/Segment.hpp @@ -123,6 +123,7 @@ namespace Opm { const SICD& spiralICD() const; const AutoICD& autoICD() const; const Valve& valve() const; + Valve& valve(); void updatePerfLength(double perf_length); void updateSpiralICD(const SICD& spiral_icd); diff --git a/opm/input/eclipse/Schedule/MSW/Valve.hpp b/opm/input/eclipse/Schedule/MSW/Valve.hpp index 8a884ddbc..4ab0ffd44 100644 --- a/opm/input/eclipse/Schedule/MSW/Valve.hpp +++ b/opm/input/eclipse/Schedule/MSW/Valve.hpp @@ -24,8 +24,11 @@ #include #include #include +#include #include +#include +#include namespace Opm { @@ -34,6 +37,19 @@ namespace Opm { class DeckKeyword; class Segment; + struct ValveUDAEval { + const SummaryState& summary_state; + const std::string& well_name; + const size_t segment_number; + const double udq_default; + + ValveUDAEval(const SummaryState& summary_state_, const std::string& well_name_, + const size_t segment_number_, const double udq_default_); + + double value(const UDAValue& value) const; + }; + + class Valve { public: @@ -58,7 +74,8 @@ namespace Opm { // parameters for constriction pressure loss double conFlowCoefficient() const; - double conCrossArea() const; + double conCrossArea(const std::optional& uda_eval = std::nullopt); + inline double conCrossAreaValue() const { return m_con_cross_area_value; } double conMaxCrossArea() const; double pipeDiameter() const; double pipeRoughness() const; @@ -84,6 +101,7 @@ namespace Opm { { serializer(m_con_flow_coeff); serializer(m_con_cross_area); + serializer(m_con_cross_area_value); serializer(m_con_max_cross_area); serializer(m_pipe_additional_length); serializer(m_pipe_diameter); @@ -94,7 +112,8 @@ namespace Opm { private: double m_con_flow_coeff; - double m_con_cross_area; + UDAValue m_con_cross_area; + double m_con_cross_area_value; double m_con_max_cross_area; double m_pipe_additional_length; diff --git a/opm/input/eclipse/Schedule/MSW/WellSegments.hpp b/opm/input/eclipse/Schedule/MSW/WellSegments.hpp index e97a30236..e4a4b1b04 100644 --- a/opm/input/eclipse/Schedule/MSW/WellSegments.hpp +++ b/opm/input/eclipse/Schedule/MSW/WellSegments.hpp @@ -93,6 +93,7 @@ namespace Opm { const Segment& getFromSegmentNumber(const int segment_number) const; const Segment& operator[](size_t idx) const; + Segment& operator[](size_t idx); void orderSegments(); void updatePerfLength(const WellConnections& connections); diff --git a/opm/input/eclipse/Schedule/Well/Well.hpp b/opm/input/eclipse/Schedule/Well/Well.hpp index 76bde9442..d65de318c 100644 --- a/opm/input/eclipse/Schedule/Well/Well.hpp +++ b/opm/input/eclipse/Schedule/Well/Well.hpp @@ -397,6 +397,7 @@ public: const std::vector getConnections(int completion) const; const WellConnections& getConnections() const; const WellSegments& getSegments() const; + WellSegments& getSegments(); int maxSegmentID() const; int maxBranchID() const; diff --git a/src/opm/input/eclipse/Schedule/MSW/Segment.cpp b/src/opm/input/eclipse/Schedule/MSW/Segment.cpp index 59a998eb9..751a9fdef 100644 --- a/src/opm/input/eclipse/Schedule/MSW/Segment.cpp +++ b/src/opm/input/eclipse/Schedule/MSW/Segment.cpp @@ -374,6 +374,10 @@ namespace Opm { return std::get(this->m_icd); } + Valve& Segment::valve() { + return std::get(this->m_icd); + } + int Segment::ecl_type_id() const { switch (this->segmentType()) { case SegmentType::REGULAR: diff --git a/src/opm/input/eclipse/Schedule/MSW/Valve.cpp b/src/opm/input/eclipse/Schedule/MSW/Valve.cpp index fae47e64a..995ce16f1 100644 --- a/src/opm/input/eclipse/Schedule/MSW/Valve.cpp +++ b/src/opm/input/eclipse/Schedule/MSW/Valve.cpp @@ -24,6 +24,31 @@ namespace Opm { + ValveUDAEval::ValveUDAEval(const SummaryState& summary_state_, const std::string& well_name_, + const size_t segment_number_, const double udq_default_) : + summary_state(summary_state_), + well_name(well_name_), + segment_number(segment_number_), + udq_default(udq_default_) + {} + + + double ValveUDAEval::value(const UDAValue& value) const { + if (value.is()) + return value.getSI(); + + const std::string& string_var = value.get(); + double output_value = udq_default; + + if (summary_state.has_segment_var(well_name, value.get(), segment_number)) + output_value = summary_state.get_segment_var(well_name, string_var, segment_number); + else if (summary_state.has(string_var)) + output_value = summary_state.get(string_var); + + return value.get_dim().convertRawToSi(output_value); + } + + Valve::Valve() : Valve(0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, ICDStatus::SHUT) { @@ -38,7 +63,8 @@ namespace Opm { double pipeCrossA, ICDStatus stat) : m_con_flow_coeff(conFlowCoeff) - , m_con_cross_area(conCrossA) + , m_con_cross_area(UDAValue(conCrossA)) + , m_con_cross_area_value(conCrossA) , m_con_max_cross_area(conMaxCrossA) , m_pipe_additional_length(pipeAddLength) , m_pipe_diameter(pipeDiam) @@ -50,7 +76,8 @@ namespace Opm { Valve::Valve(const DeckRecord& record) : m_con_flow_coeff(record.getItem("CV").get(0)) - , m_con_cross_area(record.getItem("AREA").getSIDouble(0)) + , m_con_cross_area(record.getItem("AREA").get(0)) + , m_con_cross_area_value(m_con_cross_area.is() ? m_con_cross_area.getSI() : -1.0) { // we initialize negative values for the values are defaulted const double value_for_default = -1.e100; @@ -98,7 +125,8 @@ namespace Opm { { Valve result; result.m_con_flow_coeff = 1.0; - result.m_con_cross_area = 2.0; + result.m_con_cross_area = UDAValue(2.0); + result.m_con_cross_area_value = 2.0; result.m_con_max_cross_area = 3.0; result.m_pipe_additional_length = 4.0; result.m_pipe_diameter = 5.0; @@ -134,8 +162,11 @@ namespace Opm { return m_con_flow_coeff; } - double Valve::conCrossArea() const { - return m_con_cross_area; + double Valve::conCrossArea(const std::optional& uda_eval_optional) { + m_con_cross_area_value = uda_eval_optional.has_value() ? + uda_eval_optional.value().value(m_con_cross_area) : + m_con_cross_area.getSI(); + return m_con_cross_area_value; } double Valve::pipeAdditionalLength() const { @@ -180,7 +211,8 @@ namespace Opm { bool Valve::operator==(const Valve& data) const { return this->conFlowCoefficient() == data.conFlowCoefficient() && - this->conCrossArea() == data.conCrossArea() && + this->m_con_cross_area == data.m_con_cross_area && + this->conCrossAreaValue() == data.conCrossAreaValue() && this->conMaxCrossArea() == data.conMaxCrossArea() && this->pipeAdditionalLength() == data.pipeAdditionalLength() && this->pipeDiameter() == data.pipeDiameter() && diff --git a/src/opm/input/eclipse/Schedule/MSW/WellSegments.cpp b/src/opm/input/eclipse/Schedule/MSW/WellSegments.cpp index 863671e15..aa1eae1e2 100644 --- a/src/opm/input/eclipse/Schedule/MSW/WellSegments.cpp +++ b/src/opm/input/eclipse/Schedule/MSW/WellSegments.cpp @@ -137,6 +137,10 @@ namespace Opm { return m_segments[idx]; } + Segment& WellSegments::operator[](size_t idx) { + return m_segments[idx]; + } + int WellSegments::segmentNumberToIndex(const int segment_number) const { const auto it = segment_number_to_index.find(segment_number); if (it != segment_number_to_index.end()) { diff --git a/src/opm/input/eclipse/Schedule/Well/Well.cpp b/src/opm/input/eclipse/Schedule/Well/Well.cpp index 6ce34f26d..6f96783cf 100644 --- a/src/opm/input/eclipse/Schedule/Well/Well.cpp +++ b/src/opm/input/eclipse/Schedule/Well/Well.cpp @@ -1144,6 +1144,14 @@ const WellSegments& Well::getSegments() const { throw std::logic_error("Asked for segment information in not MSW well: " + this->name()); } +WellSegments& Well::getSegments() { + if (this->segments) + return *this->segments; + else + throw std::logic_error("Asked for segment information in not MSW well: " + this->name()); +} + + int Well::maxSegmentID() const { return (this->segments == nullptr) diff --git a/src/opm/input/eclipse/share/keywords/000_Eclipse100/W/WSEGVALV b/src/opm/input/eclipse/share/keywords/000_Eclipse100/W/WSEGVALV index d0a2a5b3a..aa2df61a1 100644 --- a/src/opm/input/eclipse/share/keywords/000_Eclipse100/W/WSEGVALV +++ b/src/opm/input/eclipse/share/keywords/000_Eclipse100/W/WSEGVALV @@ -22,7 +22,7 @@ { "item": 4, "name": "AREA", - "value_type": "DOUBLE", + "value_type": "UDA", "dimension": "Length*Length" }, { diff --git a/src/opm/output/eclipse/AggregateMSWData.cpp b/src/opm/output/eclipse/AggregateMSWData.cpp index 4e7be694a..1b36ef8e4 100644 --- a/src/opm/output/eclipse/AggregateMSWData.cpp +++ b/src/opm/output/eclipse/AggregateMSWData.cpp @@ -577,7 +577,7 @@ namespace { usys.from_si(M::length, valve.pipeAdditionalLength()); rSeg[baseIndex + Ix::ValveArea] = - usys.from_si(M::length, usys.from_si(M::length, valve.conCrossArea())); + usys.from_si(M::length, usys.from_si(M::length, valve.conCrossAreaValue())); rSeg[baseIndex + Ix::ValveFlowCoeff] = valve.conFlowCoefficient(); rSeg[baseIndex + Ix::ValveMaxArea] = diff --git a/src/opm/output/eclipse/WriteRFT.cpp b/src/opm/output/eclipse/WriteRFT.cpp index db7cf4488..fab26119a 100644 --- a/src/opm/output/eclipse/WriteRFT.cpp +++ b/src/opm/output/eclipse/WriteRFT.cpp @@ -1445,7 +1445,7 @@ namespace { { assert ((segment.segmentNumber() > 1) && segment.isValve()); - const auto Ac = segment.valve().conCrossArea(); + const auto Ac = segment.valve().conCrossAreaValue(); const auto Ac_max = valveMaximumCrossSectionalArea(segment); return Ac / Ac_max;