From 0301d402a4a8e7ee4955e55df05c9f90583110b9 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Thu, 31 Jan 2019 11:25:57 +0100 Subject: [PATCH 01/20] Load Restart: Chase Semantic Update of data::Connection Commit ef606711 switched to saving data::Connection::index as the global cell index of the reservoir connection. Update restore layer accordingly. --- src/opm/output/eclipse/LoadRestart.cpp | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/src/opm/output/eclipse/LoadRestart.cpp b/src/opm/output/eclipse/LoadRestart.cpp index 34b10a07a..844e9e738 100644 --- a/src/opm/output/eclipse/LoadRestart.cpp +++ b/src/opm/output/eclipse/LoadRestart.cpp @@ -332,12 +332,10 @@ namespace { continue; } - const auto active_index = grid.activeIndex(i, j, k); - well.connections.emplace_back(); auto& connection = well.connections.back(); - connection.index = active_index; + connection.index = grid.getGlobalIndex(i, j, k); connection.pressure = *opm_xwel_data++; connection.reservoir_rate = *opm_xwel_data++; From e77bd172c5259c65aa2139e571cc44ed4e432c07 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Fri, 21 Dec 2018 13:59:21 +0100 Subject: [PATCH 02/20] Solution Output: Remove Special Treatment for TEMP It is no longer needed as this aspect is handled further up in the call chain. --- src/opm/output/eclipse/RestartIO.cpp | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/src/opm/output/eclipse/RestartIO.cpp b/src/opm/output/eclipse/RestartIO.cpp index d00b98411..72b9daaa8 100644 --- a/src/opm/output/eclipse/RestartIO.cpp +++ b/src/opm/output/eclipse/RestartIO.cpp @@ -381,8 +381,8 @@ namespace { void writeSolution(ecl_rst_file_type* rst_file, const RestartValue& value, - const bool ecl_compatible_rst, - const bool write_double_arg) + const bool ecl_compatible_rst, + const bool write_double_arg) { ecl_rst_file_start_solution(rst_file); @@ -396,8 +396,6 @@ namespace { }; for (const auto& elm : value.solution) { - if (ecl_compatible_rst && (elm.first == "TEMP")) continue; - if (elm.second.target == data::TargetType::RESTART_SOLUTION) { write(elm.first, elm.second.data, write_double_arg); @@ -415,7 +413,9 @@ namespace { ecl_rst_file_end_solution(rst_file); - if (ecl_compatible_rst) return; + if (ecl_compatible_rst) { + return; + } for (const auto& elm : value.solution) { if (elm.second.target == data::TargetType::RESTART_AUXILIARY) { From ff1b245b8f2363feb38057b8853d033605f7976a Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Fri, 2 Nov 2018 13:59:17 +0100 Subject: [PATCH 03/20] Restart Facility: Idenfity Certain Items in XGRP and RSEG This commit adds a few identifiers for items in the XGRP and RSEG restart vectors to enable restoring segment rates and cumulative group production/injection quantities. --- CMakeLists_files.cmake | 2 + opm/output/eclipse/VectorItems/group.hpp | 56 +++++++++++++++++++++ opm/output/eclipse/VectorItems/msw.hpp | 64 ++++++++++++++++++++++++ 3 files changed, 122 insertions(+) create mode 100644 opm/output/eclipse/VectorItems/group.hpp create mode 100644 opm/output/eclipse/VectorItems/msw.hpp diff --git a/CMakeLists_files.cmake b/CMakeLists_files.cmake index 6a9ccaa23..1714be2a5 100644 --- a/CMakeLists_files.cmake +++ b/CMakeLists_files.cmake @@ -528,7 +528,9 @@ if(ENABLE_ECL_OUTPUT) opm/output/data/Solution.hpp opm/output/data/Wells.hpp opm/output/eclipse/VectorItems/connection.hpp + opm/output/eclipse/VectorItems/group.hpp opm/output/eclipse/VectorItems/intehead.hpp + opm/output/eclipse/VectorItems/msw.hpp opm/output/eclipse/VectorItems/well.hpp opm/output/eclipse/AggregateGroupData.hpp opm/output/eclipse/AggregateConnectionData.hpp diff --git a/opm/output/eclipse/VectorItems/group.hpp b/opm/output/eclipse/VectorItems/group.hpp new file mode 100644 index 000000000..8c714fc06 --- /dev/null +++ b/opm/output/eclipse/VectorItems/group.hpp @@ -0,0 +1,56 @@ +/* + Copyright (c) 2018 Equinor ASA + + 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 + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + +#ifndef OPM_OUTPUT_ECLIPSE_VECTOR_GROUP_HPP +#define OPM_OUTPUT_ECLIPSE_VECTOR_GROUP_HPP + +#include + +namespace Opm { namespace RestartIO { namespace Helpers { namespace VectorItems { + + namespace XGroup { + enum index : std::vector::size_type { + OilPrRate = 0, // Group's oil production rate + WatPrRate = 1, // Group's water production rate + GasPrRate = 2, // Group's gas production rate + LiqPrRate = 3, // Group's liquid production rate + + WatInjRate = 5, // Group's water injection rate + GasInjRate = 6, // Group's gas injection rate + + WatCut = 8, // Group's producing water cut + GORatio = 9, // Group's producing gas/oil ratio + + OilPrTotal = 10, // Group's total cumulative oil production + WatPrTotal = 11, // Group's total cumulative water production + GasPrTotal = 12, // Group's total cumulative gas production + VoidPrTotal = 13, // Group's total cumulative reservoir + // voidage production + + WatInjTotal = 15, // Group's total cumulative water injection + GasInjTotal = 16, // Group's total cumulative gas injection + + OilPrPot = 22, // Group's oil production potential + WatPrPot = 23, // Group's water production potential + }; + } // XGroup + +}}}} // Opm::RestartIO::Helpers::VectorItems + +#endif // OPM_OUTPUT_ECLIPSE_VECTOR_GROUP_HPP diff --git a/opm/output/eclipse/VectorItems/msw.hpp b/opm/output/eclipse/VectorItems/msw.hpp new file mode 100644 index 000000000..e0a0bea72 --- /dev/null +++ b/opm/output/eclipse/VectorItems/msw.hpp @@ -0,0 +1,64 @@ +/* + Copyright (c) 2018 Equinor ASA + + 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 + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + +#ifndef OPM_OUTPUT_ECLIPSE_VECTOR_MSW_HPP +#define OPM_OUTPUT_ECLIPSE_VECTOR_MSW_HPP + +#include + +namespace Opm { namespace RestartIO { namespace Helpers { namespace VectorItems { + + namespace ISeg { + enum index : std::vector::size_type { + SegNo = 0, // Segment number (one-based) + OutSeg = 1, // Outlet segment (one-based) + InSegCurBranch = 2, // Inflow segment current branch (one-based) + BranchNo = 3, // Branch number (one-based) + }; + } // ISeg + + namespace RSeg { + enum index : std::vector::size_type { + DistOutlet = 0, // Segment's distance to outlet + OutletDepthDiff = 1, // Segment's depth differential to outlet + SegDiam = 2, // Internal diameter of segment + SegRough = 3, // Roughness parameter of segment + SegArea = 4, // Cross-sectional area of segment + SegVolume = 5, // Physical volume of segment + DistBHPRef = 6, // Segment's distance to BHP reference node + DepthBHPRef = 7, // Segment's depth differential to BHP ref. node + + TotFlowRate = 8, // Normalised total segment flow rate + WatFlowFract = 9, // Normalised Water flow rate fraction + GasFlowFract = 10, // Normalised Gas flow rate fraction + Pressure = 11, // Segment pressure + + item40 = 39, // Unknown + item106 = 105, // Unknown + item107 = 106, // Unknown + item108 = 107, // Unknown + item109 = 108, // Unknown + item110 = 109, // Unknown + item111 = 110, // Unknown + }; + } // RSeg + +}}}} // Opm::RestartIO::Helpers::VectorItems + +#endif // OPM_OUTPUT_ECLIPSE_VECTOR_MSW_HPP From 5736dd6e9cb1745dd259028ee47cfa4cf5070435 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Fri, 2 Nov 2018 14:14:21 +0100 Subject: [PATCH 04/20] Load Restart: Add Windowed Views into Well and Group Data This commit adds two new helper classes, {Well,Group}Vectors which simplify accessing the portions of {I,X}{CON,WEL,GRP} pertaining to any individual connection, well or group. These classes are both implemented in terms of the existing helper functions getPtr() getDataWindow() getInteHeadElem() but hide away the complexity of passing the correct maximum number of connections and number of vector elements per connection, well or group. Reimplement the vector queries in terms of these helper classes to (hopefully) reduce the mental load when reading the restart code. --- src/opm/output/eclipse/LoadRestart.cpp | 354 +++++++++++++++++-------- 1 file changed, 239 insertions(+), 115 deletions(-) diff --git a/src/opm/output/eclipse/LoadRestart.cpp b/src/opm/output/eclipse/LoadRestart.cpp index 844e9e738..30032d06d 100644 --- a/src/opm/output/eclipse/LoadRestart.cpp +++ b/src/opm/output/eclipse/LoadRestart.cpp @@ -1,4 +1,5 @@ /* + Copyright (c) 2018 Equinor ASA Copyright (c) 2016 Statoil ASA Copyright (c) 2013-2015 Andreas Lauser Copyright (c) 2013 SINTEF ICT, Applied Mathematics. @@ -145,6 +146,227 @@ RestartFileView& RestartFileView::operator=(RestartFileView&& rhs) return *this; } +// --------------------------------------------------------------------- + +namespace { + template + const T* getPtr(const ::Opm::RestartIO::ecl_kw_type* kw) + { + return (kw == nullptr) ? nullptr + : static_cast(ecl_kw_iget_ptr(kw /* <- ADL */, 0)); + } + + template + boost::iterator_range + getDataWindow(const T* 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); + + const auto* begin = arr + off; + const auto* end = begin + windowSize; + + return { begin, end }; + } + + int getInteHeadElem(const ::Opm::RestartIO::ecl_kw_type* intehead, + const std::vector::size_type i) + { + return getPtr(intehead)[i]; + } +} + +// --------------------------------------------------------------------- + +class WellVectors +{ +public: + template + using Window = boost::iterator_range; + + explicit WellVectors(const RestartFileView& rst_view, + const ::Opm::RestartIO::ecl_kw_type* intehead); + + bool hasDefinedWellValues() const; + bool hasDefinedConnectionValues() const; + + Window iwel(const std::size_t wellID) const; + Window xwel(const std::size_t wellID) const; + + Window + icon(const std::size_t wellID, const std::size_t connID) const; + + Window + xcon(const std::size_t wellID, const std::size_t connID) const; + +private: + std::size_t maxConnPerWell_; + std::size_t numIWelElem_; + std::size_t numXWelElem_; + std::size_t numIConElem_; + std::size_t numXConElem_; + + const ::Opm::RestartIO::ecl_kw_type* iwel_; + const ::Opm::RestartIO::ecl_kw_type* xwel_; + + const ::Opm::RestartIO::ecl_kw_type* icon_; + const ::Opm::RestartIO::ecl_kw_type* xcon_; +}; + +WellVectors::WellVectors(const RestartFileView& rst_view, + const ::Opm::RestartIO::ecl_kw_type* intehead) + : maxConnPerWell_(getInteHeadElem(intehead, VI::intehead::NCWMAX)) + , numIWelElem_ (getInteHeadElem(intehead, VI::intehead::NIWELZ)) + , numXWelElem_ (getInteHeadElem(intehead, VI::intehead::NXWELZ)) + , numIConElem_ (getInteHeadElem(intehead, VI::intehead::NICONZ)) + , numXConElem_ (getInteHeadElem(intehead, VI::intehead::NXCONZ)) + , iwel_ (rst_view.getKeyword("IWEL")) + , xwel_ (rst_view.getKeyword("XWEL")) + , icon_ (rst_view.getKeyword("ICON")) + , xcon_ (rst_view.getKeyword("XCON")) +{} + +bool WellVectors::hasDefinedWellValues() const +{ + return ! ((this->iwel_ == nullptr) || + (this->xwel_ == nullptr)); +} + +bool WellVectors::hasDefinedConnectionValues() const +{ + return ! ((this->icon_ == nullptr) || + (this->xcon_ == nullptr)); +} + +WellVectors::Window +WellVectors::iwel(const std::size_t wellID) const +{ + if (! this->hasDefinedWellValues()) { + throw std::logic_error { + "Cannot Request IWEL Values Unless Defined" + }; + } + + return getDataWindow(getPtr(this->iwel_), + this->numIWelElem_, wellID); +} + +WellVectors::Window +WellVectors::xwel(const std::size_t wellID) const +{ + if (! this->hasDefinedWellValues()) { + throw std::logic_error { + "Cannot Request XWEL Values Unless Defined" + }; + } + + return getDataWindow(getPtr(this->xwel_), + this->numXWelElem_, wellID); +} + +WellVectors::Window +WellVectors::icon(const std::size_t wellID, const std::size_t connID) const +{ + if (! this->hasDefinedConnectionValues()) { + throw std::logic_error { + "Cannot Request ICON Values Unless Defined" + }; + } + + return getDataWindow(getPtr(this->icon_), this->numIConElem_, + wellID, connID, this->maxConnPerWell_); +} + +WellVectors::Window +WellVectors::xcon(const std::size_t wellID, const std::size_t connID) const +{ + if (! this->hasDefinedConnectionValues()) { + throw std::logic_error { + "Cannot Request XCON Values Unless Defined" + }; + } + + return getDataWindow(getPtr(this->xcon_), this->numXConElem_, + wellID, connID, this->maxConnPerWell_); +} + +// --------------------------------------------------------------------- + +class GroupVectors +{ +public: + template + using Window = boost::iterator_range; + + explicit GroupVectors(const RestartFileView& rst_view, + const ::Opm::RestartIO::ecl_kw_type* intehead); + + bool hasDefinedValues() const; + + std::size_t maxGroups() const; + + Window igrp(const std::size_t groupID) const; + Window xgrp(const std::size_t groupID) const; + +private: + std::size_t maxNumGroups_; + std::size_t numIGrpElem_; + std::size_t numXGrpElem_; + + const ::Opm::RestartIO::ecl_kw_type* igrp_; + const ::Opm::RestartIO::ecl_kw_type* xgrp_; +}; + +GroupVectors::GroupVectors(const RestartFileView& rst_view, + const ::Opm::RestartIO::ecl_kw_type* intehead) + : maxNumGroups_(getInteHeadElem(intehead, VI::intehead::NGMAXZ) - 1) // -FIELD + , numIGrpElem_ (getInteHeadElem(intehead, VI::intehead::NIGRPZ)) + , numXGrpElem_ (getInteHeadElem(intehead, VI::intehead::NXGRPZ)) + , igrp_ (rst_view.getKeyword("IGRP")) + , xgrp_ (rst_view.getKeyword("XGRP")) +{} + +bool GroupVectors::hasDefinedValues() const +{ + return ! ((this->igrp_ == nullptr) || + (this->xgrp_ == nullptr)); +} + +std::size_t GroupVectors::maxGroups() const +{ + return this->maxNumGroups_; +} + +GroupVectors::Window +GroupVectors::igrp(const std::size_t groupID) const +{ + if (! this->hasDefinedValues()) { + throw std::logic_error { + "Cannot Request IWEL Values Unless Defined" + }; + } + + return getDataWindow(getPtr(this->igrp_), + this->numIGrpElem_, groupID); +} + +GroupVectors::Window +GroupVectors::xgrp(const std::size_t groupID) const +{ + if (! this->hasDefinedValues()) { + throw std::logic_error { + "Cannot Request XGRP Values Unless Defined" + }; + } + + return getDataWindow(getPtr(this->xgrp_), + this->numXGrpElem_, groupID); +} + namespace { void throwIfMissingRequired(const Opm::RestartKey& rst_key) { @@ -348,101 +570,16 @@ namespace { return wells; } - template - const T* getPtr(const ::Opm::RestartIO::ecl_kw_type* kw) - { - return (kw == nullptr) ? nullptr - : static_cast(ecl_kw_iget_ptr(kw /* <- ADL */, 0)); - } - - int getInteHeadElem(const ::Opm::RestartIO::ecl_kw_type* intehead, - const std::vector::size_type i) - { - return getPtr(intehead)[i]; - } - - struct WellArrayDim - { - explicit WellArrayDim(const ::Opm::RestartIO::ecl_kw_type* intehead); - - std::size_t maxConnPerWell; - std::size_t numIWelElem; - std::size_t numXWelElem; - std::size_t numIConElm; - std::size_t numXConElm; - }; - - WellArrayDim::WellArrayDim(const ::Opm::RestartIO::ecl_kw_type* intehead) - : maxConnPerWell(getInteHeadElem(intehead, VI::intehead::NCWMAX)) - , numIWelElem (getInteHeadElem(intehead, VI::intehead::NIWELZ)) - , numXWelElem (getInteHeadElem(intehead, VI::intehead::NXWELZ)) - , numIConElm (getInteHeadElem(intehead, VI::intehead::NICONZ)) - , numXConElm (getInteHeadElem(intehead, VI::intehead::NXCONZ)) - {} - - template - boost::iterator_range - getDataWindow(const T* arr, - const std::size_t windowSize, - const std::size_t well, - const std::size_t conn = 0, - const std::size_t maxConnPerWell = 1) - { - const auto off = windowSize * (conn + maxConnPerWell*well); - const auto* begin = arr + off; - const auto* end = begin + windowSize; - - return { begin, end }; - } - - boost::iterator_range - getIWelWindow(const int* iwel, - const WellArrayDim& wdim, - const std::size_t well) - { - return getDataWindow(iwel, wdim.numIWelElem, well); - } - - boost::iterator_range - getXWelWindow(const double* xwel, - const WellArrayDim& wdim, - const std::size_t well) - { - return getDataWindow(xwel, wdim.numXWelElem, well); - } - - boost::iterator_range - getIConWindow(const int* icon, - const WellArrayDim& wdim, - const std::size_t well, - const std::size_t conn) - { - return getDataWindow(icon, wdim.numIConElm, well, - conn, wdim.maxConnPerWell); - } - - boost::iterator_range - getXConWindow(const double* xcon, - const WellArrayDim& wdim, - const std::size_t well, - const std::size_t conn) - { - return getDataWindow(xcon, wdim.numXConElm, well, - conn, wdim.maxConnPerWell); - } - std::unordered_map::size_type> - seqID_to_resID(const WellArrayDim& wdim, - const std::size_t wellID, + seqID_to_resID(const std::size_t wellID, const std::size_t nConn, - const int* icon_full) + const WellVectors& wellData) { using SizeT = boost::iterator_range::size_type; auto seqToRes = std::unordered_map{}; for (auto connID = 0*nConn; connID < nConn; ++connID) { - const auto icon = - getIConWindow(icon_full, wdim, wellID, connID); + const auto icon = wellData.icon(wellID, connID); seqToRes.emplace(icon[VI::IConn::index::SeqIndex] - 1, connID); } @@ -454,23 +591,20 @@ namespace { const std::size_t wellID, const std::size_t sim_step, const Opm::EclipseGrid& grid, - const WellArrayDim& wdim, const Opm::UnitSystem& usys, const Opm::Phases& phases, - const int* iwel_full, - const int* icon_full, - const double* xcon_full, + const WellVectors& wellData, Opm::data::Well& xw) { using M = ::Opm::UnitSystem::measure; - const auto iwel = getIWelWindow(iwel_full, wdim, wellID); + const auto iwel = wellData.iwel(wellID); const auto nConn = static_cast( iwel[VI::IWell::index::NConn]); xw.connections.resize(nConn, Opm::data::Connection{}); - if ((icon_full == nullptr) || (xcon_full == nullptr)) { + if (! wellData.hasDefinedConnectionValues()) { // Result set does not provide certain pieces of // information which are needed to reconstruct // connection flow rates. Nothing to do here. @@ -482,14 +616,12 @@ namespace { const auto wat = phases.active(Opm::Phase::WATER); const auto conns = well.getActiveConnections(sim_step, grid); - const auto seq_to_res = - seqID_to_resID(wdim, wellID, nConn, icon_full); + const auto seq_to_res = seqID_to_resID(wellID, nConn, wellData); auto linConnID = std::size_t{0}; for (const auto& conn : conns) { const auto connID = seq_to_res.at(conn.getSeqIndex()); - const auto xcon = - getXConWindow(xcon_full, wdim, wellID, connID); + const auto xcon = wellData.xcon(wellID, connID); auto& xc = xw.connections[linConnID++]; @@ -518,15 +650,11 @@ namespace { const std::size_t wellID, const std::size_t sim_step, const Opm::EclipseGrid& grid, - const WellArrayDim& wdim, const Opm::UnitSystem& usys, const Opm::Phases& phases, - const int* iwel_full, - const double* xwel_full, - const int* icon_full, - const double* xcon_full) + const WellVectors& wellData) { - if ((iwel_full == nullptr) || (xwel_full == nullptr)) { + if (! wellData.hasDefinedWellValues()) { // Result set does not provide well information. // No wells? In any case, nothing to do here. return {}; @@ -534,7 +662,7 @@ namespace { using M = ::Opm::UnitSystem::measure; - const auto xwel = getXWelWindow(xwel_full, wdim, wellID); + const auto xwel = wellData.xwel(wellID); const auto oil = phases.active(Opm::Phase::OIL); const auto gas = phases.active(Opm::Phase::GAS); @@ -566,8 +694,8 @@ namespace { xw.thp = xw.temperature = 0.0; // 3) Restore connection flow rates (xw.connections[i].rates) - restoreConnRates(well, wellID, sim_step, grid, wdim, usys, phases, - iwel_full, icon_full, xcon_full, xw); + restoreConnRates(well, wellID, sim_step, + grid, usys, phases, wellData, xw); return xw; } @@ -588,15 +716,11 @@ namespace { return soln; } - const auto wdim = WellArrayDim{ intehead }; + const auto wellData = WellVectors{ rst_view, intehead }; + const auto& units = es.getUnits(); const auto& phases = es.runspec().phases(); - const auto* iwel_full = getPtr (rst_view.getKeyword("IWEL")); - const auto* xwel_full = getPtr(rst_view.getKeyword("XWEL")); - const auto* icon_full = getPtr (rst_view.getKeyword("ICON")); - const auto* xcon_full = getPtr(rst_view.getKeyword("XCON")); - const auto sim_step = rst_view.simStep(); const auto& wells = schedule.getWells(sim_step); for (auto nWells = wells.size(), wellID = 0*nWells; @@ -605,8 +729,8 @@ namespace { const auto* well = wells[wellID]; soln[well->name()] = - restore_well(*well, wellID, sim_step, grid, wdim, units, phases, - iwel_full, xwel_full, icon_full, xcon_full); + restore_well(*well, wellID, sim_step, grid, + units, phases, wellData); } return soln; From 8f2a36593d83df34278ba1f2846680f012a102f3 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Fri, 2 Nov 2018 14:55:42 +0100 Subject: [PATCH 05/20] Load Restart: Load Segment Vectors if Available This commit extends RestartIO::load() to also pick up segment related quantities (flow rates SOFR, SGFR, SWFR, and pressure SPR) if these are available in the current restart file. --- src/opm/output/eclipse/LoadRestart.cpp | 152 ++++++++++++++++++++++++- 1 file changed, 149 insertions(+), 3 deletions(-) diff --git a/src/opm/output/eclipse/LoadRestart.cpp b/src/opm/output/eclipse/LoadRestart.cpp index 30032d06d..aa1116e43 100644 --- a/src/opm/output/eclipse/LoadRestart.cpp +++ b/src/opm/output/eclipse/LoadRestart.cpp @@ -28,6 +28,7 @@ #include #include +#include #include #include @@ -367,6 +368,77 @@ GroupVectors::xgrp(const std::size_t groupID) const this->numXGrpElem_, groupID); } +// --------------------------------------------------------------------- + +class SegmentVectors +{ +public: + template + using Window = boost::iterator_range; + + explicit SegmentVectors(const RestartFileView& rst_view, + const ::Opm::RestartIO::ecl_kw_type* intehead); + + bool hasDefinedValues() const; + + Window + iseg(const std::size_t mswID, const std::size_t segID) const; + + Window + rseg(const std::size_t mswID, const std::size_t segID) const; + +private: + std::size_t maxSegPerWell_; + std::size_t numISegElm_; + std::size_t numRSegElm_; + + const ::Opm::RestartIO::ecl_kw_type* iseg_; + const ::Opm::RestartIO::ecl_kw_type* rseg_; +}; + +SegmentVectors::SegmentVectors(const RestartFileView& rst_view, + const ::Opm::RestartIO::ecl_kw_type* intehead) + : maxSegPerWell_(getInteHeadElem(intehead, VI::intehead::NSEGMX)) + , numISegElm_ (getInteHeadElem(intehead, VI::intehead::NISEGZ)) + , numRSegElm_ (getInteHeadElem(intehead, VI::intehead::NRSEGZ)) + , iseg_ (rst_view.getKeyword("ISEG")) + , rseg_ (rst_view.getKeyword("RSEG")) +{} + +bool SegmentVectors::hasDefinedValues() const +{ + return ! ((this->iseg_ == nullptr) || + (this->rseg_ == nullptr)); +} + +SegmentVectors::Window +SegmentVectors::iseg(const std::size_t mswID, const std::size_t segID) const +{ + if (! this->hasDefinedValues()) { + throw std::logic_error { + "Cannot Request ISEG Values Unless Defined" + }; + } + + return getDataWindow(getPtr(this->iseg_), this->numISegElm_, + mswID, segID, this->maxSegPerWell_); +} + +SegmentVectors::Window +SegmentVectors::rseg(const std::size_t mswID, const std::size_t segID) const +{ + if (! this->hasDefinedValues()) { + throw std::logic_error { + "Cannot Request RSEG Values Unless Defined" + }; + } + + return getDataWindow(getPtr(this->rseg_), this->numRSegElm_, + mswID, segID, this->maxSegPerWell_); +} + +// --------------------------------------------------------------------- + namespace { void throwIfMissingRequired(const Opm::RestartKey& rst_key) { @@ -645,6 +717,64 @@ namespace { } } + void restoreSegmentQuantities(const std::size_t mswID, + const std::size_t numSeg, + const Opm::UnitSystem& usys, + const Opm::Phases& phases, + const SegmentVectors& segData, + Opm::data::Well& xw) + { + // Note sign of flow rates. RSEG stores flow rates as positive from + // reservoir to well--i.e., towards the producer/platform. The Flow + // simulator uses the opposite sign convention. + + using M = ::Opm::UnitSystem::measure; + + const auto oil = phases.active(Opm::Phase::OIL); + const auto gas = phases.active(Opm::Phase::GAS); + const auto wat = phases.active(Opm::Phase::WATER); + + // Renormalisation constant for gas okay in non-field unit systems. + // A bit more questionable for field units. + const auto watRenormalisation = 10.0; + const auto gasRenormalisation = 1000.0; + + for (auto segID = 0*numSeg; segID < numSeg; ++segID) { + const auto iseg = segData.iseg(mswID, segID); + const auto rseg = segData.rseg(mswID, segID); + + const auto segNumber = iseg[VI::ISeg::index::SegNo]; // One-based + + auto& segment = xw.segments[segNumber]; + + segment.segNumber = segNumber; + segment.pressure = + usys.to_si(M::pressure, rseg[VI::RSeg::index::Pressure]); + + const auto totFlow = rseg[VI::RSeg::index::TotFlowRate]; + const auto watFraction = rseg[VI::RSeg::index::WatFlowFract]; + const auto gasFraction = rseg[VI::RSeg::index::GasFlowFract]; + + if (wat) { + const auto qW = totFlow * watFraction * watRenormalisation; + segment.rates.set(Opm::data::Rates::opt::wat, + - usys.to_si(M::liquid_surface_rate, qW)); + } + + if (oil) { + const auto qO = totFlow * (1.0 - (watFraction + gasFraction)); + segment.rates.set(Opm::data::Rates::opt::oil, + - usys.to_si(M::liquid_surface_rate, qO)); + } + + if (gas) { + const auto qG = totFlow * gasFraction * gasRenormalisation; + segment.rates.set(Opm::data::Rates::opt::gas, + - usys.to_si(M::gas_surface_rate, qG)); + } + } + } + Opm::data::Well restore_well(const Opm::Well& well, const std::size_t wellID, @@ -652,7 +782,8 @@ namespace { const Opm::EclipseGrid& grid, const Opm::UnitSystem& usys, const Opm::Phases& phases, - const WellVectors& wellData) + const WellVectors& wellData, + const SegmentVectors& segData) { if (! wellData.hasDefinedWellValues()) { // Result set does not provide well information. @@ -697,6 +828,20 @@ namespace { restoreConnRates(well, wellID, sim_step, grid, usys, phases, wellData, xw); + // 4) Restore segment quantities if applicable. + if (well.isMultiSegment(sim_step) && + segData.hasDefinedValues()) + { + const auto iwel = wellData.iwel(wellID); + const auto mswID = iwel[VI::IWell::index::MsWID]; // One-based + const auto numSeg = iwel[VI::IWell::index::NWseg]; + + if ((mswID > 0) && (numSeg > 0)) { + restoreSegmentQuantities(mswID - 1, numSeg, usys, + phases, segData, xw); + } + } + return xw; } @@ -716,7 +861,8 @@ namespace { return soln; } - const auto wellData = WellVectors{ rst_view, intehead }; + const auto wellData = WellVectors { rst_view, intehead }; + const auto segData = SegmentVectors{ rst_view, intehead }; const auto& units = es.getUnits(); const auto& phases = es.runspec().phases(); @@ -730,7 +876,7 @@ namespace { soln[well->name()] = restore_well(*well, wellID, sim_step, grid, - units, phases, wellData); + units, phases, wellData, segData); } return soln; From bed76c330b97bd440c0007406ab0b31b27745164 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Fri, 2 Nov 2018 15:45:33 +0100 Subject: [PATCH 06/20] Summary: Add Means of Resetting Cumulative Quantities This commit introduces a new mutating operation on Summary objects, Summary::reset_cumulative_quantities(const SummaryState&) which overwrites the values of cumulative ("total") summary quantities with those of the input argument. The only *intended* use case is reinitialising cumulative quantities in the case of simulation restart, but other uses may exist. Add a test case to exercise the new interface. --- opm/output/eclipse/Summary.hpp | 2 + src/opm/output/eclipse/Summary.cpp | 19 +++++ tests/test_Summary.cpp | 117 +++++++++++++++++++++++++++++ 3 files changed, 138 insertions(+) diff --git a/opm/output/eclipse/Summary.hpp b/opm/output/eclipse/Summary.hpp index 5709e2bc4..1b52c8cfa 100644 --- a/opm/output/eclipse/Summary.hpp +++ b/opm/output/eclipse/Summary.hpp @@ -63,6 +63,8 @@ class Summary { const SummaryState& get_restart_vectors() const; + void reset_cumulative_quantities(const SummaryState& rstrt); + private: class keyword_handlers; diff --git a/src/opm/output/eclipse/Summary.cpp b/src/opm/output/eclipse/Summary.cpp index 526ef53ad..0b70db34f 100644 --- a/src/opm/output/eclipse/Summary.cpp +++ b/src/opm/output/eclipse/Summary.cpp @@ -1527,4 +1527,23 @@ const SummaryState& Summary::get_restart_vectors() const return this->prev_state; } +void Summary::reset_cumulative_quantities(const SummaryState& rstrt) +{ + for (const auto& f : this->handlers->handlers) { + if (! smspec_node_is_total(f.first)) { + // Ignore quantities that are not cumulative ("total"). + continue; + } + + const auto* genkey = smspec_node_get_gen_key1(f.first); + if (rstrt.has(genkey)) { + // Assume 'rstrt' uses output units. This is satisfied if rstrt + // is constructed from information in a restart file--i.e., from + // the double precision restart vectors 'XGRP' and 'XWEL' during + // RestartIO::load(). + this->prev_state.add(genkey, rstrt.get(genkey)); + } + } +} + }} // namespace Opm::out diff --git a/tests/test_Summary.cpp b/tests/test_Summary.cpp index 9dcf2cbe3..53ea476b6 100644 --- a/tests/test_Summary.cpp +++ b/tests/test_Summary.cpp @@ -2899,3 +2899,120 @@ BOOST_AUTO_TEST_CASE(Write_Read) } BOOST_AUTO_TEST_SUITE_END() + +// ===================================================================== + +BOOST_AUTO_TEST_SUITE(Reset_Cumulative_Vectors) + +BOOST_AUTO_TEST_CASE(Reset) +{ + const auto config = setup { "test.Reset.Cumulative" }; + ::Opm::out::Summary smry { + config.es, config.config, config.grid, + config.schedule, "Ignore.This" + }; + + auto rstrt = ::Opm::SummaryState{}; + rstrt.add("WOPT:W_1", 1.0); + rstrt.add("WWPT:W_1", 2.0); + rstrt.add("WGPT:W_1", 3.0); + rstrt.add("WVPT:W_1", 4.0); + + rstrt.add("WWIT:W_1", 5.0); + rstrt.add("WGIT:W_1", 6.0); + + rstrt.add("GOPT:NoSuchGroup", 1.0); + rstrt.add("GWPT:NoSuchGroup", 2.0); + rstrt.add("GGPT:NoSuchGroup", 3.0); + rstrt.add("GVPT:NoSuchGroup", 4.0); + + rstrt.add("GWIT:NoSuchGroup", 5.0); + rstrt.add("GGIT:NoSuchGroup", 6.0); + + rstrt.add("FOPT", 10.0); + rstrt.add("FWPT", 20.0); + rstrt.add("FGPT", 30.0); + rstrt.add("FVPT", 40.0); + + rstrt.add("FWIT", 50.0); + rstrt.add("FGIT", 60.0); + + smry.reset_cumulative_quantities(rstrt); + + const auto& sumstate = smry.get_restart_vectors(); + + // Cumulatives don't affect rates, BHP, WCT, or GOR. + for (const auto* w : { "W_1", "W_2", "W_3" }) { + auto get = [w, &sumstate](const std::string& vector) { + return sumstate.get(vector + ':' + std::string(w)); + }; + + BOOST_CHECK_THROW(get("WWPR"), std::out_of_range); + BOOST_CHECK_THROW(get("WOPR"), std::out_of_range); + BOOST_CHECK_THROW(get("WGPR"), std::out_of_range); + BOOST_CHECK_THROW(get("WVPR"), std::out_of_range); + + BOOST_CHECK_THROW(get("WWIR"), std::out_of_range); + BOOST_CHECK_THROW(get("WGIR"), std::out_of_range); + + BOOST_CHECK_THROW(get("WBHP"), std::out_of_range); + BOOST_CHECK_THROW(get("WWCT"), std::out_of_range); + BOOST_CHECK_THROW(get("WGOR"), std::out_of_range); + } + + // Cumulatives reset for W_1. + BOOST_CHECK_CLOSE(sumstate.get("WOPT:W_1"), 1.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WWPT:W_1"), 2.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WGPT:W_1"), 3.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WVPT:W_1"), 4.0, 1.0e-10); + + BOOST_CHECK_CLOSE(sumstate.get("WWIT:W_1"), 5.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WGIT:W_1"), 6.0, 1.0e-10); + + // Cumulatives unset for W_2. + BOOST_CHECK_CLOSE(sumstate.get("WOPT:W_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WWPT:W_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WGPT:W_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WVPT:W_2"), 0.0, 1.0e-10); + + BOOST_CHECK_CLOSE(sumstate.get("WWIT:W_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WGIT:W_2"), 0.0, 1.0e-10); + + // Cumulatives unset for W_3. + BOOST_CHECK_CLOSE(sumstate.get("WOPT:W_3"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WWPT:W_3"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WGPT:W_3"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WVPT:W_3"), 0.0, 1.0e-10); + + BOOST_CHECK_CLOSE(sumstate.get("WWIT:W_3"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WGIT:W_3"), 0.0, 1.0e-10); + + // Cumulatives unset for G_1. + BOOST_CHECK_CLOSE(sumstate.get("GOPT:G_1"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GWPT:G_1"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GGPT:G_1"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GVPT:G_1"), 0.0, 1.0e-10); + + BOOST_CHECK_CLOSE(sumstate.get("GWIT:G_1"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GGIT:G_1"), 0.0, 1.0e-10); + + // Cumulatives unset for G_2. + BOOST_CHECK_CLOSE(sumstate.get("GOPT:G_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GWPT:G_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GGPT:G_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GVPT:G_2"), 0.0, 1.0e-10); + + BOOST_CHECK_CLOSE(sumstate.get("GWIT:G_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GGIT:G_2"), 0.0, 1.0e-10); + + // Cumulatives reset for FIELD. + BOOST_CHECK_CLOSE(sumstate.get("FOPT"), 10.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("FWPT"), 20.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("FGPT"), 30.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("FVPT"), 40.0, 1.0e-10); + + BOOST_CHECK_CLOSE(sumstate.get("FWIT"), 50.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("FGIT"), 60.0, 1.0e-10); +} + +BOOST_AUTO_TEST_SUITE_END() From b3cfb4ee78e36838f3661dee18a79e3b0bd2b80b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Fri, 2 Nov 2018 15:51:41 +0100 Subject: [PATCH 07/20] Load Restart: Add Means of Restoring Cumulative Quantities This commit introduces a new helper function restore_cumulative() which creates a SummaryState object containing a subset of known cumulative quantities--notably cumulative phase and reservoir voidage production and cumulative water and gas injection for individula wells and groups--including FIELD. This is the main facility for resetting the simulator's notion of cumulative production following simulation restart. --- src/opm/output/eclipse/LoadRestart.cpp | 112 +++++++++++++++++++++++++ 1 file changed, 112 insertions(+) diff --git a/src/opm/output/eclipse/LoadRestart.cpp b/src/opm/output/eclipse/LoadRestart.cpp index aa1116e43..62e72a615 100644 --- a/src/opm/output/eclipse/LoadRestart.cpp +++ b/src/opm/output/eclipse/LoadRestart.cpp @@ -27,6 +27,7 @@ #include #include +#include #include #include #include @@ -37,6 +38,7 @@ #include #include #include +#include #include #include @@ -881,6 +883,116 @@ namespace { return soln; } + + void assign_well_cumulatives(const std::string& well, + const std::size_t wellID, + const WellVectors& wellData, + Opm::SummaryState& smry) + { + if (! wellData.hasDefinedWellValues()) { + // Result set does not provide well information. + // No wells? In any case, nothing to do here. + return; + } + + auto key = [&well](const std::string& vector) -> std::string + { + return vector + ':' + well; + }; + + const auto xwel = wellData.xwel(wellID); + + // No unit conversion here. Smry expects output units. + smry.add(key("WOPT"), xwel[VI::XWell::index::OilPrTotal]); + smry.add(key("WWPT"), xwel[VI::XWell::index::WatPrTotal]); + smry.add(key("WGPT"), xwel[VI::XWell::index::GasPrTotal]); + smry.add(key("WVPT"), xwel[VI::XWell::index::VoidPrTotal]); + + smry.add(key("WWIT"), xwel[VI::XWell::index::WatInjTotal]); + smry.add(key("WGIT"), xwel[VI::XWell::index::GasInjTotal]); + } + + void assign_group_cumulatives(const std::string& group, + const std::size_t groupID, + const GroupVectors& groupData, + Opm::SummaryState& smry) + { + if (! groupData.hasDefinedValues()) { + // Result set does not provide group information. + // No wells? In any case, nothing to do here. + return; + } + + auto key = [&group](const std::string& vector) -> std::string + { + return (group == "FIELD") + ? 'F' + vector + : 'G' + vector + ':' + group; + }; + + const auto xgrp = groupData.xgrp(groupID); + + // No unit conversion here. Smry expects output units. + smry.add(key("OPT"), xgrp[VI::XGroup::index::OilPrTotal]); + smry.add(key("WPT"), xgrp[VI::XGroup::index::WatPrTotal]); + smry.add(key("GPT"), xgrp[VI::XGroup::index::GasPrTotal]); + smry.add(key("VPT"), xgrp[VI::XGroup::index::VoidPrTotal]); + + smry.add(key("WIT"), xgrp[VI::XGroup::index::WatInjTotal]); + smry.add(key("GIT"), xgrp[VI::XGroup::index::GasInjTotal]); + } + + Opm::SummaryState + restore_cumulative(const RestartFileView& rst_view, + const ::Opm::Schedule& schedule) + { + auto smry = Opm::SummaryState{}; + + const auto sim_step = rst_view.simStep(); + const auto* intehead = rst_view.getKeyword("INTEHEAD"); + + if (intehead == nullptr) { return smry; } + + // Well cumulatives + { + const auto wellData = WellVectors { rst_view, intehead }; + const auto& wells = schedule.getWells(sim_step); + + for (auto nWells = wells.size(), wellID = 0*nWells; + wellID < nWells; ++wellID) + { + assign_well_cumulatives(wells[wellID]->name(), + wellID, wellData, smry); + } + } + + // Group cumulatives, including FIELD. + { + const auto groupData = GroupVectors { rst_view, intehead }; + + for (const auto* group : schedule.getGroups(sim_step)) { + const auto& gname = group->name(); + + // Note: Order of group values in {I,X}GRP arrays mostly + // matches group's order of occurrence in .DATA file. + // Values pertaingin to FIELD are stored at zero-based order + // index NGMAXZ (maximum number of groups in model). The + // latter value is groupData.maxGroups(). + // + // As a final wrinkle, Flow internally stores FIELD at + // seqIndex() == 0, so subtract one to account for this + // fact. Max(seqIndex(), 1) - 1 is just a bit of future + // proofing and robustness in case that ever changes. + const auto groupOrderIx = (gname == "FIELD") + ? groupData.maxGroups() // NGMAXZ -- Item 3 of WELLDIMS + : std::max(group->seqIndex(), std::size_t{1}) - 1; + + assign_group_cumulatives(gname, groupOrderIx, groupData, smry); + } + } + + return smry; + } } // Anonymous namespace namespace Opm { namespace RestartIO { From 9e97dde46670f41d191526ccf997b1154f7e6e48 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Fri, 2 Nov 2018 15:57:48 +0100 Subject: [PATCH 08/20] Restore Cumulative Summary Quantities in Restart This commit hooks the new helper function restore_cumulative() up to the gateway restart function RestartIO::load() by changing the return type from "RestartValue" to std::pair with the pair's '.second' being restore_cumulative()'s return value. Update callers including unit tests accordingly, specifically such that the gateway function EclipseIO::loadRestart() internally resets its Summary object's cumulative quantities using 'load().second'. This is, strictly speaking, a violation of the "const" qualifier on EclipseIO::loadRestart(), but the language permits the usage because the 'impl' pointer in this case will be a constant pointer to a mutable 'Impl' object. --- opm/output/eclipse/RestartIO.hpp | 76 +++++++-------- src/opm/output/eclipse/EclipseIO.cpp | 9 +- src/opm/output/eclipse/LoadRestart.cpp | 8 +- tests/test_Restart.cpp | 122 ++++++++++++++++++++++++- 4 files changed, 170 insertions(+), 45 deletions(-) diff --git a/opm/output/eclipse/RestartIO.hpp b/opm/output/eclipse/RestartIO.hpp index 6359a1b17..741731d4a 100644 --- a/opm/output/eclipse/RestartIO.hpp +++ b/opm/output/eclipse/RestartIO.hpp @@ -1,4 +1,5 @@ /* + Copyright (c) 2018 Equinor ASA Copyright (c) 2016 Statoil ASA Copyright (c) 2013-2015 Andreas Lauser Copyright (c) 2013 SINTEF ICT, Applied Mathematics. @@ -23,9 +24,6 @@ #ifndef RESTART_IO_HPP #define RESTART_IO_HPP -#include -#include - #include #include #include @@ -35,24 +33,29 @@ #include #include +#include + #include #include #include #include +#include +#include +#include + namespace Opm { -class EclipseGrid; -class EclipseState; -class Phases; -class Schedule; -class SummaryState; + class EclipseGrid; + class EclipseState; + class Phases; + class Schedule; -namespace RestartIO { +} // namespace Opm /* - The two loose functions RestartIO::save() and RestartIO::load() can + The two free functions RestartIO::save() and RestartIO::load() can be used to save and load reservoir and well state from restart files. Observe that these functions 'just do it', i.e. the checking of which report step to load from, if output is enabled at all and @@ -69,39 +72,30 @@ namespace RestartIO { load("CASE.X0010" , 99 , ...) save("CASE.X0010" , 99 , ...) - will read and write to the file "CASE.X0010" - completely ignoring + will read from and write to the file "CASE.X0010" - completely ignoring the report step argument '99'. */ +namespace Opm { namespace RestartIO { -/*void save(const std::string& filename, - int report_step, - double seconds_elapsed, - data::Solution cells, - data::Wells wells, - const EclipseState& es, - const EclipseGrid& grid, - const Schedule& schedule, - std::map> extra_data = {}, - bool write_double = false); -*/ -void save(const std::string& filename, - int report_step, - double seconds_elapsed, - RestartValue value, - const EclipseState& es, - const EclipseGrid& grid, - const Schedule& schedule, - const SummaryState& sumState, - bool write_double = false); + void save(const std::string& filename, + int report_step, + double seconds_elapsed, + RestartValue value, + const EclipseState& es, + const EclipseGrid& grid, + const Schedule& schedule, + const SummaryState& sumState, + bool write_double = false); -RestartValue load( const std::string& filename, - int report_step, - const std::vector& solution_keys, - const EclipseState& es, - const EclipseGrid& grid, - const Schedule& schedule, - const std::vector& extra_keys = {}); + std::pair + load(const std::string& filename, + int report_step, + const std::vector& solution_keys, + const EclipseState& es, + const EclipseGrid& grid, + const Schedule& schedule, + const std::vector& extra_keys = {}); -} -} -#endif +}} // namespace Opm::RestartIO + +#endif // RESTART_IO_HPP diff --git a/src/opm/output/eclipse/EclipseIO.cpp b/src/opm/output/eclipse/EclipseIO.cpp index 9d5cb12bd..f89ee5f7b 100644 --- a/src/opm/output/eclipse/EclipseIO.cpp +++ b/src/opm/output/eclipse/EclipseIO.cpp @@ -500,7 +500,14 @@ RestartValue EclipseIO::loadRestart(const std::vector& solution_keys report_step, false ); - return RestartIO::load( filename , report_step , solution_keys , es, grid , schedule, extra_keys); + auto rst = RestartIO::load(filename, report_step, solution_keys, + es, grid, schedule, extra_keys); + + // Technically a violation of 'const'. Allowed because 'impl' is + // constant pointer to mutable Impl. + this->impl->summary.reset_cumulative_quantities(rst.second); + + return std::move(rst.first); } EclipseIO::EclipseIO( const EclipseState& es, diff --git a/src/opm/output/eclipse/LoadRestart.cpp b/src/opm/output/eclipse/LoadRestart.cpp index 62e72a615..2028d1de9 100644 --- a/src/opm/output/eclipse/LoadRestart.cpp +++ b/src/opm/output/eclipse/LoadRestart.cpp @@ -997,7 +997,7 @@ namespace { namespace Opm { namespace RestartIO { - RestartValue + std::pair load(const std::string& filename, int report_step, const std::vector& solution_keys, @@ -1023,6 +1023,10 @@ namespace Opm { namespace RestartIO { restoreExtra(rst_view, extra_keys, es.getUnits(), rst_value); } - return rst_value; + return { + std::move(rst_value), + restore_cumulative(rst_view, schedule) + }; } + }} // Opm::RestartIO diff --git a/tests/test_Restart.cpp b/tests/test_Restart.cpp index ec4c7df53..0fb26a47d 100644 --- a/tests/test_Restart.cpp +++ b/tests/test_Restart.cpp @@ -405,6 +405,36 @@ Opm::SummaryState sim_state() state.add("WGVIR:OP_3", 0.0); state.add("WWVIR:OP_3", 43.21); + state.add("GOPR:OP" , 110.0); + state.add("GWPR:OP" , 120.0); + state.add("GGPR:OP" , 130.0); + state.add("GVPR:OP" , 140.0); + state.add("GOPT:OP" , 1100.0); + state.add("GWPT:OP" , 1200.0); + state.add("GGPT:OP" , 1300.0); + state.add("GVPT:OP" , 1400.0); + state.add("GWIR:OP" , - 256.0); + state.add("GGIR:OP" , - 65536.0); + state.add("GWIT:OP" , 31415.9); + state.add("GGIT:OP" , 27182.8); + state.add("GWCT:OP" , 0.625); + state.add("GGOR:OP" , 1234.5); + + state.add("FOPR" , 1100.0); + state.add("FWPR" , 1200.0); + state.add("FGPR" , 1300.0); + state.add("FVPR" , 1400.0); + state.add("FOPT" , 11000.0); + state.add("FWPT" , 12000.0); + state.add("FGPT" , 13000.0); + state.add("FVPT" , 14000.0); + state.add("FWIR" , - 2560.0); + state.add("FGIR" , - 655360.0); + state.add("FWIT" , 314159.2); + state.add("FGIT" , 271828.1); + state.add("FWCT" , 0.625); + state.add("FGOR" , 1234.5); + return state; } @@ -688,7 +718,7 @@ BOOST_AUTO_TEST_CASE(ExtraData_content) { /* extra_keys = */ { {"EXTRA" , UnitSystem::measure::pressure, true} , {"EXTRA2", UnitSystem::measure::identity, false} - }); + }).first; BOOST_CHECK(!rst_value.hasExtra("EXTRA2")); BOOST_CHECK( rst_value.hasExtra("EXTRA")); @@ -777,4 +807,94 @@ BOOST_AUTO_TEST_CASE(STORE_THPRES) { +BOOST_AUTO_TEST_CASE(Restore_Cumulatives) +{ + Setup setup("FIRST_SIM.DATA"); + + // Write fully ECLIPSE compatible output. This also saves cumulatives. + setup.es.getIOConfig().setEclCompatibleRST(true); + + const auto restart_value = RestartValue { + mkSolution(setup.grid.getNumActive()), + mkWells() + }; + const auto sumState = sim_state(); + + RestartIO::save("FILE.UNRST", 1, 100, restart_value, + setup.es, setup.grid, setup.schedule, sumState); + + const auto rst_value = RestartIO::load("FILE.UNRST", 1, + /* solution_keys = */ { + RestartKey("SWAT", UnitSystem::measure::identity), + }, + setup.es, setup.grid, setup.schedule, + /* extra_keys = */ {}); + + const auto& rstSumState = rst_value.second; + + // Verify that the restored summary state has all of its requisite + // cumulative summary vectors. + + // Producer => W*IT saved/restored as zero (0.0) + BOOST_CHECK(rstSumState.has("WOPT:OP_1")); + BOOST_CHECK(rstSumState.has("WGPT:OP_1")); + BOOST_CHECK(rstSumState.has("WWPT:OP_1")); + BOOST_CHECK(rstSumState.has("WVPT:OP_1")); + BOOST_CHECK(rstSumState.has("WWIT:OP_1")); + BOOST_CHECK(rstSumState.has("WGIT:OP_1")); + + BOOST_CHECK_CLOSE(rstSumState.get("WOPT:OP_1"), 10.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WGPT:OP_1"), 30.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WWPT:OP_1"), 20.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WVPT:OP_1"), 40.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WWIT:OP_1"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WGIT:OP_1"), 0.0, 1.0e-10); + + // Gas injector => W*PT and WWIT saved/restored as zero (0.0) + BOOST_CHECK(rstSumState.has("WOPT:OP_2")); + BOOST_CHECK(rstSumState.has("WGPT:OP_2")); + BOOST_CHECK(rstSumState.has("WWPT:OP_2")); + BOOST_CHECK(rstSumState.has("WVPT:OP_2")); + BOOST_CHECK(rstSumState.has("WWIT:OP_2")); + BOOST_CHECK(rstSumState.has("WGIT:OP_2")); + + BOOST_CHECK_CLOSE(rstSumState.get("WOPT:OP_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WGPT:OP_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WWPT:OP_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WVPT:OP_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WWIT:OP_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WGIT:OP_2"), 2000.0, 1.0e-10); + + // Group cumulatives saved/restored for all phases + BOOST_CHECK(rstSumState.has("GOPT:OP")); + BOOST_CHECK(rstSumState.has("GGPT:OP")); + BOOST_CHECK(rstSumState.has("GWPT:OP")); + BOOST_CHECK(rstSumState.has("GVPT:OP")); + BOOST_CHECK(rstSumState.has("GWIT:OP")); + BOOST_CHECK(rstSumState.has("GGIT:OP")); + + BOOST_CHECK_CLOSE(rstSumState.get("GOPT:OP"), 1100.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("GWPT:OP"), 1200.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("GGPT:OP"), 1300.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("GVPT:OP"), 1400.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("GWIT:OP"), 31415.9, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("GGIT:OP"), 27182.8, 1.0e-10); + + // Field cumulatives saved/restored for all phases + BOOST_CHECK(rstSumState.has("FOPT")); + BOOST_CHECK(rstSumState.has("FGPT")); + BOOST_CHECK(rstSumState.has("FWPT")); + BOOST_CHECK(rstSumState.has("FVPT")); + BOOST_CHECK(rstSumState.has("FWIT")); + BOOST_CHECK(rstSumState.has("FGIT")); + + BOOST_CHECK_CLOSE(rstSumState.get("FOPT"), 11000.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("FWPT"), 12000.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("FGPT"), 13000.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("FVPT"), 14000.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("FWIT"), 314159.2, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("FGIT"), 271828.1, 1.0e-10); +} + + } From 551f1f8c76ba516802110d0454d6f4031d95a063 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Fri, 16 Nov 2018 12:31:46 +0100 Subject: [PATCH 09/20] Connection Rates: Lookup Connections by (I,J,K) Index This is more robust than the linear index. --- src/opm/output/eclipse/LoadRestart.cpp | 111 +++++++++++++++++-------- 1 file changed, 76 insertions(+), 35 deletions(-) diff --git a/src/opm/output/eclipse/LoadRestart.cpp b/src/opm/output/eclipse/LoadRestart.cpp index 2028d1de9..3b57f80a7 100644 --- a/src/opm/output/eclipse/LoadRestart.cpp +++ b/src/opm/output/eclipse/LoadRestart.cpp @@ -45,9 +45,10 @@ #include #include #include +#include #include #include -#include +#include #include #include @@ -644,21 +645,66 @@ namespace { return wells; } - std::unordered_map::size_type> - seqID_to_resID(const std::size_t wellID, - const std::size_t nConn, - const WellVectors& wellData) + std::map< + std::tuple, + boost::iterator_range::size_type + > + ijk_to_resID(const std::size_t wellID, + const std::size_t nConn, + const WellVectors& wellData) { using SizeT = boost::iterator_range::size_type; - auto seqToRes = std::unordered_map{}; + auto ijkToRes = std::map, SizeT>{}; for (auto connID = 0*nConn; connID < nConn; ++connID) { const auto icon = wellData.icon(wellID, connID); - seqToRes.emplace(icon[VI::IConn::index::SeqIndex] - 1, connID); + const auto i = icon[VI::IConn::index::CellI] - 1; + const auto j = icon[VI::IConn::index::CellJ] - 1; + const auto k = icon[VI::IConn::index::CellK] - 1; + + ijkToRes.emplace(std::make_tuple(i, j, k), connID); } - return seqToRes; + return ijkToRes; + } + + void restoreConnRates(const WellVectors::Window& xcon, + const Opm::UnitSystem& usys, + const bool oil, + const bool gas, + const bool wat, + Opm::data::Connection& xc) + { + using M = ::Opm::UnitSystem::measure; + + if (wat) { + xc.rates.set(Opm::data::Rates::opt::wat, + - usys.to_si(M::liquid_surface_rate, + xcon[VI::XConn::index::WaterRate])); + } + + if (oil) { + xc.rates.set(Opm::data::Rates::opt::oil, + - usys.to_si(M::liquid_surface_rate, + xcon[VI::XConn::index::OilRate])); + } + + if (gas) { + xc.rates.set(Opm::data::Rates::opt::gas, + - usys.to_si(M::gas_surface_rate, + xcon[VI::XConn::index::GasRate])); + } + } + + void zeroConnRates(const bool oil, + const bool gas, + const bool wat, + Opm::data::Connection& xc) + { + if (wat) { xc.rates.set(Opm::data::Rates::opt::wat, 0.0); } + if (oil) { xc.rates.set(Opm::data::Rates::opt::oil, 0.0); } + if (gas) { xc.rates.set(Opm::data::Rates::opt::gas, 0.0); } } void restoreConnRates(const Opm::Well& well, @@ -670,51 +716,46 @@ namespace { const WellVectors& wellData, Opm::data::Well& xw) { - using M = ::Opm::UnitSystem::measure; - const auto iwel = wellData.iwel(wellID); const auto nConn = static_cast( iwel[VI::IWell::index::NConn]); xw.connections.resize(nConn, Opm::data::Connection{}); - if (! wellData.hasDefinedConnectionValues()) { - // Result set does not provide certain pieces of - // information which are needed to reconstruct - // connection flow rates. Nothing to do here. - return; - } - const auto oil = phases.active(Opm::Phase::OIL); const auto gas = phases.active(Opm::Phase::GAS); const auto wat = phases.active(Opm::Phase::WATER); + for (auto& xc : xw.connections) { + zeroConnRates(oil, gas, wat, xc); + } + + if (! wellData.hasDefinedConnectionValues()) { + // Result set does not provide certain pieces of information + // which are needed to reconstruct connection flow rates. + // Nothing to do except to return all zero rates. + return; + } + const auto conns = well.getActiveConnections(sim_step, grid); - const auto seq_to_res = seqID_to_resID(wellID, nConn, wellData); + const auto ijk_to_res = ijk_to_resID(wellID, nConn, wellData); auto linConnID = std::size_t{0}; for (const auto& conn : conns) { - const auto connID = seq_to_res.at(conn.getSeqIndex()); - const auto xcon = wellData.xcon(wellID, connID); + if (++linConnID > nConn) { continue; } - auto& xc = xw.connections[linConnID++]; + auto& xc = xw.connections[linConnID - 1]; - if (wat) { - xc.rates.set(Opm::data::Rates::opt::wat, - - usys.to_si(M::liquid_surface_rate, - xcon[VI::XConn::index::WaterRate])); - } + const auto ijk = + std::make_tuple(conn.getI(), conn.getJ(), conn.getK()); - if (oil) { - xc.rates.set(Opm::data::Rates::opt::oil, - - usys.to_si(M::liquid_surface_rate, - xcon[VI::XConn::index::OilRate])); - } + auto resPos = ijk_to_res.find(ijk); - if (gas) { - xc.rates.set(Opm::data::Rates::opt::gas, - - usys.to_si(M::gas_surface_rate, - xcon[VI::XConn::index::GasRate])); + if (resPos != std::end(ijk_to_res)) { + const auto connID = resPos->second; + const auto xcon = wellData.xcon(wellID, connID); + + restoreConnRates(xcon, usys, oil, gas, wat, xc); } } } From ecca78af851e7928fbdcae64674be42c96d7eb72 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Tue, 4 Dec 2018 14:15:13 +0100 Subject: [PATCH 10/20] DOUBHEAD: Add OPM-Specific Method for Storing Next Timestep This commit creates an OPM-specific extension of the DOUBHEAD vector of a restart step. We reuse the TSINIT item (zero-based index 1) to store the next timestep. Local testing suggests that ECLIPSE does not use this value as part of restarting a simulation so this item is a reasonable compromise for creating a mostly ECLIPSE-compatible restart file that still enables communicating the suggested next timestep if we're restarting Flow from a result set created by Flow. --- opm/output/eclipse/DoubHEAD.hpp | 1 + opm/output/eclipse/VectorItems/doubhead.hpp | 37 +++++++++++++++++++++ src/opm/output/eclipse/DoubHEAD.cpp | 23 ++++++++++--- 3 files changed, 56 insertions(+), 5 deletions(-) create mode 100644 opm/output/eclipse/VectorItems/doubhead.hpp diff --git a/opm/output/eclipse/DoubHEAD.hpp b/opm/output/eclipse/DoubHEAD.hpp index 7f6b51643..c4348a12f 100755 --- a/opm/output/eclipse/DoubHEAD.hpp +++ b/opm/output/eclipse/DoubHEAD.hpp @@ -53,6 +53,7 @@ namespace Opm { namespace RestartIO { const double cnvT); DoubHEAD& timeStamp(const TimeStamp& ts); + DoubHEAD& nextStep(const double nextTimeStep); DoubHEAD& drsdt(const Schedule& sched, const std::size_t lookup_step, diff --git a/opm/output/eclipse/VectorItems/doubhead.hpp b/opm/output/eclipse/VectorItems/doubhead.hpp new file mode 100644 index 000000000..36a70b8d2 --- /dev/null +++ b/opm/output/eclipse/VectorItems/doubhead.hpp @@ -0,0 +1,37 @@ +/* + Copyright (c) 2018 Equinor ASA + + 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 + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + +#ifndef OPM_OUTPUT_ECLIPSE_VECTOR_DOUBHEAD_HPP +#define OPM_OUTPUT_ECLIPSE_VECTOR_DOUBHEAD_HPP + +#include + +namespace Opm { namespace RestartIO { namespace Helpers { namespace VectorItems { + + // This is a subset of the items in src/opm/output/eclipse/DoubHEAD.cpp . + // Promote items from that list to this in order to make them public. + enum doubhead : std::vector::size_type { + TsInit = 1, // Maximum Length of Next Timestep + TsMaxz = 2, // Maximum Length of Timestep After Next + TsMinz = 3, // Minumum Length of All Timesteps + }; + +}}}} // Opm::RestartIO::Helpers::VectorItems + +#endif // OPM_OUTPUT_ECLIPSE_VECTOR_DOUBHEAD_HPP diff --git a/src/opm/output/eclipse/DoubHEAD.cpp b/src/opm/output/eclipse/DoubHEAD.cpp index 8868902d8..b8e2ccb1c 100755 --- a/src/opm/output/eclipse/DoubHEAD.cpp +++ b/src/opm/output/eclipse/DoubHEAD.cpp @@ -28,6 +28,8 @@ #include // Opm::RestartIO::makeUTCTime() +#include + #include #include #include @@ -38,12 +40,14 @@ #include #include +namespace VI = Opm::RestartIO::Helpers::VectorItems; + enum Index : std::vector::size_type { // 0..9 SimTime = 0, - TsInit = 1, - TsMaxz = 2, - TsMinz = 3, + TsInit = VI::doubhead::TsInit, + TsMaxz = VI::doubhead::TsMaxz, + TsMinz = VI::doubhead::TsMinz, TsMchp = 4, TsFMax = 5, TsFMin = 6, @@ -303,7 +307,7 @@ enum Index : std::vector::size_type { dh_218 = 218, dh_219 = 219, - // 220..227 + // 220..228 dh_220 = 220, dh_221 = 221, dh_222 = 222, @@ -376,7 +380,6 @@ namespace { Opm::RestartIO::DoubHEAD::DoubHEAD() : data_(Index::NUMBER_OF_ITEMS, 0.0) - //: data_(Index::NUMBER_OF_ITEMS, -1.0e20) { // Numbers below have unknown usage, values have been determined by // experiments to be constant across a range of reference cases. @@ -588,6 +591,16 @@ Opm::RestartIO::DoubHEAD::timeStamp(const TimeStamp& ts) return *this; } +Opm::RestartIO::DoubHEAD& +Opm::RestartIO::DoubHEAD::nextStep(const double nextTimeStep) +{ + if (nextTimeStep > 0.0) { + this->data_[Index::TsInit] = nextTimeStep; + } + + return *this; +} + Opm::RestartIO::DoubHEAD& Opm::RestartIO::DoubHEAD::drsdt(const Schedule& sched, const std::size_t lookup_step, From 466693f1d06514add01a1ba0ae6a8215481d047e Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Tue, 4 Dec 2018 16:44:11 +0100 Subject: [PATCH 11/20] Load Restart: Add Means of Retrieving Next Timestep From DOUBHEAD This commit restores Flow's suggested next timestep size from the TSINIT item of the DOUBHEAD vector pertaining to a particular restart/report step. If the item is defaulted, which typically happens if the result set is created by ECLIPSE, then we don't try to restore the stepsize and RestartValue::hasExtra("OPMEXTRA") will be false. --- src/opm/output/eclipse/LoadRestart.cpp | 58 ++++++++++++++++++++++---- 1 file changed, 50 insertions(+), 8 deletions(-) diff --git a/src/opm/output/eclipse/LoadRestart.cpp b/src/opm/output/eclipse/LoadRestart.cpp index 3b57f80a7..3650a950c 100644 --- a/src/opm/output/eclipse/LoadRestart.cpp +++ b/src/opm/output/eclipse/LoadRestart.cpp @@ -27,6 +27,7 @@ #include #include +#include #include #include #include @@ -474,6 +475,25 @@ namespace { } } + std::vector + getOpmExtraFromDoubHEAD(const RestartFileView& rst_view, + const bool required, + const Opm::UnitSystem& usys) + { + using M = Opm::UnitSystem::measure; + + const auto* doubhead = + getPtr(rst_view.getKeyword("DOUBHEAD")); + + const auto TsInit = doubhead[VI::doubhead::TsInit]; + + if (TsInit < 0.0) { + throwIfMissingRequired({ "OPMEXTRA", M::identity, required }); + } + + return { usys.to_si(M::time, TsInit) }; + } + Opm::data::Solution restoreSOLUTION(const RestartFileView& rst_view, const std::vector& solution_keys, @@ -518,17 +538,39 @@ namespace { const auto& vector = extra.key; const auto* kw = rst_view.getKeyword(vector.c_str()); - if (kw == nullptr) { - throwIfMissingRequired(extra); + auto kwdata = std::vector{}; - // Requested vector not available, but caller does not - // actually require the vector for restart purposes. - // Skip this. - continue; + if (kw == nullptr) { + // Requested vector not available in result set. Take + // appropriate action depending on specific vector and + // 'extra.required'. + + if (vector != "OPMEXTRA") { + throwIfMissingRequired(extra); + + // Requested vector not available, but caller does not + // actually require the vector for restart purposes. + // Skip this. + continue; + } + else { + // Special case handling of OPMEXTRA. Single item + // possibly stored in TSINIT item of DOUBHEAD. Try to + // recover this. Function throws if item is defaulted + // and caller requires that item be present through the + // 'extra.required' mechanism. + + kwdata = getOpmExtraFromDoubHEAD(rst_view, + extra.required, + usys); + } + } + else { + // Requisite vector available in result set. Recover data. + kwdata = double_vector(kw); } - // Requisite vector available in result set. Recover data. - rst_value.addExtra(vector, extra.dim, double_vector(kw)); + rst_value.addExtra(vector, extra.dim, std::move(kwdata)); } for (auto& extra_value : rst_value.extra) { From bfe25071843ff72acf5efdc91bb2074618e7f554 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Wed, 5 Dec 2018 10:23:35 +0100 Subject: [PATCH 12/20] Restart Facility: Output Suggested Next Timestep as TSINIT This commit activates the support for storing Flow's suggested next timestep size as the TSINIT item (zero-based index 1) of the restart file's DOUBHEAD vector. Local testing suggests that this value is essential to even being able to restart Flow from a result set generated by Flow itself. --- opm/output/eclipse/WriteRestartHelpers.hpp | 3 ++- src/opm/output/eclipse/CreateDoubHead.cpp | 22 +++++++++++++++------ src/opm/output/eclipse/RestartIO.cpp | 23 ++++++++++++++++------ 3 files changed, 35 insertions(+), 13 deletions(-) diff --git a/opm/output/eclipse/WriteRestartHelpers.hpp b/opm/output/eclipse/WriteRestartHelpers.hpp index fae90be5b..ba0eabcc0 100755 --- a/opm/output/eclipse/WriteRestartHelpers.hpp +++ b/opm/output/eclipse/WriteRestartHelpers.hpp @@ -47,7 +47,8 @@ namespace Opm { namespace RestartIO { namespace Helpers { createDoubHead(const EclipseState& es, const Schedule& sched, const std::size_t lookup_step, - const double simTime); + const double simTime, + const double nextTimeStep); diff --git a/src/opm/output/eclipse/CreateDoubHead.cpp b/src/opm/output/eclipse/CreateDoubHead.cpp index abddcd7df..165dcd4f7 100755 --- a/src/opm/output/eclipse/CreateDoubHead.cpp +++ b/src/opm/output/eclipse/CreateDoubHead.cpp @@ -1,4 +1,5 @@ /* + Copyright (c) 2018 Equinor ASA Copyright (c) 2018 Statoil ASA This file is part of the Open Porous Media project (OPM). @@ -24,6 +25,7 @@ #include #include +#include #include #include @@ -76,15 +78,23 @@ Opm::RestartIO::Helpers:: createDoubHead(const EclipseState& es, const Schedule& sched, const std::size_t lookup_step, - const double simTime) + const double simTime, + const double nextTimeStep) { - const auto& usys = es.getDeckUnitSystem(); - const auto dh = DoubHEAD{} - .tuningParameters(sched.getTuning(), lookup_step, - getTimeConv(usys)) + const auto& usys = es.getDeckUnitSystem(); + const auto tconv = getTimeConv(usys); + + auto dh = DoubHEAD{} + .tuningParameters(sched.getTuning(), lookup_step, tconv) .timeStamp (computeTimeStamp(sched, simTime)) - .drsdt (sched, lookup_step, getTimeConv(usys)) + .drsdt (sched, lookup_step, tconv) ; + if (nextTimeStep > 0.0) { + using M = ::Opm::UnitSystem::measure; + + dh.nextStep(usys.from_si(M::time, nextTimeStep)); + } + return dh.data(); } diff --git a/src/opm/output/eclipse/RestartIO.cpp b/src/opm/output/eclipse/RestartIO.cpp index 72b9daaa8..401e2dbae 100644 --- a/src/opm/output/eclipse/RestartIO.cpp +++ b/src/opm/output/eclipse/RestartIO.cpp @@ -74,6 +74,13 @@ namespace { return extra_solution.count(vector) > 0; } + double nextStepSize(const Opm::RestartValue& rst_value) + { + return rst_value.hasExtra("OPMEXTRA") + ? rst_value.getExtra("OPMEXTRA")[0] + : 0.0; + } + std::vector serialize_OPM_IWEL(const data::Wells& wells, const std::vector& sched_wells) @@ -262,9 +269,10 @@ namespace { std::vector writeHeader(::Opm::RestartIO::ecl_rst_file_type* rst_file, - int sim_step, - int report_step, - double simTime, + const int sim_step, + const int report_step, + const double next_step_size, + const double simTime, const Schedule& schedule, const EclipseGrid& grid, const EclipseState& es) @@ -282,7 +290,8 @@ namespace { write_kw(rst_file, "LOGIHEAD", lh); // write DOUBHEAD to restart file - const auto dh = Helpers::createDoubHead(es, schedule, sim_step, simTime); + const auto dh = Helpers::createDoubHead(es, schedule, sim_step, + simTime, next_step_size); write_kw(rst_file, "DOUBHEAD", dh); // return the inteHead vector @@ -469,8 +478,10 @@ void save(const std::string& filename, units.from_si(restart_key.dim, data); } - const auto inteHD = writeHeader(rst_file.get(), sim_step, report_step, - seconds_elapsed, schedule, grid, es); + const auto inteHD = + writeHeader(rst_file.get(), sim_step, report_step, + nextStepSize(value), seconds_elapsed, + schedule, grid, es); writeGroup(rst_file.get(), sim_step, ecl_compatible_rst, schedule, sumState, inteHD); From b311fbb9d018d7468681f9f51108b0f1e8c35dc6 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Fri, 14 Dec 2018 20:34:32 +0100 Subject: [PATCH 13/20] Save/Restore Connection Index and Pressure Values Further analysis suggests Item 35 (index 34) of XCON is the connection pressure. Use this information to save/restore the values of data::Connection::pressure. While here, also fix a spelling error in LoadRestart.cpp. --- opm/output/eclipse/VectorItems/connection.hpp | 2 ++ .../eclipse/AggregateConnectionData.cpp | 2 ++ src/opm/output/eclipse/LoadRestart.cpp | 32 ++++++++++++------- 3 files changed, 25 insertions(+), 11 deletions(-) diff --git a/opm/output/eclipse/VectorItems/connection.hpp b/opm/output/eclipse/VectorItems/connection.hpp index ae725de1c..85e0efe29 100644 --- a/opm/output/eclipse/VectorItems/connection.hpp +++ b/opm/output/eclipse/VectorItems/connection.hpp @@ -66,6 +66,8 @@ namespace Opm { namespace RestartIO { namespace Helpers { namespace VectorItems WaterRate = 1, // Surface flow rate (water) GasRate = 2, // Surface Flow rate (gas) + Pressure = 34, // Connection pressure value + ResVRate = 49, // Reservoir voidage rate }; } // XConn diff --git a/src/opm/output/eclipse/AggregateConnectionData.cpp b/src/opm/output/eclipse/AggregateConnectionData.cpp index b67121be8..e336b0ff4 100755 --- a/src/opm/output/eclipse/AggregateConnectionData.cpp +++ b/src/opm/output/eclipse/AggregateConnectionData.cpp @@ -246,6 +246,8 @@ namespace { using Ix = ::Opm::RestartIO::Helpers::VectorItems::XConn::index; using R = ::Opm::data::Rates::opt; + xConn[Ix::Pressure] = units.from_si(M::pressure, x.pressure); + // Note flow rate sign. Treat production rates as positive. const auto& Q = x.rates; diff --git a/src/opm/output/eclipse/LoadRestart.cpp b/src/opm/output/eclipse/LoadRestart.cpp index 3650a950c..24f9cd7b8 100644 --- a/src/opm/output/eclipse/LoadRestart.cpp +++ b/src/opm/output/eclipse/LoadRestart.cpp @@ -749,15 +749,18 @@ namespace { if (gas) { xc.rates.set(Opm::data::Rates::opt::gas, 0.0); } } - void restoreConnRates(const Opm::Well& well, - const std::size_t wellID, - const std::size_t sim_step, - const Opm::EclipseGrid& grid, - const Opm::UnitSystem& usys, - const Opm::Phases& phases, - const WellVectors& wellData, - Opm::data::Well& xw) + void restoreConnResults(const Opm::Well& well, + const std::size_t wellID, + const std::size_t sim_step, + const Opm::EclipseGrid& grid, + const Opm::UnitSystem& usys, + const Opm::Phases& phases, + const WellVectors& wellData, + Opm::data::Well& xw) { + using M = ::Opm::UnitSystem::measure; + using Ix = ::Opm::RestartIO::Helpers::VectorItems::XConn::index; + const auto iwel = wellData.iwel(wellID); const auto nConn = static_cast( iwel[VI::IWell::index::NConn]); @@ -798,6 +801,12 @@ namespace { const auto xcon = wellData.xcon(wellID, connID); restoreConnRates(xcon, usys, oil, gas, wat, xc); + + xc.index = grid.getGlobalIndex(std::get<0>(ijk), + std::get<1>(ijk), + std::get<2>(ijk)); + + xc.pressure = usys.to_si(M::pressure, xcon[Ix::Pressure]); } } } @@ -910,8 +919,9 @@ namespace { xw.thp = xw.temperature = 0.0; // 3) Restore connection flow rates (xw.connections[i].rates) - restoreConnRates(well, wellID, sim_step, - grid, usys, phases, wellData, xw); + // and pressure values (xw.connections[i].pressure). + restoreConnResults(well, wellID, sim_step, + grid, usys, phases, wellData, xw); // 4) Restore segment quantities if applicable. if (well.isMultiSegment(sim_step) && @@ -1058,7 +1068,7 @@ namespace { // Note: Order of group values in {I,X}GRP arrays mostly // matches group's order of occurrence in .DATA file. - // Values pertaingin to FIELD are stored at zero-based order + // Values pertaining to FIELD are stored at zero-based order // index NGMAXZ (maximum number of groups in model). The // latter value is groupData.maxGroups(). // From 169ec766057c29db3ef38da0353d6710eb2e9d06 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Fri, 11 Jan 2019 08:22:48 +0100 Subject: [PATCH 14/20] Save/Restore: Support ECLIPSE Compatible Hysteresis Vectors This commit adds an intermediate layer to the save/restore code that will translate between Flow's hysteresis parameters and a compatible subset of ECLIPSE's hysteresis model parameters when creating a stricly ECLIPSE compatible restart file. In particular we save and restore Flow's KRNSW_OW and PCSWM_OW parameters in terms of the SOMAX parameter using the relations KRNSW_OW = 1 - SOMAX PCSWM_OW = 1 - SOMAX Similarly, we save and restore Flow's KRNSW_GO and PCSWM_GO parameters in terms of ECLIPSE's SGMAX parameter using the relations KRNSW_GO = 1 - SGMAX PCSWM_GO = 1 - SGMAX This does implicitly assume that KRNSW_OW = PCSWM_OW and that KRNSW_GO = PCSWM_GO and will likely need refinement later. On the other hand, the current relations are sufficient for the Norne model. --- src/opm/output/eclipse/LoadRestart.cpp | 123 +++++++++++++++++++++---- src/opm/output/eclipse/RestartIO.cpp | 68 ++++++++++++++ 2 files changed, 172 insertions(+), 19 deletions(-) diff --git a/src/opm/output/eclipse/LoadRestart.cpp b/src/opm/output/eclipse/LoadRestart.cpp index 24f9cd7b8..92da30c2d 100644 --- a/src/opm/output/eclipse/LoadRestart.cpp +++ b/src/opm/output/eclipse/LoadRestart.cpp @@ -46,6 +46,7 @@ #include #include #include +#include #include #include #include @@ -475,6 +476,102 @@ namespace { } } + void insertSolutionVector(const std::vector& vector, + const Opm::RestartKey& value, + const std::vector::size_type numcells, + Opm::data::Solution& sol) + { + if (vector.size() != numcells) { + throw std::runtime_error { + "Restart file: Could not restore '" + + value.key + + "', mismatched number of cells" + }; + } + + sol.insert(value.key, value.dim, vector, + Opm::data::TargetType::RESTART_SOLUTION); + } + + void loadIfAvailable(const RestartFileView& rst_view, + const Opm::RestartKey& value, + const std::vector::size_type numcells, + Opm::data::Solution& sol) + { + const auto* kw = rst_view.getKeyword(value.key.c_str()); + + if (kw == nullptr) { + throwIfMissingRequired(value); + + // If we get here, the requested value was not available in the + // result set. However, the client does not actually require + // the value for restart purposes so we can safely skip this. + return; + } + + insertSolutionVector(double_vector(kw), value, numcells, sol); + } + + void loadHysteresisIfAvailable(const RestartFileView& rst_view, + const std::string& primary, + const Opm::RestartKey& fallback_key, + const std::vector::size_type numcells, + Opm::data::Solution& sol) + { + const auto* kw = rst_view.getKeyword(primary.c_str()); + + if (kw == nullptr) { + // Primary key does not exist in rst_view. Attempt to load + // fallback keys directly. + + loadIfAvailable(rst_view, fallback_key, numcells, sol); + } + else { + // Primary exists in rst_view. Translate to Flow's hysteresis + // parameter. + auto smax = double_vector(kw); + + std::transform(std::begin(smax), std::end(smax), std::begin(smax), + [](const double s) { return 1.0 - s; }); + + insertSolutionVector(smax, fallback_key, numcells, sol); + } + } + + bool isHysteresis(const std::string& vector) + { + for (const auto* flow_hyst_key : { "KRNSW_OW", "PCSWM_OW", + "KRNSW_GO", "PCSWM_GO", }) + { + if (vector == flow_hyst_key) { return true; } + } + + return false; + } + + void restoreHysteresisVector(const Opm::RestartKey& value, + const RestartFileView& rst_view, + const int numcells, + Opm::data::Solution& sol) + { + const auto& key = value.key; + + if ((key == "KRNSW_OW") || (key == "PCSWM_OW")) + { + // Attempt to load from SOMAX, fall back to value.key if + // unavailable--typically in OPM Extended restart file. + loadHysteresisIfAvailable(rst_view, "SOMAX", + value, numcells, sol); + } + else if ((key == "KRNSW_GO") || (key == "PCSWM_GO")) + { + // Attempt to load from SGMAX, fall back to value.key if + // unavailable--typically in OPM Extended restart file. + loadHysteresisIfAvailable(rst_view, "SGMAX", + value, numcells, sol); + } + } + std::vector getOpmExtraFromDoubHEAD(const RestartFileView& rst_view, const bool required, @@ -502,28 +599,16 @@ namespace { Opm::data::Solution sol(/* init_si = */ false); for (const auto& value : solution_keys) { - const auto& vector = value.key; - const auto* kw = rst_view.getKeyword(vector.c_str()); - - if (kw == nullptr) { - throwIfMissingRequired(value); - - // Requested vector not available, but caller does not - // actually require the vector for restart purposes. - // Skip this. + if (isHysteresis(value.key)) { + // Special case handling of hysteresis data. Possibly needs + // translation from ECLIPSE-compatible set to Flow's known + // set of hysteresis vectors. + restoreHysteresisVector(value, rst_view, numcells, sol); continue; } - if (Opm::RestartIO::ecl_kw_get_size(kw) != numcells) { - throw std::runtime_error { - "Restart file: Could not restore " - + std::string(Opm::RestartIO::ecl_kw_get_header(kw)) - + ", mismatched number of cells" - }; - } - - sol.insert(vector, value.dim, double_vector(kw), - Opm::data::TargetType::RESTART_SOLUTION); + // Load regular (non-hysteresis) vector if available. + loadIfAvailable(rst_view, value, numcells, sol); } return sol; diff --git a/src/opm/output/eclipse/RestartIO.cpp b/src/opm/output/eclipse/RestartIO.cpp index 401e2dbae..530190ad6 100644 --- a/src/opm/output/eclipse/RestartIO.cpp +++ b/src/opm/output/eclipse/RestartIO.cpp @@ -42,6 +42,7 @@ #include #include +#include #include #include #include @@ -388,6 +389,69 @@ namespace { write_kw(rst_file, "XCON", connectionData.getXConn()); } + bool haveHysteresis(const RestartValue& value) + { + for (const auto* key : { "KRNSW_OW", "PCSWM_OW", + "KRNSW_GO", "PCSWM_GO", }) + { + if (value.solution.has(key)) { return true; } + } + + return false; + } + + std::vector + convertedHysteresisSat(const RestartValue& value, + const std::string& primary, + const std::string& fallback) + { + auto smax = std::vector{}; + + if (value.solution.has(primary)) { + smax = value.solution.data(primary); + } + else if (value.solution.has(fallback)) { + smax = value.solution.data(fallback); + } + + if (! smax.empty()) { + std::transform(std::begin(smax), std::end(smax), std::begin(smax), + [](const double s) { return 1.0 - s; }); + } + + return smax; + } + + template + void writeEclipseCompatHysteresis(const RestartValue& value, + const bool write_double, + OutputVector&& writeVector) + { + // Convert Flow-specific vectors {KRNSW,PCSWM}_OW to ECLIPSE's + // requisite SOMAX vector. Only partially characterised. + // Sufficient for Norne. + { + const auto somax = + convertedHysteresisSat(value, "KRNSW_OW", "PCSWM_OW"); + + if (! somax.empty()) { + writeVector("SOMAX", somax, write_double); + } + } + + // Convert Flow-specific vectors {KRNSW,PCSWM}_GO to ECLIPSE's + // requisite SGMAX vector. Only partially characterised. + // Sufficient for Norne. + { + const auto sgmax = + convertedHysteresisSat(value, "KRNSW_GO", "PCSWM_GO"); + + if (! sgmax.empty()) { + writeVector("SGMAX", sgmax, write_double); + } + } + } + void writeSolution(ecl_rst_file_type* rst_file, const RestartValue& value, const bool ecl_compatible_rst, @@ -420,6 +484,10 @@ namespace { } } + if (ecl_compatible_rst && haveHysteresis(value)) { + writeEclipseCompatHysteresis(value, write_double_arg, write); + } + ecl_rst_file_end_solution(rst_file); if (ecl_compatible_rst) { From e685eccdc902c67361f0d497afcfa77739d369f7 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Tue, 22 Jan 2019 19:00:46 +0100 Subject: [PATCH 15/20] Restart: Save/Restore Historical Cumulative Production This commit adds support for saving and restoring cumulative production quantities like WOPTH, GWPTH, FGPTH, GGITH, and WWITH. While here, also cater to the case of wells alternating between injecting gas and injecting water. This means that we'll save and restore cumulative gas/water injection for all injectors, irrespective of current injecting phase. Update unit tests accordingly. --- opm/output/eclipse/AggregateGroupData.hpp | 18 +++- opm/output/eclipse/VectorItems/group.hpp | 11 ++ opm/output/eclipse/VectorItems/well.hpp | 32 ++++-- src/opm/output/eclipse/AggregateGroupData.cpp | 7 +- src/opm/output/eclipse/AggregateWellData.cpp | 101 ++++++++++-------- src/opm/output/eclipse/LoadRestart.cpp | 13 +++ src/opm/output/eclipse/Summary.cpp | 2 + tests/test_AggregateWellData.cpp | 64 ++++++++--- tests/test_Restart.cpp | 79 +++++++++++++- tests/test_Summary.cpp | 56 ++++++++++ 10 files changed, 310 insertions(+), 73 deletions(-) diff --git a/opm/output/eclipse/AggregateGroupData.hpp b/opm/output/eclipse/AggregateGroupData.hpp index bfaa61677..3ac4620fe 100644 --- a/opm/output/eclipse/AggregateGroupData.hpp +++ b/opm/output/eclipse/AggregateGroupData.hpp @@ -89,12 +89,16 @@ namespace Opm { namespace RestartIO { namespace Helpers { const std::vector restart_group_keys = {"GOPP", "GWPP", "GOPR", "GWPR", "GGPR", "GVPR", "GWIR", "GGIR", "GWCT", "GGOR", "GOPT", "GWPT", "GGPT", "GVPT", "GWIT", - "GGIT"}; + "GGIT", + "GOPTH", "GWPTH", "GGPTH", + "GWITH", "GGITH"}; const std::vector restart_field_keys = {"FOPP", "FWPP", "FOPR", "FWPR", "FGPR", "FVPR", "FWIR", "FGIR", "FWCT", "FGOR", "FOPT", "FWPT", "FGPT", "FVPT", "FWIT", - "FGIT"}; + "FGIT", + "FOPTH", "FWPTH", "FGPTH", + "FWITH", "FGITH"}; const std::map groupKeyToIndex = { {"GOPR", 0}, @@ -113,6 +117,11 @@ namespace Opm { namespace RestartIO { namespace Helpers { {"GGIT", 16}, {"GOPP", 22}, {"GWPP", 23}, + {"GOPTH", 135}, + {"GWPTH", 139}, + {"GWITH", 140}, + {"GGPTH", 143}, + {"GGITH", 144}, }; const std::map fieldKeyToIndex = { @@ -132,6 +141,11 @@ namespace Opm { namespace RestartIO { namespace Helpers { {"FGIT", 16}, {"FOPP", 22}, {"FWPP", 23}, + {"FOPTH", 135}, + {"FWPTH", 139}, + {"FWITH", 140}, + {"FGPTH", 143}, + {"FGITH", 144}, }; private: diff --git a/opm/output/eclipse/VectorItems/group.hpp b/opm/output/eclipse/VectorItems/group.hpp index 8c714fc06..46853a081 100644 --- a/opm/output/eclipse/VectorItems/group.hpp +++ b/opm/output/eclipse/VectorItems/group.hpp @@ -48,6 +48,17 @@ namespace Opm { namespace RestartIO { namespace Helpers { namespace VectorItems OilPrPot = 22, // Group's oil production potential WatPrPot = 23, // Group's water production potential + + HistOilPrTotal = 135, // Group's total cumulative oil + // production (observed/historical rates) + HistWatPrTotal = 139, // Group's total cumulative water + // production (observed/historical rates) + HistWatInjTotal = 140, // Group's total cumulative water + // injection (observed/historical rates) + HistGasPrTotal = 143, // Group's total cumulative gas + // production (observed/historical rates) + HistGasInjTotal = 144, // Group's total cumulative gas injection + // (observed/historical rates) }; } // XGroup diff --git a/opm/output/eclipse/VectorItems/well.hpp b/opm/output/eclipse/VectorItems/well.hpp index 1203a20c4..bbde343f3 100644 --- a/opm/output/eclipse/VectorItems/well.hpp +++ b/opm/output/eclipse/VectorItems/well.hpp @@ -116,10 +116,15 @@ namespace Opm { namespace RestartIO { namespace Helpers { namespace VectorItems BHPTarget = 6, // Well's bottom hole pressure target DatumDepth = 9, // Well's reference depth for BHP - LiqRateTarget_2 = 33, //Well's liquid rate target/limit for a well on WCONINJH control or for a producer - GasRateTarget_2 = 54, //Well's gas rate target/limit for a well on WCONINJH control or for producer - BHPTarget_2 = 55, //Well's bottom hole pressure target/limit - + + HistLiqRateTarget = 33, // Well's historical/observed liquid + // rate target/limit + + HistGasRateTarget = 54, // Well's historical/observed gas rate + // target/limit + + HistBHPTarget = 55, // Well's historical/observed bottom + // hole pressure target/limit }; } // SWell @@ -146,13 +151,22 @@ namespace Opm { namespace RestartIO { namespace Helpers { namespace VectorItems GasFVF = 34, // Well's producing gas formation volume factor. - item37 = 36, // Unknown - item38 = 37, // Unknown + item37 = 36, // Unknown + item38 = 37, // Unknown - BHPTarget = 41, // Well's current BHP Target/Limit + BHPTarget = 41, // Well's current BHP Target/Limit - item82 = 81, // Unknown - item83 = 82, // Unknown + HistOilPrTotal = 75, // Well's total cumulative oil production + // (observed/historical rates) + HistWatPrTotal = 76, // Well's total cumulative water + // production (observed/historical rates) + HistGasPrTotal = 77, // Well's total cumulative gas production + // (observed(historical rates) + + HistWatInjTotal = 81, // Well's total cumulative water injection + // (observed/historical rates) + HistGasInjTotal = 82, // Well's total cumulative gas injection + // (observed/historical rates) WatVoidPrRate = 122, // Well's voidage production rate GasVoidPrRate = 123, // Well's voidage production rate diff --git a/src/opm/output/eclipse/AggregateGroupData.cpp b/src/opm/output/eclipse/AggregateGroupData.cpp index eb20e4e01..f55bc5d8d 100644 --- a/src/opm/output/eclipse/AggregateGroupData.cpp +++ b/src/opm/output/eclipse/AggregateGroupData.cpp @@ -381,7 +381,8 @@ namespace { for (const auto& key : keys) { if ((key[3] == 'T') && ((key[2] == 'I') || (key[2] == 'P'))) { - // Don't write cumulative quantities in case of + // Don't write cumulative quantities in case of OPM + // extended restart files. continue; } @@ -424,7 +425,7 @@ namespace { }*/ } } - } + } // XGrp namespace ZGrp { std::size_t entriesPerGroup(const std::vector& inteHead) @@ -447,7 +448,7 @@ namespace { WV::WindowSize{ entriesPerGroup(inteHead) } }; } - } + } // ZGrp } // Anonymous void diff --git a/src/opm/output/eclipse/AggregateWellData.cpp b/src/opm/output/eclipse/AggregateWellData.cpp index ca0aee54e..d2339afe4 100644 --- a/src/opm/output/eclipse/AggregateWellData.cpp +++ b/src/opm/output/eclipse/AggregateWellData.cpp @@ -30,7 +30,6 @@ #include #include #include -//#include #include #include @@ -511,13 +510,13 @@ namespace { if ((pp.GasRate != 0.0) || (!predMode)) { sWell[Ix::GasRateTarget] = swprop(M::gas_surface_rate, pp.GasRate); - sWell[Ix::GasRateTarget_2] = sWell[Ix::GasRateTarget]; + sWell[Ix::HistGasRateTarget] = sWell[Ix::GasRateTarget]; } if (pp.LiquidRate != 0.0 || (!predMode)) { sWell[Ix::LiqRateTarget] = swprop(M::liquid_surface_rate, pp.LiquidRate); - sWell[Ix::LiqRateTarget_2] = sWell[Ix::LiqRateTarget]; + sWell[Ix::HistLiqRateTarget] = sWell[Ix::LiqRateTarget]; } else { sWell[Ix::LiqRateTarget] = @@ -542,26 +541,29 @@ namespace { sWell[Ix::BHPTarget] = pp.BHPLimit != 0.0 ? swprop(M::pressure, pp.BHPLimit) : swprop(M::pressure, 1.0*::Opm::unit::atm); - sWell[Ix::BHPTarget_2] = sWell[Ix::BHPTarget]; + sWell[Ix::HistBHPTarget] = sWell[Ix::BHPTarget]; } else if (well.isInjector(sim_step)) { const auto& ip = well.getInjectionProperties(sim_step); using IP = ::Opm::WellInjector::ControlModeEnum; - using IT = ::Opm::WellInjector::TypeEnum; + using IT = ::Opm::WellInjector::TypeEnum; - if (ip.hasInjectionControl(IP::RATE)) { - if (ip.injectorType == IT::OIL) { - sWell[Ix::OilRateTarget] = swprop(M::liquid_surface_rate, ip.surfaceInjectionRate); - } - if (ip.injectorType == IT::WATER) { - sWell[Ix::WatRateTarget] = swprop(M::liquid_surface_rate, ip.surfaceInjectionRate); - sWell[Ix::LiqRateTarget_2] = sWell[Ix::WatRateTarget]; - } - if (ip.injectorType == IT::GAS) { - sWell[Ix::GasRateTarget] = swprop(M::gas_surface_rate, ip.surfaceInjectionRate); - sWell[Ix::GasRateTarget_2] = sWell[Ix::GasRateTarget]; - } + if (ip.hasInjectionControl(IP::RATE)) { + if (ip.injectorType == IT::OIL) { + sWell[Ix::OilRateTarget] = + swprop(M::liquid_surface_rate, ip.surfaceInjectionRate); + } + if (ip.injectorType == IT::WATER) { + sWell[Ix::WatRateTarget] = + swprop(M::liquid_surface_rate, ip.surfaceInjectionRate); + sWell[Ix::HistLiqRateTarget] = sWell[Ix::WatRateTarget]; + } + if (ip.injectorType == IT::GAS) { + sWell[Ix::GasRateTarget] = + swprop(M::gas_surface_rate, ip.surfaceInjectionRate); + sWell[Ix::HistGasRateTarget] = sWell[Ix::GasRateTarget]; + } } if (ip.hasInjectionControl(IP::RESV)) { @@ -575,7 +577,7 @@ namespace { sWell[Ix::BHPTarget] = ip.hasInjectionControl(IP::BHP) ? swprop(M::pressure, ip.BHPLimit) : swprop(M::pressure, 1.0E05*::Opm::unit::psia); - sWell[Ix::BHPTarget_2] = sWell[Ix::BHPTarget]; + sWell[Ix::HistBHPTarget] = sWell[Ix::BHPTarget]; } sWell[Ix::DatumDepth] = @@ -647,15 +649,41 @@ namespace { xWell[Ix::WatCut] = get("WWCT"); xWell[Ix::GORatio] = get("WGOR"); - if (ecl_compatible_rst) { - xWell[Ix::OilPrTotal] = get("WOPT"); - xWell[Ix::WatPrTotal] = get("WWPT"); - xWell[Ix::GasPrTotal] = get("WGPT"); - xWell[Ix::VoidPrTotal] = get("WVPT"); - } + if (ecl_compatible_rst) { + xWell[Ix::OilPrTotal] = get("WOPT"); + xWell[Ix::WatPrTotal] = get("WWPT"); + xWell[Ix::GasPrTotal] = get("WGPT"); + xWell[Ix::VoidPrTotal] = get("WVPT"); + } + // Not fully characterised. xWell[Ix::item37] = xWell[Ix::WatPrRate]; xWell[Ix::item38] = xWell[Ix::GasPrRate]; + + if (ecl_compatible_rst) { + xWell[Ix::HistOilPrTotal] = get("WOPTH"); + xWell[Ix::HistWatPrTotal] = get("WWPTH"); + xWell[Ix::HistGasPrTotal] = get("WGPTH"); + } + } + + template + void assignCommonInjector(GetSummaryVector& get, + const bool ecl_compatible_rst, + XWellArray& xWell) + { + using Ix = ::Opm::RestartIO::Helpers::VectorItems::XWell::index; + + xWell[Ix::FlowBHP] = get("WBHP"); + + if (ecl_compatible_rst) { + // Note: Assign both water and gas cumulatives to support + // case of well alternating between injecting water and gas. + xWell[Ix::WatInjTotal] = get("WWIT"); + xWell[Ix::GasInjTotal] = get("WGIT"); + xWell[Ix::HistWatInjTotal] = get("WWITH"); + xWell[Ix::HistGasInjTotal] = get("WGITH"); + } } template @@ -673,19 +701,14 @@ namespace { return smry.has(key) ? smry.get(key) : 0.0; }; - // Injection rates reported as negative, cumulative - // totals as positive. + assignCommonInjector(get, ecl_compatible_rst, xWell); + + // Injection rates reported as negative. xWell[Ix::WatPrRate] = -get("WWIR"); xWell[Ix::LiqPrRate] = xWell[Ix::WatPrRate]; - xWell[Ix::FlowBHP] = get("WBHP"); - - if (ecl_compatible_rst) { - xWell[Ix::WatInjTotal] = get("WWIT"); - } - + // Not fully characterised. xWell[Ix::item37] = xWell[Ix::WatPrRate]; - xWell[Ix::item82] = xWell[Ix::WatInjTotal]; xWell[Ix::WatVoidPrRate] = -get("WWVIR"); } @@ -705,17 +728,12 @@ namespace { return smry.has(key) ? smry.get(key) : 0.0; }; - // Injection rates reported as negative production rates, - // cumulative injection totals as positive. + assignCommonInjector(get, ecl_compatible_rst, xWell); + + // Injection rates reported as negative production rates. xWell[Ix::GasPrRate] = -get("WGIR"); xWell[Ix::VoidPrRate] = -get("WGVIR"); - xWell[Ix::FlowBHP] = get("WBHP"); - - if (ecl_compatible_rst) { - xWell[Ix::GasInjTotal] = get("WGIT"); - } - xWell[Ix::GasFVF] = (std::abs(xWell[Ix::GasPrRate]) > 0.0) ? xWell[Ix::VoidPrRate] / xWell[Ix::GasPrRate] : 0.0; @@ -724,7 +742,6 @@ namespace { // Not fully characterised. xWell[Ix::item38] = xWell[Ix::GasPrRate]; - xWell[Ix::item83] = xWell[Ix::GasInjTotal]; xWell[Ix::GasVoidPrRate] = xWell[Ix::VoidPrRate]; } diff --git a/src/opm/output/eclipse/LoadRestart.cpp b/src/opm/output/eclipse/LoadRestart.cpp index 92da30c2d..07e5d35ae 100644 --- a/src/opm/output/eclipse/LoadRestart.cpp +++ b/src/opm/output/eclipse/LoadRestart.cpp @@ -1088,6 +1088,13 @@ namespace { smry.add(key("WWIT"), xwel[VI::XWell::index::WatInjTotal]); smry.add(key("WGIT"), xwel[VI::XWell::index::GasInjTotal]); + + smry.add(key("WOPTH"), xwel[VI::XWell::index::HistOilPrTotal]); + smry.add(key("WWPTH"), xwel[VI::XWell::index::HistWatPrTotal]); + smry.add(key("WGPTH"), xwel[VI::XWell::index::HistGasPrTotal]); + + smry.add(key("WWITH"), xwel[VI::XWell::index::HistWatInjTotal]); + smry.add(key("WGITH"), xwel[VI::XWell::index::HistGasInjTotal]); } void assign_group_cumulatives(const std::string& group, @@ -1118,6 +1125,12 @@ namespace { smry.add(key("WIT"), xgrp[VI::XGroup::index::WatInjTotal]); smry.add(key("GIT"), xgrp[VI::XGroup::index::GasInjTotal]); + + smry.add(key("OPTH"), xgrp[VI::XGroup::index::HistOilPrTotal]); + smry.add(key("WPTH"), xgrp[VI::XGroup::index::HistWatPrTotal]); + smry.add(key("GPTH"), xgrp[VI::XGroup::index::HistGasPrTotal]); + smry.add(key("WITH"), xgrp[VI::XGroup::index::HistWatInjTotal]); + smry.add(key("GITH"), xgrp[VI::XGroup::index::HistGasInjTotal]); } Opm::SummaryState diff --git a/src/opm/output/eclipse/Summary.cpp b/src/opm/output/eclipse/Summary.cpp index 0b70db34f..1ef8d5abe 100644 --- a/src/opm/output/eclipse/Summary.cpp +++ b/src/opm/output/eclipse/Summary.cpp @@ -65,6 +65,8 @@ namespace { "WIR", "GIR", "WIT", "GIT", "WCT", "GOR", + "OPTH", "WPTH", "GPTH", + "WITH", "GITH", }; } diff --git a/tests/test_AggregateWellData.cpp b/tests/test_AggregateWellData.cpp index 51c1062fd..3bc0d498a 100644 --- a/tests/test_AggregateWellData.cpp +++ b/tests/test_AggregateWellData.cpp @@ -1,4 +1,5 @@ /* + Copyright 2019 Equinor Copyright 2018 Statoil ASA This file is part of the Open Porous Media project (OPM). @@ -232,6 +233,11 @@ TSTEP -- 8 state.add("WWCT:OP_1" , 0.625); state.add("WGOR:OP_1" , 234.5); state.add("WBHP:OP_1" , 314.15); + state.add("WOPTH:OP_1", 345.6); + state.add("WWPTH:OP_1", 456.7); + state.add("WGPTH:OP_1", 567.8); + state.add("WWITH:OP_1", 0.0); + state.add("WGITH:OP_1", 0.0); state.add("WGVIR:OP_1", 0.0); state.add("WWVIR:OP_1", 0.0); @@ -250,6 +256,11 @@ TSTEP -- 8 state.add("WWCT:OP_2" , 0.0); state.add("WGOR:OP_2" , 0.0); state.add("WBHP:OP_2" , 400.6); + state.add("WOPTH:OP_2", 0.0); + state.add("WWPTH:OP_2", 0.0); + state.add("WGPTH:OP_2", 0.0); + state.add("WWITH:OP_2", 1515.0); + state.add("WGITH:OP_2", 3030.0); state.add("WGVIR:OP_2", 1234.0); state.add("WWVIR:OP_2", 4321.0); @@ -268,6 +279,11 @@ TSTEP -- 8 state.add("WWCT:OP_3" , 0.0625); state.add("WGOR:OP_3" , 1234.5); state.add("WBHP:OP_3" , 314.15); + state.add("WOPTH:OP_3", 2345.6); + state.add("WWPTH:OP_3", 3456.7); + state.add("WGPTH:OP_3", 4567.8); + state.add("WWITH:OP_3", 0.0); + state.add("WGITH:OP_3", 0.0); state.add("WGVIR:OP_3", 0.0); state.add("WWVIR:OP_3", 43.21); @@ -561,17 +577,27 @@ BOOST_AUTO_TEST_CASE (Dynamic_Well_Data_Step1) BOOST_CHECK_CLOSE(xwell[i0 + Ix::LiqPrRate], 1.0 + 2.0, 1.0e-10); BOOST_CHECK_CLOSE(xwell[i0 + Ix::VoidPrRate], 4.0, 1.0e-10); - BOOST_CHECK_CLOSE(xwell[i0 + Ix::FlowBHP], 314.15, 1.0e-10); - BOOST_CHECK_CLOSE(xwell[i0 + Ix::WatCut] , 0.625, 1.0e-10); - BOOST_CHECK_CLOSE(xwell[i0 + Ix::GORatio], 234.5, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i0 + Ix::FlowBHP], 314.15 , 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i0 + Ix::WatCut] , 0.625, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i0 + Ix::GORatio], 234.5 , 1.0e-10); BOOST_CHECK_CLOSE(xwell[i0 + Ix::OilPrTotal], 10.0, 1.0e-10); BOOST_CHECK_CLOSE(xwell[i0 + Ix::WatPrTotal], 20.0, 1.0e-10); BOOST_CHECK_CLOSE(xwell[i0 + Ix::GasPrTotal], 30.0, 1.0e-10); BOOST_CHECK_CLOSE(xwell[i0 + Ix::VoidPrTotal], 40.0, 1.0e-10); - BOOST_CHECK_CLOSE(xwell[i0 + Ix::item37], xwell[i0 + Ix::WatPrRate], 1.0e-10); - BOOST_CHECK_CLOSE(xwell[i0 + Ix::item38], xwell[i0 + Ix::GasPrRate], 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i0 + Ix::item37], + xwell[i0 + Ix::WatPrRate], 1.0e-10); + + BOOST_CHECK_CLOSE(xwell[i0 + Ix::item38], + xwell[i0 + Ix::GasPrRate], 1.0e-10); + + BOOST_CHECK_CLOSE(xwell[i0 + Ix::HistOilPrTotal], 345.6, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i0 + Ix::HistWatPrTotal], 456.7, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i0 + Ix::HistGasPrTotal], 567.8, 1.0e-10); + + BOOST_CHECK_CLOSE(xwell[i0 + Ix::HistWatInjTotal], 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i0 + Ix::HistGasInjTotal], 0.0, 1.0e-10); } // XWEL (OP_2) @@ -585,6 +611,7 @@ BOOST_AUTO_TEST_CASE (Dynamic_Well_Data_Step1) BOOST_CHECK_CLOSE(xwell[i1 + Ix::VoidPrRate], -1234.0, 1.0e-10); BOOST_CHECK_CLOSE(xwell[i1 + Ix::FlowBHP], 400.6, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i1 + Ix::WatInjTotal], 1000.0, 1.0e-10); BOOST_CHECK_CLOSE(xwell[i1 + Ix::GasInjTotal], 2000.0, 1.0e-10); // Bg = VGIR / GIR = 1234.0 / 200.0 @@ -593,11 +620,11 @@ BOOST_AUTO_TEST_CASE (Dynamic_Well_Data_Step1) BOOST_CHECK_CLOSE(xwell[i1 + Ix::item38], xwell[i1 + Ix::GasPrRate], 1.0e-10); - BOOST_CHECK_CLOSE(xwell[i1 + Ix::item83], - xwell[i1 + Ix::GasInjTotal], 1.0e-10); - - BOOST_CHECK_CLOSE(xwell[i1 + Ix::GasVoidPrRate], - xwell[i1 + Ix::VoidPrRate], 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i1 + Ix::HistOilPrTotal] , 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i1 + Ix::HistWatPrTotal] , 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i1 + Ix::HistGasPrTotal] , 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i1 + Ix::HistWatInjTotal], 1515.0, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i1 + Ix::HistGasInjTotal], 3030.0, 1.0e-10); } } @@ -672,6 +699,10 @@ BOOST_AUTO_TEST_CASE (Dynamic_Well_Data_Step2) BOOST_CHECK_CLOSE(xwell[i0 + Ix::item38], xwell[i0 + Ix::GasPrRate], 1.0e-10); + + BOOST_CHECK_CLOSE(xwell[i0 + Ix::HistOilPrTotal], 345.6, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i0 + Ix::HistWatPrTotal], 456.7, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i0 + Ix::HistGasPrTotal], 567.8, 1.0e-10); } // XWEL (OP_2) -- water injector @@ -690,14 +721,17 @@ BOOST_AUTO_TEST_CASE (Dynamic_Well_Data_Step2) BOOST_CHECK_CLOSE(xwell[i1 + Ix::FlowBHP], 400.6, 1.0e-10); BOOST_CHECK_CLOSE(xwell[i1 + Ix::WatInjTotal], 1000.0, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i1 + Ix::GasInjTotal], 2000.0, 1.0e-10); // Copy of WWIR BOOST_CHECK_CLOSE(xwell[i1 + Ix::item37], xwell[i1 + Ix::WatPrRate], 1.0e-10); - // Copy of WWIT - BOOST_CHECK_CLOSE(xwell[i1 + Ix::item82], - xwell[i1 + Ix::WatInjTotal], 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i1 + Ix::HistOilPrTotal] , 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i1 + Ix::HistWatPrTotal] , 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i1 + Ix::HistGasPrTotal] , 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i1 + Ix::HistWatInjTotal], 1515.0, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i1 + Ix::HistGasInjTotal], 3030.0, 1.0e-10); // WWVIR BOOST_CHECK_CLOSE(xwell[i1 + Ix::WatVoidPrRate], @@ -733,6 +767,10 @@ BOOST_AUTO_TEST_CASE (Dynamic_Well_Data_Step2) // Copy of WGPR BOOST_CHECK_CLOSE(xwell[i2 + Ix::item38], xwell[i2 + Ix::GasPrRate], 1.0e-10); + + BOOST_CHECK_CLOSE(xwell[i2 + Ix::HistOilPrTotal], 2345.6, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i2 + Ix::HistWatPrTotal], 3456.7, 1.0e-10); + BOOST_CHECK_CLOSE(xwell[i2 + Ix::HistGasPrTotal], 4567.8, 1.0e-10); } } diff --git a/tests/test_Restart.cpp b/tests/test_Restart.cpp index 0fb26a47d..8c3edaedb 100644 --- a/tests/test_Restart.cpp +++ b/tests/test_Restart.cpp @@ -366,6 +366,11 @@ Opm::SummaryState sim_state() state.add("WWCT:OP_1" , 0.625); state.add("WGOR:OP_1" , 234.5); state.add("WBHP:OP_1" , 314.15); + state.add("WOPTH:OP_1", 345.6); + state.add("WWPTH:OP_1", 456.7); + state.add("WGPTH:OP_1", 567.8); + state.add("WWITH:OP_1", 0.0); + state.add("WGITH:OP_1", 0.0); state.add("WGVIR:OP_1", 0.0); state.add("WWVIR:OP_1", 0.0); @@ -384,6 +389,11 @@ Opm::SummaryState sim_state() state.add("WWCT:OP_2" , 0.0); state.add("WGOR:OP_2" , 0.0); state.add("WBHP:OP_2" , 400.6); + state.add("WOPTH:OP_2", 0.0); + state.add("WWPTH:OP_2", 0.0); + state.add("WGPTH:OP_2", 0.0); + state.add("WWITH:OP_2", 1515.0); + state.add("WGITH:OP_2", 3030.0); state.add("WGVIR:OP_2", 1234.0); state.add("WWVIR:OP_2", 4321.0); @@ -402,6 +412,11 @@ Opm::SummaryState sim_state() state.add("WWCT:OP_3" , 0.0625); state.add("WGOR:OP_3" , 1234.5); state.add("WBHP:OP_3" , 314.15); + state.add("WOPTH:OP_3", 2345.6); + state.add("WWPTH:OP_3", 3456.7); + state.add("WGPTH:OP_3", 4567.8); + state.add("WWITH:OP_3", 0.0); + state.add("WGITH:OP_3", 0.0); state.add("WGVIR:OP_3", 0.0); state.add("WWVIR:OP_3", 43.21); @@ -419,6 +434,13 @@ Opm::SummaryState sim_state() state.add("GGIT:OP" , 27182.8); state.add("GWCT:OP" , 0.625); state.add("GGOR:OP" , 1234.5); + state.add("GGVIR:OP", 123.45); + state.add("GWVIR:OP", 1234.56); + state.add("GOPTH:OP", 5678.90); + state.add("GWPTH:OP", 6789.01); + state.add("GGPTH:OP", 7890.12); + state.add("GWITH:OP", 8901.23); + state.add("GGITH:OP", 9012.34); state.add("FOPR" , 1100.0); state.add("FWPR" , 1200.0); @@ -433,7 +455,14 @@ Opm::SummaryState sim_state() state.add("FWIT" , 314159.2); state.add("FGIT" , 271828.1); state.add("FWCT" , 0.625); - state.add("FGOR" , 1234.5); + state.add("FGOR" , 1234.5); + state.add("FOPTH", 56789.01); + state.add("FWPTH", 67890.12); + state.add("FGPTH", 78901.23); + state.add("FWITH", 89012.34); + state.add("FGITH", 90123.45); + state.add("FGVIR", 1234.56); + state.add("FWVIR", 12345.67); return state; } @@ -835,13 +864,18 @@ BOOST_AUTO_TEST_CASE(Restore_Cumulatives) // Verify that the restored summary state has all of its requisite // cumulative summary vectors. - // Producer => W*IT saved/restored as zero (0.0) + // Producer => W*IT{,H} saved/restored as zero (0.0) BOOST_CHECK(rstSumState.has("WOPT:OP_1")); BOOST_CHECK(rstSumState.has("WGPT:OP_1")); BOOST_CHECK(rstSumState.has("WWPT:OP_1")); BOOST_CHECK(rstSumState.has("WVPT:OP_1")); BOOST_CHECK(rstSumState.has("WWIT:OP_1")); BOOST_CHECK(rstSumState.has("WGIT:OP_1")); + BOOST_CHECK(rstSumState.has("WOPTH:OP_1")); + BOOST_CHECK(rstSumState.has("WGPTH:OP_1")); + BOOST_CHECK(rstSumState.has("WWPTH:OP_1")); + BOOST_CHECK(rstSumState.has("WWITH:OP_1")); + BOOST_CHECK(rstSumState.has("WGITH:OP_1")); BOOST_CHECK_CLOSE(rstSumState.get("WOPT:OP_1"), 10.0, 1.0e-10); BOOST_CHECK_CLOSE(rstSumState.get("WGPT:OP_1"), 30.0, 1.0e-10); @@ -849,21 +883,36 @@ BOOST_AUTO_TEST_CASE(Restore_Cumulatives) BOOST_CHECK_CLOSE(rstSumState.get("WVPT:OP_1"), 40.0, 1.0e-10); BOOST_CHECK_CLOSE(rstSumState.get("WWIT:OP_1"), 0.0, 1.0e-10); BOOST_CHECK_CLOSE(rstSumState.get("WGIT:OP_1"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WOPTH:OP_1"), 345.6, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WWPTH:OP_1"), 456.7, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WGPTH:OP_1"), 567.8, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WWITH:OP_1"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WGITH:OP_1"), 0.0, 1.0e-10); - // Gas injector => W*PT and WWIT saved/restored as zero (0.0) + // Gas injector => W*PT{,H} saved/restored as zero (0.0) BOOST_CHECK(rstSumState.has("WOPT:OP_2")); BOOST_CHECK(rstSumState.has("WGPT:OP_2")); BOOST_CHECK(rstSumState.has("WWPT:OP_2")); BOOST_CHECK(rstSumState.has("WVPT:OP_2")); BOOST_CHECK(rstSumState.has("WWIT:OP_2")); BOOST_CHECK(rstSumState.has("WGIT:OP_2")); + BOOST_CHECK(rstSumState.has("WOPTH:OP_2")); + BOOST_CHECK(rstSumState.has("WGPTH:OP_2")); + BOOST_CHECK(rstSumState.has("WWPTH:OP_2")); + BOOST_CHECK(rstSumState.has("WWITH:OP_2")); + BOOST_CHECK(rstSumState.has("WGITH:OP_2")); BOOST_CHECK_CLOSE(rstSumState.get("WOPT:OP_2"), 0.0, 1.0e-10); BOOST_CHECK_CLOSE(rstSumState.get("WGPT:OP_2"), 0.0, 1.0e-10); BOOST_CHECK_CLOSE(rstSumState.get("WWPT:OP_2"), 0.0, 1.0e-10); BOOST_CHECK_CLOSE(rstSumState.get("WVPT:OP_2"), 0.0, 1.0e-10); - BOOST_CHECK_CLOSE(rstSumState.get("WWIT:OP_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WWIT:OP_2"), 1000.0, 1.0e-10); BOOST_CHECK_CLOSE(rstSumState.get("WGIT:OP_2"), 2000.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WOPTH:OP_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WGPTH:OP_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WWPTH:OP_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WWITH:OP_2"), 1515.0, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("WGITH:OP_2"), 3030.0, 1.0e-10); // Group cumulatives saved/restored for all phases BOOST_CHECK(rstSumState.has("GOPT:OP")); @@ -872,6 +921,11 @@ BOOST_AUTO_TEST_CASE(Restore_Cumulatives) BOOST_CHECK(rstSumState.has("GVPT:OP")); BOOST_CHECK(rstSumState.has("GWIT:OP")); BOOST_CHECK(rstSumState.has("GGIT:OP")); + BOOST_CHECK(rstSumState.has("GOPTH:OP")); + BOOST_CHECK(rstSumState.has("GGPTH:OP")); + BOOST_CHECK(rstSumState.has("GWPTH:OP")); + BOOST_CHECK(rstSumState.has("GWITH:OP")); + BOOST_CHECK(rstSumState.has("GGITH:OP")); BOOST_CHECK_CLOSE(rstSumState.get("GOPT:OP"), 1100.0, 1.0e-10); BOOST_CHECK_CLOSE(rstSumState.get("GWPT:OP"), 1200.0, 1.0e-10); @@ -880,6 +934,12 @@ BOOST_AUTO_TEST_CASE(Restore_Cumulatives) BOOST_CHECK_CLOSE(rstSumState.get("GWIT:OP"), 31415.9, 1.0e-10); BOOST_CHECK_CLOSE(rstSumState.get("GGIT:OP"), 27182.8, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("GOPTH:OP"), 5678.90, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("GGPTH:OP"), 7890.12, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("GWPTH:OP"), 6789.01, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("GWITH:OP"), 8901.23, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("GGITH:OP"), 9012.34, 1.0e-10); + // Field cumulatives saved/restored for all phases BOOST_CHECK(rstSumState.has("FOPT")); BOOST_CHECK(rstSumState.has("FGPT")); @@ -887,6 +947,11 @@ BOOST_AUTO_TEST_CASE(Restore_Cumulatives) BOOST_CHECK(rstSumState.has("FVPT")); BOOST_CHECK(rstSumState.has("FWIT")); BOOST_CHECK(rstSumState.has("FGIT")); + BOOST_CHECK(rstSumState.has("FOPTH")); + BOOST_CHECK(rstSumState.has("FGPTH")); + BOOST_CHECK(rstSumState.has("FWPTH")); + BOOST_CHECK(rstSumState.has("FWITH")); + BOOST_CHECK(rstSumState.has("FGITH")); BOOST_CHECK_CLOSE(rstSumState.get("FOPT"), 11000.0, 1.0e-10); BOOST_CHECK_CLOSE(rstSumState.get("FWPT"), 12000.0, 1.0e-10); @@ -894,6 +959,12 @@ BOOST_AUTO_TEST_CASE(Restore_Cumulatives) BOOST_CHECK_CLOSE(rstSumState.get("FVPT"), 14000.0, 1.0e-10); BOOST_CHECK_CLOSE(rstSumState.get("FWIT"), 314159.2, 1.0e-10); BOOST_CHECK_CLOSE(rstSumState.get("FGIT"), 271828.1, 1.0e-10); + + BOOST_CHECK_CLOSE(rstSumState.get("FOPTH"), 56789.01, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("FGPTH"), 78901.23, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("FWPTH"), 67890.12, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("FWITH"), 89012.34, 1.0e-10); + BOOST_CHECK_CLOSE(rstSumState.get("FGITH"), 90123.45, 1.0e-10); } diff --git a/tests/test_Summary.cpp b/tests/test_Summary.cpp index 53ea476b6..cd5ccdd0a 100644 --- a/tests/test_Summary.cpp +++ b/tests/test_Summary.cpp @@ -2921,6 +2921,13 @@ BOOST_AUTO_TEST_CASE(Reset) rstrt.add("WWIT:W_1", 5.0); rstrt.add("WGIT:W_1", 6.0); + rstrt.add("WOPTH:W_1", 0.1); + rstrt.add("WWPTH:W_1", 0.2); + rstrt.add("WGPTH:W_1", 0.3); + + rstrt.add("WWITH:W_1", 0.5); + rstrt.add("WGITH:W_1", 0.6); + rstrt.add("GOPT:NoSuchGroup", 1.0); rstrt.add("GWPT:NoSuchGroup", 2.0); rstrt.add("GGPT:NoSuchGroup", 3.0); @@ -2937,6 +2944,13 @@ BOOST_AUTO_TEST_CASE(Reset) rstrt.add("FWIT", 50.0); rstrt.add("FGIT", 60.0); + rstrt.add("FOPTH", 0.01); + rstrt.add("FWPTH", 0.02); + rstrt.add("FGPTH", 0.03); + + rstrt.add("FWITH", 0.05); + rstrt.add("FGITH", 0.06); + smry.reset_cumulative_quantities(rstrt); const auto& sumstate = smry.get_restart_vectors(); @@ -2969,6 +2983,13 @@ BOOST_AUTO_TEST_CASE(Reset) BOOST_CHECK_CLOSE(sumstate.get("WWIT:W_1"), 5.0, 1.0e-10); BOOST_CHECK_CLOSE(sumstate.get("WGIT:W_1"), 6.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WOPTH:W_1"), 0.1, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WWPTH:W_1"), 0.2, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WGPTH:W_1"), 0.3, 1.0e-10); + + BOOST_CHECK_CLOSE(sumstate.get("WWITH:W_1"), 0.5, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WGITH:W_1"), 0.6, 1.0e-10); + // Cumulatives unset for W_2. BOOST_CHECK_CLOSE(sumstate.get("WOPT:W_2"), 0.0, 1.0e-10); BOOST_CHECK_CLOSE(sumstate.get("WWPT:W_2"), 0.0, 1.0e-10); @@ -2978,6 +2999,13 @@ BOOST_AUTO_TEST_CASE(Reset) BOOST_CHECK_CLOSE(sumstate.get("WWIT:W_2"), 0.0, 1.0e-10); BOOST_CHECK_CLOSE(sumstate.get("WGIT:W_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WOPTH:W_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WWPTH:W_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WGPTH:W_2"), 0.0, 1.0e-10); + + BOOST_CHECK_CLOSE(sumstate.get("WWITH:W_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WGITH:W_2"), 0.0, 1.0e-10); + // Cumulatives unset for W_3. BOOST_CHECK_CLOSE(sumstate.get("WOPT:W_3"), 0.0, 1.0e-10); BOOST_CHECK_CLOSE(sumstate.get("WWPT:W_3"), 0.0, 1.0e-10); @@ -2987,6 +3015,13 @@ BOOST_AUTO_TEST_CASE(Reset) BOOST_CHECK_CLOSE(sumstate.get("WWIT:W_3"), 0.0, 1.0e-10); BOOST_CHECK_CLOSE(sumstate.get("WGIT:W_3"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WOPTH:W_3"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WWPTH:W_3"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WGPTH:W_3"), 0.0, 1.0e-10); + + BOOST_CHECK_CLOSE(sumstate.get("WWITH:W_3"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("WGITH:W_3"), 0.0, 1.0e-10); + // Cumulatives unset for G_1. BOOST_CHECK_CLOSE(sumstate.get("GOPT:G_1"), 0.0, 1.0e-10); BOOST_CHECK_CLOSE(sumstate.get("GWPT:G_1"), 0.0, 1.0e-10); @@ -2996,6 +3031,13 @@ BOOST_AUTO_TEST_CASE(Reset) BOOST_CHECK_CLOSE(sumstate.get("GWIT:G_1"), 0.0, 1.0e-10); BOOST_CHECK_CLOSE(sumstate.get("GGIT:G_1"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GOPTH:G_1"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GWPTH:G_1"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GGPTH:G_1"), 0.0, 1.0e-10); + + BOOST_CHECK_CLOSE(sumstate.get("GWITH:G_1"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GGITH:G_1"), 0.0, 1.0e-10); + // Cumulatives unset for G_2. BOOST_CHECK_CLOSE(sumstate.get("GOPT:G_2"), 0.0, 1.0e-10); BOOST_CHECK_CLOSE(sumstate.get("GWPT:G_2"), 0.0, 1.0e-10); @@ -3005,6 +3047,13 @@ BOOST_AUTO_TEST_CASE(Reset) BOOST_CHECK_CLOSE(sumstate.get("GWIT:G_2"), 0.0, 1.0e-10); BOOST_CHECK_CLOSE(sumstate.get("GGIT:G_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GOPTH:G_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GWPTH:G_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GGPTH:G_2"), 0.0, 1.0e-10); + + BOOST_CHECK_CLOSE(sumstate.get("GWITH:G_2"), 0.0, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("GGITH:G_2"), 0.0, 1.0e-10); + // Cumulatives reset for FIELD. BOOST_CHECK_CLOSE(sumstate.get("FOPT"), 10.0, 1.0e-10); BOOST_CHECK_CLOSE(sumstate.get("FWPT"), 20.0, 1.0e-10); @@ -3013,6 +3062,13 @@ BOOST_AUTO_TEST_CASE(Reset) BOOST_CHECK_CLOSE(sumstate.get("FWIT"), 50.0, 1.0e-10); BOOST_CHECK_CLOSE(sumstate.get("FGIT"), 60.0, 1.0e-10); + + BOOST_CHECK_CLOSE(sumstate.get("FOPTH"), 0.01, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("FWPTH"), 0.02, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("FGPTH"), 0.03, 1.0e-10); + + BOOST_CHECK_CLOSE(sumstate.get("FWITH"), 0.05, 1.0e-10); + BOOST_CHECK_CLOSE(sumstate.get("FGITH"), 0.06, 1.0e-10); } BOOST_AUTO_TEST_SUITE_END() From 849eb80e412305d3cf13b3294efd2200bfa77162 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Tue, 22 Jan 2019 20:11:43 +0100 Subject: [PATCH 16/20] IWEL: Identify Items for Active and Requested Well Controls This commit refines our understanding of IWEL items relating to well constraints. In particular Item 8 is the well's active control mode as determined by the simulator from dynamic state variables. On the other hand, Item 16 is the constraint that is requested in the simulation run input for prediction wells (keywords WCONINJE, WCONPROD) while Item 50 is the requested constraint for wells controlled by observed rates (WCONINJH, WCONHIST). Special Note: This commit outputs the requested control mode to Item 8 and will need an update later once the simulator becomes aware of the distinction. --- opm/output/eclipse/VectorItems/well.hpp | 9 ++++-- src/opm/output/eclipse/AggregateWellData.cpp | 29 ++++++++++++++++++-- tests/test_AggregateWellData.cpp | 4 +-- 3 files changed, 35 insertions(+), 7 deletions(-) diff --git a/opm/output/eclipse/VectorItems/well.hpp b/opm/output/eclipse/VectorItems/well.hpp index bbde343f3..4a7d1564b 100644 --- a/opm/output/eclipse/VectorItems/well.hpp +++ b/opm/output/eclipse/VectorItems/well.hpp @@ -33,19 +33,24 @@ namespace Opm { namespace RestartIO { namespace Helpers { namespace VectorItems NConn = 4, // Number of active cells connected to well Group = 5, // Index (one-based) of well's current group WType = 6, // Well type - WCtrl = 7, // Well control + ActWCtrl = 7, // Well's active target control mode (constraint). item9 = 8, // Unknown item11 = 10, // Unknown VFPTab = 11, // ID (one-based) of well's current VFP table + PredReqWCtrl = 15, // Well's requested control mode from + // simulation deck (WCONINJE, WCONPROD). + item18 = 17, // Unknown XFlow = 22, item25 = 24, // Unknown item32 = 31, // Unknown item48 = 47, // Unknown - item50 = 49, // Unknown + + HistReqWCtrl = 49, // Well's requested control mode from + // simulation deck (WCONHIST, WCONINJH) MsWID = 70, // Multisegment well ID // Value 0 for regular wells diff --git a/src/opm/output/eclipse/AggregateWellData.cpp b/src/opm/output/eclipse/AggregateWellData.cpp index d2339afe4..1e9a4594c 100644 --- a/src/opm/output/eclipse/AggregateWellData.cpp +++ b/src/opm/output/eclipse/AggregateWellData.cpp @@ -337,7 +337,6 @@ namespace { GroupMapNameInd); iWell[Ix::WType] = wellType (well, sim_step); - iWell[Ix::WCtrl] = ctrlMode (well, sim_step); iWell[Ix::VFPTab] = wellVFPTab(well, sim_step); iWell[Ix::XFlow] = well.getAllowCrossFlow() ? 1 : 0; @@ -348,7 +347,31 @@ namespace { iWell[Ix::item32] = 7; iWell[Ix::item48] = - 1; - iWell[Ix::item50] = iWell[Ix::WCtrl]; + // Deliberate misrepresentation. Function 'ctrlMode()' returns + // the target control mode requested in the simulation deck. + // This item is supposed to be the well's actual, active target + // control mode in the simulator. + iWell[Ix::ActWCtrl] = ctrlMode(well, sim_step); + + const auto isPred = + (well.isProducer(sim_step) && + well.getProductionProperties(sim_step).predictionMode) + || + (well.isInjector(sim_step) && + well.getInjectionProperties(sim_step).predictionMode); + + if (isPred) { + // Well in prediction mode (WCONPROD, WCONINJE). Assign + // requested control mode for prediction. + iWell[Ix::PredReqWCtrl] = iWell[Ix::ActWCtrl]; + iWell[Ix::HistReqWCtrl] = 0; + } + else { + // Well controlled by observed rates/BHP (WCONHIST, + // WCONINJH). Assign requested control mode for history. + iWell[Ix::PredReqWCtrl] = 0; // Possibly =1 instead. + iWell[Ix::HistReqWCtrl] = iWell[Ix::ActWCtrl]; + } // Multi-segmented well information iWell[Ix::MsWID] = 0; // MS Well ID (0 or 1..#MS wells) @@ -386,7 +409,7 @@ namespace { }); iWell[Ix::item9] = any_flowing_conn - ? iWell[Ix::WCtrl] : -1; + ? iWell[Ix::ActWCtrl] : -1; iWell[Ix::item11] = 1; } diff --git a/tests/test_AggregateWellData.cpp b/tests/test_AggregateWellData.cpp index 3bc0d498a..3a8be3c56 100644 --- a/tests/test_AggregateWellData.cpp +++ b/tests/test_AggregateWellData.cpp @@ -548,7 +548,7 @@ BOOST_AUTO_TEST_CASE (Dynamic_Well_Data_Step1) const auto& iwell = awd.getIWell(); - BOOST_CHECK_EQUAL(iwell[i0 + Ix::item9 ], iwell[i0 + Ix::WCtrl]); + BOOST_CHECK_EQUAL(iwell[i0 + Ix::item9 ], iwell[i0 + Ix::ActWCtrl]); BOOST_CHECK_EQUAL(iwell[i0 + Ix::item11], 1); } @@ -668,7 +668,7 @@ BOOST_AUTO_TEST_CASE (Dynamic_Well_Data_Step2) const auto& iwell = awd.getIWell(); BOOST_CHECK_EQUAL(iwell[i1 + Ix::item9], - iwell[i1 + Ix::WCtrl]); + iwell[i1 + Ix::ActWCtrl]); BOOST_CHECK_EQUAL(iwell[i1 + Ix::item11], 1); } From c9331d7fedba8ae73f72fa6cf74cf4a843fa17af Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Tue, 22 Jan 2019 20:23:12 +0100 Subject: [PATCH 17/20] Save: Simplify Main Gateway Function In particular, defer unit conversion to the member function RestartValue::convertFromSI() and detect presence of multisegment wells using std::any_of(). We don't need to know the exact number of MS wells to activate segment output. --- src/opm/output/eclipse/RestartIO.cpp | 14 ++++---------- 1 file changed, 4 insertions(+), 10 deletions(-) diff --git a/src/opm/output/eclipse/RestartIO.cpp b/src/opm/output/eclipse/RestartIO.cpp index 530190ad6..05b33b81b 100644 --- a/src/opm/output/eclipse/RestartIO.cpp +++ b/src/opm/output/eclipse/RestartIO.cpp @@ -538,13 +538,7 @@ void save(const std::string& filename, write_double = false; // Convert solution fields and extra values from SI to user units. - value.solution.convertFromSI(units); - for (auto& extra_value : value.extra) { - const auto& restart_key = extra_value.first; - auto& data = extra_value.second; - - units.from_si(restart_key.dim, data); - } + value.convertFromSI(units); const auto inteHD = writeHeader(rst_file.get(), sim_step, report_step, @@ -559,14 +553,14 @@ void save(const std::string& filename, const auto& wells = schedule.getWells(sim_step); if (! wells.empty()) { - const auto numMSW = - std::count_if(std::begin(wells), std::end(wells), + const auto haveMSW = + std::any_of(std::begin(wells), std::end(wells), [sim_step](const Well* well) { return well->isMultiSegment(sim_step); }); - if (numMSW > 0) { + if (haveMSW) { writeMSWData(rst_file.get(), sim_step, units, schedule, grid, sumState, value.wells, inteHD); } From 0d0b8b257f47580d70f7209bf581ea4b4f6b7986 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Tue, 22 Jan 2019 20:28:19 +0100 Subject: [PATCH 18/20] Well Data: Prune Disabled Code Blocks These have not been used for a long time. --- src/opm/output/eclipse/AggregateWellData.cpp | 33 -------------------- 1 file changed, 33 deletions(-) diff --git a/src/opm/output/eclipse/AggregateWellData.cpp b/src/opm/output/eclipse/AggregateWellData.cpp index 1e9a4594c..cc41da631 100644 --- a/src/opm/output/eclipse/AggregateWellData.cpp +++ b/src/opm/output/eclipse/AggregateWellData.cpp @@ -135,29 +135,6 @@ namespace { } return ind; } - /*int groupIndex(const std::string& grpName, - const std::vector& groupNames, - const int maxGroups) - { - if (grpName == "FIELD") { - // Not really supposed to happen since wells are - // not supposed to be parented dirctly to FIELD. - return maxGroups + 1; - } - - auto b = std::begin(groupNames); - auto e = std::end (groupNames); - auto i = std::find(b, e, grpName); - - if (i == e) { - // Not really supposed to happen since wells are - // not supposed to be parented dirctly to FIELD. - return maxGroups + 1; - } - - // One-based indexing. - return std::distance(b, i) + 1; - } */ int wellType(const Opm::Well& well, const std::size_t sim_step) @@ -198,16 +175,6 @@ namespace { using WMCtrlVal = ::Opm::RestartIO::Helpers:: VectorItems::IWell::Value::WellCtrlMode; - /*{ - const auto stat = well.getStatus(sim_step); - - using WStat = ::Opm::WellCommon::StatusEnum; - - if (stat == WStat::SHUT) { - return WMCtrlVal::Shut; - } - }*/ - if (well.isInjector(sim_step)) { const auto& prop = well .getInjectionProperties(sim_step); From a51d18cab0b23f570633d251a0fbab3f8ee06709 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Tue, 22 Jan 2019 20:45:28 +0100 Subject: [PATCH 19/20] Restart Output: Adjust Whitespace Mostly removing leading 'tab' characters, splitting long lines and realigning '=' characters where appropriate. Also make a few more objects 'const'. Finally, don't copy large objects (e.g., the Schedule) when we only need to read from it. --- opm/output/eclipse/VectorItems/well.hpp | 21 ++- src/opm/output/eclipse/AggregateGroupData.cpp | 37 +++-- src/opm/output/eclipse/AggregateWellData.cpp | 150 +++++++++--------- src/opm/output/eclipse/RestartIO.cpp | 50 +++--- 4 files changed, 136 insertions(+), 122 deletions(-) diff --git a/opm/output/eclipse/VectorItems/well.hpp b/opm/output/eclipse/VectorItems/well.hpp index 4a7d1564b..1e3ec5dbb 100644 --- a/opm/output/eclipse/VectorItems/well.hpp +++ b/opm/output/eclipse/VectorItems/well.hpp @@ -26,19 +26,18 @@ namespace Opm { namespace RestartIO { namespace Helpers { namespace VectorItems namespace IWell { enum index : std::vector::size_type { - IHead = 0, // I-location (one-based) of well head - JHead = 1, // J-location (one-based) of well head - FirstK = 2, // Layer ID (one-based) of top/first connection - LastK = 3, // Layer ID (one-based) of bottom/last connection - NConn = 4, // Number of active cells connected to well - Group = 5, // Index (one-based) of well's current group - WType = 6, // Well type + IHead = 0, // I-location (one-based) of well head + JHead = 1, // J-location (one-based) of well head + FirstK = 2, // Layer ID (one-based) of top/first connection + LastK = 3, // Layer ID (one-based) of bottom/last connection + NConn = 4, // Number of active cells connected to well + Group = 5, // Index (one-based) of well's current group + WType = 6, // Well type (producer vs. injector) ActWCtrl = 7, // Well's active target control mode (constraint). - item9 = 8, // Unknown - item11 = 10, // Unknown - - VFPTab = 11, // ID (one-based) of well's current VFP table + item9 = 8, // Unknown + item11 = 10, // Unknown + VFPTab = 11, // ID (one-based) of well's current VFP table. PredReqWCtrl = 15, // Well's requested control mode from // simulation deck (WCONINJE, WCONPROD). diff --git a/src/opm/output/eclipse/AggregateGroupData.cpp b/src/opm/output/eclipse/AggregateGroupData.cpp index f55bc5d8d..c6c41eb99 100644 --- a/src/opm/output/eclipse/AggregateGroupData.cpp +++ b/src/opm/output/eclipse/AggregateGroupData.cpp @@ -521,19 +521,20 @@ captureDeclaredGroupData(const Opm::Schedule& sched, auto it = indexGroupMap.begin(); while (it != indexGroupMap.end()) - { - curGroups[static_cast(it->first)] = it->second; - it++; - } { - groupLoop(curGroups, [sched, simStep, inteHead, this] - (const Group& group, const std::size_t groupID) -> void - { - auto ig = this->iGroup_[groupID]; - IGrp::staticContrib(sched, group, this->nWGMax_, this->nGMaxz_, simStep, ig, inteHead); - }); + curGroups[static_cast(it->first)] = it->second; + it++; } + groupLoop(curGroups, [&sched, simStep, &inteHead, this] + (const Group& group, const std::size_t groupID) -> void + { + auto ig = this->iGroup_[groupID]; + + IGrp::staticContrib(sched, group, this->nWGMax_, this->nGMaxz_, + simStep, ig, inteHead); + }); + // Define Static Contributions to SGrp Array. groupLoop(curGroups, [this](const Group& group, const std::size_t groupID) -> void @@ -542,15 +543,19 @@ captureDeclaredGroupData(const Opm::Schedule& sched, SGrp::staticContrib(sw); }); - // Define DynamicContributions to XGrp Array. - groupLoop(curGroups, - [restart_group_keys, restart_field_keys, groupKeyToIndex, fieldKeyToIndex, ecl_compatible_rst, sumState, this] + // Define Dynamic Contributions to XGrp Array. + groupLoop(curGroups, [&restart_group_keys, &restart_field_keys, + &groupKeyToIndex, &fieldKeyToIndex, + ecl_compatible_rst, &sumState, this] (const Group& group, const std::size_t groupID) -> void { auto xg = this->xGroup_[groupID]; - XGrp::dynamicContrib( restart_group_keys, restart_field_keys, groupKeyToIndex, fieldKeyToIndex, group, sumState, ecl_compatible_rst, xg); + + XGrp::dynamicContrib(restart_group_keys, restart_field_keys, + groupKeyToIndex, fieldKeyToIndex, group, + sumState, ecl_compatible_rst, xg); }); - + // Define Static Contributions to ZGrp Array. groupLoop(curGroups, [this](const Group& group, const std::size_t groupID) -> void @@ -560,5 +565,3 @@ captureDeclaredGroupData(const Opm::Schedule& sched, }); } - -// --------------------------------------------------------------------- diff --git a/src/opm/output/eclipse/AggregateWellData.cpp b/src/opm/output/eclipse/AggregateWellData.cpp index cc41da631..9168f40f4 100644 --- a/src/opm/output/eclipse/AggregateWellData.cpp +++ b/src/opm/output/eclipse/AggregateWellData.cpp @@ -117,24 +117,22 @@ namespace { groupIndexMap.insert(groupPair); } return groupIndexMap; - } - - int groupIndex(const std::string& grpName, - const std::map & currentGroupMapNameIndex) - { + } + + int groupIndex(const std::string& grpName, + const std::map & currentGroupMapNameIndex) + { int ind = 0; - auto searchGTName = currentGroupMapNameIndex.find(grpName); - if (searchGTName != currentGroupMapNameIndex.end()) - { - ind = searchGTName->second + 1; + auto searchGTName = currentGroupMapNameIndex.find(grpName); + if (searchGTName != currentGroupMapNameIndex.end()) { + ind = searchGTName->second + 1; } - else - { + else { std::cout << "group Name: " << grpName << std::endl; throw std::invalid_argument( "Invalid group name" ); } - return ind; - } + return ind; + } int wellType(const Opm::Well& well, const std::size_t sim_step) @@ -202,16 +200,16 @@ namespace { case CMode::GRUP: return WMCtrlVal::Group; default: - { - const auto stat = well.getStatus(sim_step); + { + const auto stat = well.getStatus(sim_step); - using WStat = ::Opm::WellCommon::StatusEnum; + using WStat = ::Opm::WellCommon::StatusEnum; - if (stat == WStat::SHUT) { - return WMCtrlVal::Shut; - } - } - return WMCtrlVal::WMCtlUnk; + if (stat == WStat::SHUT) { + return WMCtrlVal::Shut; + } + } + return WMCtrlVal::WMCtlUnk; } } else if (well.isProducer(sim_step)) { @@ -229,19 +227,19 @@ namespace { case CMode::THP: return WMCtrlVal::THP; case CMode::BHP: return WMCtrlVal::BHP; case CMode::CRAT: return WMCtrlVal::CombRate; - case CMode::GRUP: return WMCtrlVal::Group; - - default: - { - const auto stat = well.getStatus(sim_step); + case CMode::GRUP: return WMCtrlVal::Group; - using WStat = ::Opm::WellCommon::StatusEnum; + default: + { + const auto stat = well.getStatus(sim_step); - if (stat == WStat::SHUT) { - return WMCtrlVal::Shut; - } - } - return WMCtrlVal::WMCtlUnk; + using WStat = ::Opm::WellCommon::StatusEnum; + + if (stat == WStat::SHUT) { + return WMCtrlVal::Shut; + } + } + return WMCtrlVal::WMCtlUnk; } } @@ -282,30 +280,29 @@ namespace { const auto& conn = well.getConnections(sim_step); iWell[Ix::NConn] = static_cast(conn.size()); - - if (well.isMultiSegment(sim_step)) - { - // Set top and bottom connections to zero for multi segment wells - iWell[Ix::FirstK] = 0; - iWell[Ix::LastK] = 0; - } - else - { - iWell[Ix::FirstK] = (iWell[Ix::NConn] == 0) - ? 0 : conn.get(0).getK() + 1; - iWell[Ix::LastK] = (iWell[Ix::NConn] == 0) - ? 0 : conn.get(conn.size() - 1).getK() + 1; - } + if (well.isMultiSegment(sim_step)) { + // Set top and bottom connections to zero for multi + // segment wells + iWell[Ix::FirstK] = 0; + iWell[Ix::LastK] = 0; + } + else { + iWell[Ix::FirstK] = (iWell[Ix::NConn] == 0) + ? 0 : conn.get(0).getK() + 1; + + iWell[Ix::LastK] = (iWell[Ix::NConn] == 0) + ? 0 : conn.get(conn.size() - 1).getK() + 1; + } } iWell[Ix::Group] = - groupIndex(trim(well.getGroupName(sim_step)), - GroupMapNameInd); + groupIndex(trim(well.getGroupName(sim_step)), + GroupMapNameInd); iWell[Ix::WType] = wellType (well, sim_step); iWell[Ix::VFPTab] = wellVFPTab(well, sim_step); - iWell[Ix::XFlow] = well.getAllowCrossFlow() ? 1 : 0; + iWell[Ix::XFlow] = well.getAllowCrossFlow() ? 1 : 0; // The following items aren't fully characterised yet, but // needed for restart of M2. Will need further refinement. @@ -470,7 +467,7 @@ namespace { void staticContrib(const Opm::Well& well, const Opm::UnitSystem& units, const std::size_t sim_step, - const ::Opm::SummaryState& smry, + const ::Opm::SummaryState& smry, SWellArray& sWell) { using Ix = ::Opm::RestartIO::Helpers::VectorItems::SWell::index; @@ -480,17 +477,17 @@ namespace { { return static_cast(units.from_si(u, x)); }; - + assignDefaultSWell(sWell); if (well.isProducer(sim_step)) { const auto& pp = well.getProductionProperties(sim_step); - const auto& predMode = pp.predictionMode; + const auto& predMode = pp.predictionMode; if ((pp.OilRate != 0.0) || (!predMode)) { sWell[Ix::OilRateTarget] = swprop(M::liquid_surface_rate, pp.OilRate); - } + } if ((pp.WaterRate != 0.0) || (!predMode)) { sWell[Ix::WatRateTarget] = @@ -500,13 +497,13 @@ namespace { if ((pp.GasRate != 0.0) || (!predMode)) { sWell[Ix::GasRateTarget] = swprop(M::gas_surface_rate, pp.GasRate); - sWell[Ix::HistGasRateTarget] = sWell[Ix::GasRateTarget]; + sWell[Ix::HistGasRateTarget] = sWell[Ix::GasRateTarget]; } if (pp.LiquidRate != 0.0 || (!predMode)) { sWell[Ix::LiqRateTarget] = swprop(M::liquid_surface_rate, pp.LiquidRate); - sWell[Ix::HistLiqRateTarget] = sWell[Ix::LiqRateTarget]; + sWell[Ix::HistLiqRateTarget] = sWell[Ix::LiqRateTarget]; } else { sWell[Ix::LiqRateTarget] = @@ -517,21 +514,23 @@ namespace { sWell[Ix::ResVRateTarget] = swprop(M::rate, pp.ResVRate); } - else if ((smry.has("WVPR:" + well.name())) && (!predMode)) { + else if ((smry.has("WVPR:" + well.name())) && (!predMode)) { // Write out summary voidage production rate if // target/limit is not set - auto vr = static_cast(smry.get("WVPR:" + well.name())); - if (vr != 0.0) sWell[Ix::ResVRateTarget] = vr; + auto vr = static_cast(smry.get("WVPR:" + well.name())); + if (vr != 0.0) { + sWell[Ix::ResVRateTarget] = vr; + } } - sWell[Ix::THPTarget] = pp.THPLimit != 0.0 + sWell[Ix::THPTarget] = pp.THPLimit != 0.0 ? swprop(M::pressure, pp.THPLimit) - : 0.; - - sWell[Ix::BHPTarget] = pp.BHPLimit != 0.0 + : 0.0; + + sWell[Ix::BHPTarget] = pp.BHPLimit != 0.0 ? swprop(M::pressure, pp.BHPLimit) : swprop(M::pressure, 1.0*::Opm::unit::atm); - sWell[Ix::HistBHPTarget] = sWell[Ix::BHPTarget]; + sWell[Ix::HistBHPTarget] = sWell[Ix::BHPTarget]; } else if (well.isInjector(sim_step)) { const auto& ip = well.getInjectionProperties(sim_step); @@ -547,18 +546,19 @@ namespace { if (ip.injectorType == IT::WATER) { sWell[Ix::WatRateTarget] = swprop(M::liquid_surface_rate, ip.surfaceInjectionRate); - sWell[Ix::HistLiqRateTarget] = sWell[Ix::WatRateTarget]; + sWell[Ix::HistLiqRateTarget] = sWell[Ix::WatRateTarget]; } if (ip.injectorType == IT::GAS) { sWell[Ix::GasRateTarget] = swprop(M::gas_surface_rate, ip.surfaceInjectionRate); - sWell[Ix::HistGasRateTarget] = sWell[Ix::GasRateTarget]; + sWell[Ix::HistGasRateTarget] = sWell[Ix::GasRateTarget]; } } - if (ip.hasInjectionControl(IP::RESV)) { - sWell[Ix::ResVRateTarget] = swprop(M::rate, ip.reservoirInjectionRate); - } + if (ip.hasInjectionControl(IP::RESV)) { + sWell[Ix::ResVRateTarget] = + swprop(M::rate, ip.reservoirInjectionRate); + } if (ip.hasInjectionControl(IP::THP)) { sWell[Ix::THPTarget] = swprop(M::pressure, ip.THPLimit); @@ -567,7 +567,7 @@ namespace { sWell[Ix::BHPTarget] = ip.hasInjectionControl(IP::BHP) ? swprop(M::pressure, ip.BHPLimit) : swprop(M::pressure, 1.0E05*::Opm::unit::psia); - sWell[Ix::HistBHPTarget] = sWell[Ix::BHPTarget]; + sWell[Ix::HistBHPTarget] = sWell[Ix::BHPTarget]; } sWell[Ix::DatumDepth] = @@ -666,14 +666,14 @@ namespace { xWell[Ix::FlowBHP] = get("WBHP"); - if (ecl_compatible_rst) { + if (ecl_compatible_rst) { // Note: Assign both water and gas cumulatives to support // case of well alternating between injecting water and gas. xWell[Ix::WatInjTotal] = get("WWIT"); xWell[Ix::GasInjTotal] = get("WGIT"); xWell[Ix::HistWatInjTotal] = get("WWITH"); xWell[Ix::HistGasInjTotal] = get("WGITH"); - } + } } template @@ -727,9 +727,11 @@ namespace { xWell[Ix::GasFVF] = (std::abs(xWell[Ix::GasPrRate]) > 0.0) ? xWell[Ix::VoidPrRate] / xWell[Ix::GasPrRate] : 0.0; - - if (std::isnan(xWell[Ix::GasFVF])) xWell[Ix::GasFVF] = 0.; - + + if (std::isnan(xWell[Ix::GasFVF])) { + xWell[Ix::GasFVF] = 0.0; + } + // Not fully characterised. xWell[Ix::item38] = xWell[Ix::GasPrRate]; diff --git a/src/opm/output/eclipse/RestartIO.cpp b/src/opm/output/eclipse/RestartIO.cpp index 05b33b81b..eb877883f 100644 --- a/src/opm/output/eclipse/RestartIO.cpp +++ b/src/opm/output/eclipse/RestartIO.cpp @@ -283,7 +283,8 @@ namespace { } // write INTEHEAD to restart file - const auto ih = Helpers::createInteHead(es, grid, schedule, simTime, sim_step, sim_step); + const auto ih = Helpers::createInteHead(es, grid, schedule, + simTime, sim_step, sim_step); write_kw(rst_file, "INTEHEAD", ih); // write LOGIHEAD to restart file @@ -301,7 +302,7 @@ namespace { void writeGroup(::Opm::RestartIO::ecl_rst_file_type* rst_file, int sim_step, - const bool ecl_compatible_rst, + const bool ecl_compatible_rst, const Schedule& schedule, const Opm::SummaryState& sumState, const std::vector& ih) @@ -311,15 +312,15 @@ namespace { auto groupData = Helpers::AggregateGroupData(ih); - auto & rst_g_keys = groupData.restart_group_keys; - auto & rst_f_keys = groupData.restart_field_keys; - auto & grpKeyToInd = groupData.groupKeyToIndex; - auto & fldKeyToInd = groupData.fieldKeyToIndex; + const auto& rst_g_keys = groupData.restart_group_keys; + const auto& rst_f_keys = groupData.restart_field_keys; + const auto& grpKeyToInd = groupData.groupKeyToIndex; + const auto& fldKeyToInd = groupData.fieldKeyToIndex; groupData.captureDeclaredGroupData(schedule, rst_g_keys, rst_f_keys, grpKeyToInd, fldKeyToInd, - ecl_compatible_rst, + ecl_compatible_rst, simStep, sumState, ih); write_kw(rst_file, "IGRP", groupData.getIGroup()); write_kw(rst_file, "SGRP", groupData.getSGroup()); @@ -332,14 +333,16 @@ namespace { const UnitSystem& units, const Schedule& schedule, const EclipseGrid& grid, - const Opm::SummaryState& sumState, - const Opm::data::Wells& wells, + const Opm::SummaryState& sumState, + const Opm::data::Wells& wells, const std::vector& ih) { // write ISEG, RSEG, ILBS and ILBR to restart file - const size_t simStep = static_cast (sim_step); + const auto simStep = static_cast (sim_step); + auto MSWData = Helpers::AggregateMSWData(ih); - MSWData.captureDeclaredMSWData(schedule, simStep, units, ih, grid, sumState, wells); + MSWData.captureDeclaredMSWData(schedule, simStep, units, + ih, grid, sumState, wells); write_kw(rst_file, "ISEG", MSWData.getISeg()); write_kw(rst_file, "ILBS", MSWData.getILBs()); @@ -349,7 +352,7 @@ namespace { void writeWell(::Opm::RestartIO::ecl_rst_file_type* rst_file, int sim_step, - const bool ecl_compatible_rst, + const bool ecl_compatible_rst, const Phases& phases, const UnitSystem& units, const EclipseGrid& grid, @@ -360,7 +363,9 @@ namespace { { auto wellData = Helpers::AggregateWellData(ih); wellData.captureDeclaredWellData(schedule, units, sim_step, sumState, ih); - wellData.captureDynamicWellData(schedule, sim_step, ecl_compatible_rst, wells, sumState); + wellData.captureDynamicWellData(schedule, sim_step, + ecl_compatible_rst, + wells, sumState); write_kw(rst_file, "IWEL", wellData.getIWell()); write_kw(rst_file, "SWEL", wellData.getSWell()); @@ -368,7 +373,7 @@ namespace { write_kw(rst_file, "ZWEL", serialize_ZWEL(wellData.getZWell())); // Extended set of OPM well vectors - if (!ecl_compatible_rst) + if (!ecl_compatible_rst) { const auto sched_wells = schedule.getWells(sim_step); @@ -382,7 +387,8 @@ namespace { } auto connectionData = Helpers::AggregateConnectionData(ih); - connectionData.captureDeclaredConnData(schedule, grid, units, wells, sim_step); + connectionData.captureDeclaredConnData(schedule, grid, units, + wells, sim_step); write_kw(rst_file, "ICON", connectionData.getIConn()); write_kw(rst_file, "SCON", connectionData.getSConn()); @@ -529,13 +535,16 @@ void save(const std::string& filename, bool write_double) { ::Opm::RestartIO::checkSaveArguments(es, value, grid); - bool ecl_compatible_rst = es.getIOConfig().getEclCompatibleRST(); + + const auto ecl_compatible_rst = es.getIOConfig().getEclCompatibleRST(); + const auto sim_step = std::max(report_step - 1, 0); const auto& units = es.getUnits(); auto rst_file = openRestartFile(filename, report_step); - if (ecl_compatible_rst) - write_double = false; + if (ecl_compatible_rst) { + write_double = false; + } // Convert solution fields and extra values from SI to user units. value.convertFromSI(units); @@ -573,8 +582,9 @@ void save(const std::string& filename, writeSolution(rst_file.get(), value, ecl_compatible_rst, write_double); - if (!ecl_compatible_rst) - ::Opm::RestartIO::writeExtraData(rst_file.get(), value.extra); + if (! ecl_compatible_rst) { + writeExtraData(rst_file.get(), value.extra); + } } }} // Opm::RestartIO From c345aa3f4eaf86251c5b9901b8820cb6cdf4e28b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Tue, 22 Jan 2019 20:48:25 +0100 Subject: [PATCH 20/20] Restart Output: Move Group Name Assignment to Separate Helper Mostly for consistency. Also, look up window index by group name. --- src/opm/output/eclipse/AggregateGroupData.cpp | 16 +++++++++++----- 1 file changed, 11 insertions(+), 5 deletions(-) diff --git a/src/opm/output/eclipse/AggregateGroupData.cpp b/src/opm/output/eclipse/AggregateGroupData.cpp index c6c41eb99..80f5c75d9 100644 --- a/src/opm/output/eclipse/AggregateGroupData.cpp +++ b/src/opm/output/eclipse/AggregateGroupData.cpp @@ -448,6 +448,12 @@ namespace { WV::WindowSize{ entriesPerGroup(inteHead) } }; } + + template + void staticContrib(const Opm::Group& group, ZGroupArray& zGroup) + { + zGroup[0] = group.name(); + } } // ZGrp } // Anonymous @@ -557,11 +563,11 @@ captureDeclaredGroupData(const Opm::Schedule& sched, }); // Define Static Contributions to ZGrp Array. - groupLoop(curGroups, - [this](const Group& group, const std::size_t groupID) -> void + groupLoop(curGroups, [this, &nameIndexMap] + (const Group& group, const std::size_t /* groupID */) -> void { - auto zw = this->zGroup_[groupID]; - zw[0] = group.name(); - }); + auto zg = this->zGroup_[ nameIndexMap.at(group.name()) ]; + ZGrp::staticContrib(group, zg); + }); }