diff --git a/src/opm/output/eclipse/AggregateMSWData.cpp b/src/opm/output/eclipse/AggregateMSWData.cpp index 68e315228..e07a33d39 100644 --- a/src/opm/output/eclipse/AggregateMSWData.cpp +++ b/src/opm/output/eclipse/AggregateMSWData.cpp @@ -180,30 +180,37 @@ namespace { return segmentOrder(segSet, 0); } + /// Accumulate connection flow rates (surface conditions) to their connecting segment. Opm::RestartIO::Helpers::SegmentSetSourceSinkTerms - getSegmentSetSSTerms(const std::string& wname, const Opm::WellSegments& segSet, const std::vector& rateConns, - const Opm::WellConnections& welConns, const Opm::UnitSystem& units) + getSegmentSetSSTerms(const Opm::WellSegments& segSet, + const std::vector& rateConns, + const Opm::WellConnections& welConns, + const Opm::UnitSystem& units) { std::vector qosc (segSet.size(), 0.); std::vector qwsc (segSet.size(), 0.); std::vector qgsc (segSet.size(), 0.); - std::vector openConnections; - using M = ::Opm::UnitSystem::measure; - using R = ::Opm::data::Rates::opt; - for (auto nConn = welConns.size(), connID = 0*nConn; connID < nConn; connID++) { - if (welConns[connID].state() == Opm::Connection::State::OPEN) openConnections.push_back(&welConns[connID]); - } - if (openConnections.size() != rateConns.size()) { - throw std::invalid_argument { - "Inconsistent number of open connections I in Opm::WellConnections (" + - std::to_string(welConns.size()) + ") and vector (" + - std::to_string(rateConns.size()) + ") in Well " + wname - }; - } - for (auto nConn = openConnections.size(), connID = 0*nConn; connID < nConn; connID++) { - const auto& segNo = openConnections[connID]->segment(); - const auto& segInd = segSet.segmentNumberToIndex(segNo); - const auto& Q = rateConns[connID].rates; + + using M = ::Opm::UnitSystem::measure; + using R = ::Opm::data::Rates::opt; + + for (const auto& conn : welConns) { + if (conn.state() != Opm::Connection::State::OPEN) { + continue; + } + + auto xcPos = std::find_if(rateConns.begin(), rateConns.end(), + [&conn](const Opm::data::Connection& xconn) -> bool + { + return xconn.index == conn.global_index(); + }); + + if (xcPos == rateConns.end()) { + continue; + } + + const auto segInd = segSet.segmentNumberToIndex(conn.segment()); + const auto& Q = xcPos->rates; auto get = [&units, &Q](const M u, const R q) -> double { @@ -225,15 +232,17 @@ namespace { } Opm::RestartIO::Helpers::SegmentSetFlowRates - getSegmentSetFlowRates(const std::string& wname, const Opm::WellSegments& segSet, const std::vector& rateConns, - const Opm::WellConnections& welConns, const Opm::UnitSystem& units) + getSegmentSetFlowRates(const Opm::WellSegments& segSet, + const std::vector& rateConns, + const Opm::WellConnections& welConns, + const Opm::UnitSystem& units) { std::vector sofr (segSet.size(), 0.); std::vector swfr (segSet.size(), 0.); std::vector sgfr (segSet.size(), 0.); // //call function to calculate the individual segment source/sink terms - auto sSSST = getSegmentSetSSTerms(wname, segSet, rateConns, welConns, units); + auto segmentSources = getSegmentSetSSTerms(segSet, rateConns, welConns, units); // find an ordered list of segments auto orderedSegmentInd = segmentOrder(segSet); @@ -244,9 +253,9 @@ namespace { const auto& segInd = sNFOSN[indOSN]; // the segment flow rates is the sum of the the source/sink terms for each segment plus the flow rates from the inflow segments // add source sink terms - sofr[segInd] += sSSST.qosc[segInd]; - swfr[segInd] += sSSST.qwsc[segInd]; - sgfr[segInd] += sSSST.qgsc[segInd]; + sofr[segInd] += segmentSources.qosc[segInd]; + swfr[segInd] += segmentSources.qwsc[segInd]; + sgfr[segInd] += segmentSources.qgsc[segInd]; // add flow from all inflow segments for (const auto& segNo : segSet[segInd].inletSegments()) { const auto & ifSegInd = segSet.segmentNumberToIndex(segNo); @@ -255,6 +264,7 @@ namespace { sgfr[segInd] += sgfr[ifSegInd]; } } + return { sofr, swfr, @@ -655,7 +665,7 @@ namespace { // find well connections and calculate segment rates based on well connection production/injection terms auto sSFR = Opm::RestartIO::Helpers::SegmentSetFlowRates{}; if (haveWellRes) { - sSFR = getSegmentSetFlowRates(well.name(), welSegSet, wRatesIt->second.connections, welConns, units); + sSFR = getSegmentSetFlowRates(welSegSet, wRatesIt->second.connections, welConns, units); } std::string stringSegNum = std::to_string(segment0.segmentNumber()); diff --git a/tests/test_AggregateMSWData.cpp b/tests/test_AggregateMSWData.cpp index b495febda..c94cad37b 100644 --- a/tests/test_AggregateMSWData.cpp +++ b/tests/test_AggregateMSWData.cpp @@ -551,6 +551,7 @@ END double qo = 5.; double qw = 4.; double qg = 50.; + int firstConnectedCell = 90; // zero-based linear index of (1,5,2) for (int i = 0; i < 5; i++) { xw["PROD"].connections.emplace_back(); auto& c = xw["PROD"].connections.back(); @@ -558,6 +559,8 @@ END c.rates.set(o::wat, qw*(float(i)+1.)) .set(o::oil, qo*(float(i)+1.)) .set(o::gas, qg*(float(i)+1.)); + + c.index = firstConnectedCell + i; } auto seg = Opm::data::Segment{}; for (std::size_t i = 1; i < 5; i++) { @@ -569,6 +572,7 @@ END xw["WINJ"].rates.set(o::oil, 0.0); xw["WINJ"].rates.set(o::gas, 0.0); qw = 7.; + firstConnectedCell = 409; // zero-based linear index of (10,1,9) for (int i = 0; i < 5; i++) { xw["WINJ"].connections.emplace_back(); auto& c = xw["WINJ"].connections.back(); @@ -576,6 +580,8 @@ END c.rates.set(o::wat, qw*(float(i)+1.)) .set(o::oil, 0.) .set(o::gas, 0.); + + c.index = firstConnectedCell - i; } } return xw;