From 9b034978bdb05a3a882128ce552a438cb7d7b47a Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Tue, 14 Nov 2023 09:20:24 +0100 Subject: [PATCH] Load D-Factor Correlation Parameters From Restart File This commit loads the connection level Peaceman denominator value, static D-factor correlation coefficient, and connection interval length from the restart file and forwards these values from the RstConnection type to the regular Connection type. Furthermore, we load the well-level D-factor correlation parameters 'A', 'B', and 'C' of the WDFACCOR keyword into the WDFAC type at restart time. Collectively, these changes enable restarting simulation runs using the WDFACCOR keyword. --- opm/input/eclipse/Schedule/Well/WDFAC.hpp | 11 ++ opm/io/eclipse/rst/connection.hpp | 35 ++-- opm/io/eclipse/rst/well.hpp | 3 + .../eclipse/Schedule/Well/Connection.cpp | 9 +- src/opm/input/eclipse/Schedule/Well/WDFAC.cpp | 12 ++ src/opm/input/eclipse/Schedule/Well/Well.cpp | 2 +- src/opm/io/eclipse/rst/connection.cpp | 169 +++++++++++------- src/opm/io/eclipse/rst/state.cpp | 3 +- src/opm/io/eclipse/rst/well.cpp | 159 +++++++++------- 9 files changed, 252 insertions(+), 151 deletions(-) diff --git a/opm/input/eclipse/Schedule/Well/WDFAC.hpp b/opm/input/eclipse/Schedule/Well/WDFAC.hpp index 22c5cc3b6..cfa1e1a62 100644 --- a/opm/input/eclipse/Schedule/Well/WDFAC.hpp +++ b/opm/input/eclipse/Schedule/Well/WDFAC.hpp @@ -101,6 +101,17 @@ namespace Opm { } }; + /// Default constructor + WDFAC() = default; + + /// Constructor + /// + /// Creates an object from restart information + /// + /// \param[in] rst_well Linearised well-level restart information, + /// including D-factor parameters. + explicit WDFAC(const RestartIO::RstWell& rst_well); + /// Serialisation test object static WDFAC serializationTestObject(); diff --git a/opm/io/eclipse/rst/connection.hpp b/opm/io/eclipse/rst/connection.hpp index 0b5191475..9bd54dfa7 100644 --- a/opm/io/eclipse/rst/connection.hpp +++ b/opm/io/eclipse/rst/connection.hpp @@ -3,7 +3,8 @@ This file is part of the Open Porous Media project (OPM). - OPM is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. @@ -19,19 +20,27 @@ #ifndef RST_CONNECTION #define RST_CONNECTION -#include - #include +#include + namespace Opm { + class UnitSystem; -namespace RestartIO { +} // namespace Opm -class Header; +namespace Opm { namespace RestartIO { + +struct RstConnection +{ + RstConnection(const UnitSystem& unit_system, + std::size_t rst_index, + int nsconz, + const int* icon, + const float* scon, + const double *xcon); -struct RstConnection { - RstConnection(const ::Opm::UnitSystem& unit_system, std::size_t rst_index, int nsconz, const int* icon, const float* scon, const double *xcon); static double inverse_peaceman(double cf, double kh, double rw, double skin); std::size_t rst_index; @@ -49,6 +58,9 @@ struct RstConnection { float depth; float diameter; float kh; + float denom; + float length; + float static_dfac_corr_coeff; float segdist_end; float segdist_start; @@ -60,11 +72,6 @@ struct RstConnection { double r0; }; +}} // namespace Opm::RestartIO -} -} - - - - -#endif +#endif // RST_CONNECTION diff --git a/opm/io/eclipse/rst/well.hpp b/opm/io/eclipse/rst/well.hpp index f38bb333c..fb3a826e1 100644 --- a/opm/io/eclipse/rst/well.hpp +++ b/opm/io/eclipse/rst/well.hpp @@ -126,6 +126,9 @@ struct RstWell float glift_min_rate; float glift_weight_factor; float glift_inc_weight_factor; + float dfac_corr_coeff_a{}; + float dfac_corr_exponent_b{}; + float dfac_corr_exponent_c{}; std::vector tracer_concentration_injection; double oil_rate; diff --git a/src/opm/input/eclipse/Schedule/Well/Connection.cpp b/src/opm/input/eclipse/Schedule/Well/Connection.cpp index be7712392..5acab796f 100644 --- a/src/opm/input/eclipse/Schedule/Well/Connection.cpp +++ b/src/opm/input/eclipse/Schedule/Well/Connection.cpp @@ -31,12 +31,15 @@ #include #include +#include #include +#include #include #include #include #include #include +#include #include namespace { @@ -53,9 +56,11 @@ namespace { props.rw = rst_conn.diameter / 2; props.r0 = rst_conn.r0; props.re = 0.0; - props.connection_length = 0.0; + props.connection_length = rst_conn.length; props.skin_factor = rst_conn.skin_factor; props.d_factor = 0.0; + props.static_dfac_corr_coeff = rst_conn.static_dfac_corr_coeff; + props.peaceman_denom = rst_conn.denom; return props; } @@ -155,7 +160,7 @@ namespace Opm rst_connection.segdist_end); } - //TODO recompute re and perf_length from the grid + // TODO: recompute re from the grid } Connection Connection::serializationTestObject() diff --git a/src/opm/input/eclipse/Schedule/Well/WDFAC.cpp b/src/opm/input/eclipse/Schedule/Well/WDFAC.cpp index 5c0295779..194e9f1f3 100644 --- a/src/opm/input/eclipse/Schedule/Well/WDFAC.cpp +++ b/src/opm/input/eclipse/Schedule/Well/WDFAC.cpp @@ -72,6 +72,18 @@ namespace Opm { // ------------------------------------------------------------------------- + WDFAC::WDFAC(const RestartIO::RstWell& rst_well) + : m_type { WDFacType::NONE } + , m_d { 0.0 } + , m_corr { rst_well.dfac_corr_coeff_a , + rst_well.dfac_corr_exponent_b , + rst_well.dfac_corr_exponent_c } + { + if (m_corr.coeff_a > 0.0) { + this->m_type = WDFacType::DAKEMODEL; + } + } + WDFAC WDFAC::serializationTestObject() { WDFAC result; diff --git a/src/opm/input/eclipse/Schedule/Well/Well.cpp b/src/opm/input/eclipse/Schedule/Well/Well.cpp index 84ef0354d..b98fa5a8c 100644 --- a/src/opm/input/eclipse/Schedule/Well/Well.cpp +++ b/src/opm/input/eclipse/Schedule/Well/Well.cpp @@ -313,7 +313,7 @@ Well::Well(const RestartIO::RstWell& rst_well, injection(std::make_shared(unit_system_arg, wname)), wvfpdp(std::make_shared()), wvfpexp(explicitTHPOptions(rst_well)), - wdfac(std::make_shared()), + wdfac(std::make_shared(rst_well)), status(status_from_int(rst_well.well_status)), well_temperature(Metric::TemperatureOffset + ParserKeywords::STCOND::TEMPERATURE::defaultValue), well_inj_mult(std::nullopt) diff --git a/src/opm/io/eclipse/rst/connection.cpp b/src/opm/io/eclipse/rst/connection.cpp index 92db47826..2b21522a7 100644 --- a/src/opm/io/eclipse/rst/connection.cpp +++ b/src/opm/io/eclipse/rst/connection.cpp @@ -3,7 +3,8 @@ This file is part of the Open Porous Media project (OPM). - OPM is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. @@ -16,99 +17,129 @@ along with OPM. If not, see . */ -#include -#include - -#include #include + #include + #include +#include +#include +#include +#include +#include namespace VI = ::Opm::RestartIO::Helpers::VectorItems; -namespace Opm { -namespace RestartIO { - namespace { template T from_int(int); -template <> -Connection::State from_int(int int_state) { - if (int_state == 1) - return Connection::State::OPEN; - - return Connection::State::SHUT; +template<> Opm::Connection::State from_int(int int_state) +{ + return (int_state == 1) + ? Opm::Connection::State::OPEN + : Opm::Connection::State::SHUT; } -template <> -Connection::Direction from_int(int int_dir) { +template<> Opm::Connection::Direction from_int(int int_dir) +{ switch (int_dir) { - case 1: - return Connection::Direction::X; - case 2: - return Connection::Direction::Y; - case 3: - return Connection::Direction::Z; - default: - throw std::invalid_argument("Can not convert: " + std::to_string(int_dir) + " to string"); + case 1: return Opm::Connection::Direction::X; + case 2: return Opm::Connection::Direction::Y; + case 3: return Opm::Connection::Direction::Z; } + + throw std::invalid_argument { + "Unable to convert direction value: " + + std::to_string(int_dir) + " to Direction category" + }; } -/* - That the CTFKind variable comes from a float looks extremely suspicious; but - it has been double checked ... -*/ -Connection::CTFKind from_float(float float_kind) { - if (float_kind == 0) - return Connection::CTFKind::Defaulted; - - return Connection::CTFKind::DeckValue; -} +// Note: CTFKind originates in SCON and is indeed a float. +Opm::Connection::CTFKind from_float(float float_kind) +{ + return (float_kind == 0.0f) + ? Opm::Connection::CTFKind::Defaulted + : Opm::Connection::CTFKind::DeckValue; } -double RstConnection::inverse_peaceman(double cf, double kh, double rw, double skin) { +float as_float(const double x) +{ + return static_cast(x); +} + +float staticDFactorCorrCoeff(const Opm::UnitSystem& usys, + const float coeff) +{ + using M = ::Opm::UnitSystem::measure; + + // Coefficient's units are [D] * [viscosity] + return as_float(usys.to_si(M::viscosity, usys.to_si(M::dfactor, coeff))); +} + +double pressEquivRadius(const float denom, + const float skin, + const float rw) +{ + // Recall: denom = log(r0 / rw) + skin + return rw * std::exp(denom - skin); +} + +} // Anonymous namespace + +double Opm::RestartIO::RstConnection::inverse_peaceman(double cf, double kh, double rw, double skin) +{ auto alpha = 3.14159265 * 2 * kh / cf - skin; return rw * std::exp(alpha); } -using M = ::Opm::UnitSystem::measure; +using M = ::Opm::UnitSystem::measure; -RstConnection::RstConnection(const ::Opm::UnitSystem& unit_system, std::size_t rst_index_, int nsconz, const int* icon, const float* scon, const double* xcon) : - rst_index( rst_index_), - ijk( {icon[VI::IConn::CellI] - 1, icon[VI::IConn::CellJ] - 1, icon[VI::IConn::CellK] - 1}), - state( from_int(icon[VI::IConn::ConnStat])), - drain_sat_table( icon[VI::IConn::Drainage]), - imb_sat_table( icon[VI::IConn::Imbibition]), - completion( icon[VI::IConn::ComplNum]), - dir( from_int(icon[VI::IConn::ConnDir])), - segment( icon[VI::IConn::Segment]), - cf_kind( from_float(1.0)), - skin_factor( scon[VI::SConn::SkinFactor]), - cf( unit_system.to_si(M::transmissibility, scon[VI::SConn::ConnTrans])), - depth( unit_system.to_si(M::length, scon[VI::SConn::Depth])), - diameter( unit_system.to_si(M::length, scon[VI::SConn::Diameter])), - kh( unit_system.to_si(M::effective_Kh, scon[VI::SConn::EffectiveKH])), - segdist_end( unit_system.to_si(M::length, scon[VI::SConn::SegDistEnd])), - segdist_start( unit_system.to_si(M::length, scon[VI::SConn::SegDistStart])), - oil_rate( unit_system.to_si(M::liquid_surface_rate, xcon[VI::XConn::OilRate])), - water_rate( unit_system.to_si(M::liquid_surface_rate, xcon[VI::XConn::WaterRate])), - gas_rate( unit_system.to_si(M::gas_surface_rate, xcon[VI::XConn::GasRate])), - pressure( unit_system.to_si(M::pressure, xcon[VI::XConn::Pressure])), - resv_rate( unit_system.to_si(M::rate, xcon[VI::XConn::ResVRate])), - r0( RstConnection::inverse_peaceman(this->cf, this->kh, this->diameter/2, this->skin_factor) ) - /* - r0: The r0 quantity is currently not written or read from the restart - file. If the r0 value is given explicitly in the deck it is possible - to give a value which is not consistent with the Peaceman formula - - that value will be lost when loading back from a restart file. - */ +Opm::RestartIO::RstConnection::RstConnection(const UnitSystem& unit_system, + const std::size_t rst_index_, + const int nsconz, + const int* icon, + const float* scon, + const double* xcon) + : rst_index { rst_index_ } + // ----------------------------------------------------------------- + // Integer values (ICON) + , ijk { icon[VI::IConn::CellI] - 1 , + icon[VI::IConn::CellJ] - 1 , + icon[VI::IConn::CellK] - 1 } + , state { from_int(icon[VI::IConn::ConnStat]) } + , drain_sat_table { icon[VI::IConn::Drainage] } + , imb_sat_table { icon[VI::IConn::Imbibition] } + , completion { icon[VI::IConn::ComplNum] } + , dir { from_int(icon[VI::IConn::ConnDir]) } + , segment { icon[VI::IConn::Segment] } + // ----------------------------------------------------------------- + // Float values (SCON) + , cf_kind { from_float(1.0f) } + , skin_factor { scon[VI::SConn::SkinFactor] } + , cf { as_float(unit_system.to_si(M::transmissibility, scon[VI::SConn::ConnTrans])) } + , depth { as_float(unit_system.to_si(M::length, scon[VI::SConn::Depth])) } + , diameter { as_float(unit_system.to_si(M::length, scon[VI::SConn::Diameter])) } + , kh { as_float(unit_system.to_si(M::effective_Kh, scon[VI::SConn::EffectiveKH])) } + , denom { scon[VI::SConn::CFDenom] } + , length { as_float(unit_system.to_si(M::length, scon[VI::SConn::EffectiveLength])) } + , static_dfac_corr_coeff { staticDFactorCorrCoeff(unit_system, scon[VI::SConn::StaticDFacCorrCoeff]) } + , segdist_end { as_float(unit_system.to_si(M::length, scon[VI::SConn::SegDistEnd])) } + , segdist_start { as_float(unit_system.to_si(M::length, scon[VI::SConn::SegDistStart])) } + // ----------------------------------------------------------------- + // Double values (XCON) + , oil_rate { unit_system.to_si(M::liquid_surface_rate, xcon[VI::XConn::OilRate]) } + , water_rate { unit_system.to_si(M::liquid_surface_rate, xcon[VI::XConn::WaterRate]) } + , gas_rate { unit_system.to_si(M::gas_surface_rate, xcon[VI::XConn::GasRate]) } + , pressure { unit_system.to_si(M::pressure, xcon[VI::XConn::Pressure]) } + , resv_rate { unit_system.to_si(M::rate, xcon[VI::XConn::ResVRate]) } + // ----------------------------------------------------------------- + // Derived quantities + , r0 { pressEquivRadius(this->denom, this->skin_factor, this->diameter / 2) } { - if (static_cast(nsconz) > VI::SConn::CFInDeck) + if (static_cast(nsconz) > VI::SConn::CFInDeck) { this->cf_kind = from_float(scon[VI::SConn::CFInDeck]); -} - -} + } } diff --git a/src/opm/io/eclipse/rst/state.cpp b/src/opm/io/eclipse/rst/state.cpp index 27549704c..f767764f6 100644 --- a/src/opm/io/eclipse/rst/state.cpp +++ b/src/opm/io/eclipse/rst/state.cpp @@ -3,7 +3,8 @@ This file is part of the Open Porous Media project (OPM). - OPM is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. diff --git a/src/opm/io/eclipse/rst/well.cpp b/src/opm/io/eclipse/rst/well.cpp index bf9c86403..89406003f 100644 --- a/src/opm/io/eclipse/rst/well.cpp +++ b/src/opm/io/eclipse/rst/well.cpp @@ -17,17 +17,24 @@ along with OPM. If not, see . */ -#include +#include #include #include -#include #include #include #include + #include +#include + +#include +#include +#include +#include + #include #include #include @@ -54,28 +61,34 @@ namespace { { return is_sentinel(raw_value) ? raw_value : convert(raw_value); } -} + + float dfactor_correlation_coefficient_a(const Opm::UnitSystem& unit_system, + const float coeff_a) + { + const auto dimension = Opm::ParserKeywords::WDFACCOR{} + .getRecord(0).get(Opm::ParserKeywords::WDFACCOR::A::itemName) + .dimensions().front(); + + return static_cast(unit_system.to_si(dimension, coeff_a)); + } + + constexpr int def_ecl_phase = 1; + constexpr int def_pvt_table = 0; +} // Anonymous namespace namespace VI = ::Opm::RestartIO::Helpers::VectorItems; +using M = ::Opm::UnitSystem::measure; -namespace Opm { -namespace RestartIO { - -constexpr int def_ecl_phase = 1; -constexpr int def_pvt_table = 0; - -using M = ::Opm::UnitSystem::measure; - -RstWell::RstWell(const ::Opm::UnitSystem& unit_system, - const RstHeader& header, - const std::string& group_arg, - const std::string* zwel, - const int * iwel, - const float * swel, - const double * xwel, - const int * icon, - const float * scon, - const double * xcon) : +Opm::RestartIO::RstWell::RstWell(const UnitSystem& unit_system, + const RstHeader& header, + const std::string& group_arg, + const std::string* zwel, + const int* iwel, + const float* swel, + const double* xwel, + const int* icon, + const float* scon, + const double* xcon) : name(rtrim_copy(zwel[0])), group(group_arg), ij( {iwel[VI::IWell::IHead] - 1, iwel[VI::IWell::JHead] - 1}), @@ -139,6 +152,9 @@ RstWell::RstWell(const ::Opm::UnitSystem& unit_system, glift_min_rate( unit_system.to_si(M::gas_surface_rate, swel[VI::SWell::LOminRate])), glift_weight_factor( swel[VI::SWell::LOweightFac]), glift_inc_weight_factor( swel[VI::SWell::LOincFac]), + dfac_corr_coeff_a(dfactor_correlation_coefficient_a(unit_system, swel[VI::SWell::DFacCorrCoeffA])), + dfac_corr_exponent_b( swel[VI::SWell::DFacCorrExpB]), + dfac_corr_exponent_c( swel[VI::SWell::DFacCorrExpC]), // oil_rate( unit_system.to_si(M::liquid_surface_rate, xwel[VI::XWell::OilPrRate])), water_rate( unit_system.to_si(M::liquid_surface_rate, xwel[VI::XWell::WatPrRate])), @@ -167,66 +183,81 @@ RstWell::RstWell(const ::Opm::UnitSystem& unit_system, gas_void_rate( unit_system.to_si(M::gas_surface_volume, xwel[VI::XWell::GasVoidPrRate])) { - for (std::size_t tracer_index = 0; tracer_index < static_cast(header.runspec.tracers().water_tracers()); tracer_index++) - this->tracer_concentration_injection.push_back( swel[VI::SWell::TracerOffset + tracer_index] ); + for (std::size_t tracer_index = 0; + tracer_index < static_cast(header.runspec.tracers().water_tracers()); + ++tracer_index) + { + this->tracer_concentration_injection.push_back(swel[VI::SWell::TracerOffset + tracer_index]); + } + for (int ic = 0; ic < iwel[VI::IWell::NConn]; ++ic) { + const std::size_t icon_offset = ic * header.niconz; + const std::size_t scon_offset = ic * header.nsconz; + const std::size_t xcon_offset = ic * header.nxconz; - for (int ic = 0; ic < iwel[VI::IWell::NConn]; ic++) { - std::size_t icon_offset = ic * header.niconz; - std::size_t scon_offset = ic * header.nsconz; - std::size_t xcon_offset = ic * header.nxconz; - this->connections.emplace_back( unit_system, ic, header.nsconz, icon + icon_offset, scon + scon_offset, xcon + xcon_offset); + this->connections.emplace_back(unit_system, ic, header.nsconz, + icon + icon_offset, + scon + scon_offset, + xcon + xcon_offset); } } -RstWell::RstWell(const ::Opm::UnitSystem& unit_system, - const RstHeader& header, - const std::string& group_arg, - const std::string* zwel, - const int * iwel, - const float * swel, - const double * xwel, - const int * icon, - const float * scon, - const double * xcon, - const std::vector& iseg, - const std::vector& rseg) : - RstWell(unit_system, header, group_arg, zwel, iwel, swel, xwel, icon, scon, xcon) +Opm::RestartIO::RstWell::RstWell(const UnitSystem& unit_system, + const RstHeader& header, + const std::string& group_arg, + const std::string* zwel, + const int* iwel, + const float* swel, + const double* xwel, + const int* icon, + const float* scon, + const double* xcon, + const std::vector& iseg, + const std::vector& rseg) + : RstWell { unit_system, header, group_arg, + zwel, iwel, swel, xwel, + icon, scon, xcon } { + if (this->msw_index == 0) { + return; + } - if (this->msw_index) { - std::unordered_map segment_map; - for (int is=0; is < header.nsegmx; is++) { - std::size_t iseg_offset = header.nisegz * (is + (this->msw_index - 1) * header.nsegmx); - std::size_t rseg_offset = header.nrsegz * (is + (this->msw_index - 1) * header.nsegmx); - auto other_segment_number = iseg[iseg_offset + VI::ISeg::SegNo]; - if (other_segment_number != 0) { - auto segment_number = is + 1; - segment_map.insert({segment_number, this->segments.size()}); - this->segments.emplace_back( unit_system, segment_number, iseg.data() + iseg_offset, rseg.data() + rseg_offset); - } + std::unordered_map segment_map; + for (int is = 0; is < header.nsegmx; ++is) { + const std::size_t iseg_offset = header.nisegz * (is + (this->msw_index - 1)*header.nsegmx); + const std::size_t rseg_offset = header.nrsegz * (is + (this->msw_index - 1)*header.nsegmx); + const auto other_segment_number = iseg[iseg_offset + VI::ISeg::SegNo]; + + if (other_segment_number != 0) { + auto segment_number = is + 1; + segment_map.insert({segment_number, this->segments.size()}); + this->segments.emplace_back(unit_system, segment_number, + iseg.data() + iseg_offset, + rseg.data() + rseg_offset); } + } - for (auto& segment : this->segments) { - if (segment.outlet_segment != 0) { - auto& outlet_segment = this->segments[ segment_map[segment.outlet_segment] ]; - outlet_segment.inflow_segments.push_back(segment.segment); - } + for (auto& segment : this->segments) { + if (segment.outlet_segment != 0) { + auto& outlet_segment = this->segments[ segment_map[segment.outlet_segment] ]; + outlet_segment.inflow_segments.push_back(segment.segment); } } } -const RstSegment& -RstWell::segment(int segment_number) const -{ - const auto& iter = std::find_if(this->segments.begin(), this->segments.end(), [segment_number](const RstSegment& segment) { return segment.segment == segment_number; }); - if (iter == this->segments.end()) +const Opm::RestartIO::RstSegment& +Opm::RestartIO::RstWell::segment(int segment_number) const +{ + auto iter = std::find_if(this->segments.begin(), this->segments.end(), + [segment_number](const RstSegment& segment) + { return segment.segment == segment_number; }); + + if (iter == this->segments.end()) { throw std::invalid_argument("No such segment"); + } return *iter; } - -}} // namepace Opm::RestartIO