diff --git a/src/opm/output/eclipse/AggregateMSWData.cpp b/src/opm/output/eclipse/AggregateMSWData.cpp index 4b581a4d7..61a1cfcf2 100644 --- a/src/opm/output/eclipse/AggregateMSWData.cpp +++ b/src/opm/output/eclipse/AggregateMSWData.cpp @@ -925,9 +925,9 @@ namespace { /// 1, 2, 3, 4, 5, 6 -- Branch (1) /// 11, 12, 13, 14, 15, 16 -- Branch (2) /// 7, 8, 9, 10 -- Branch (3) + /// 17, 18, 19 -- Branch (4) /// 20, 22, 23, 24 -- Branch (5) /// 21, -- Branch (6) - /// 17, 18, 19 -- Branch (4) /// /// +------------------------------------------------------------+ /// | | @@ -1007,6 +1007,49 @@ namespace { void traverseStructure(); private: + /// Representation of a branch kick-off point + struct KickOffPoint + { + /// Branch start segment + int segment; + + /// Branch outlet segment. Segment from which the branch + /// kicks off. + int outlet; + + /// ID of branch starting at this kick-off point. + int branch; + }; + + /// Priority queue ordering operation + struct KickOffPointPriorityOrder + { + /// Priority comparison operation + /// + /// \param[in] p1 Left-hand side of comparison + /// \param[in] p2 Right-hand side of comparison + /// + /// \return Whether \p p1 has a lower priority than \p p2. + bool operator()(const KickOffPoint& p1, + const KickOffPoint& p2) const + { + // We use '>' here to have low outlet segment IDs move + // to the front of the priority queue. If multiple + // branches kick off from the same outlet segment, we + // sort these by their branch IDs whence lower-ID + // branches are created before higher-ID branches. + return std::tie(p1.outlet, p1.branch) + > std::tie(p2.outlet, p2.branch); + } + }; + + /// Priority queue container type for branch kick-off points. + using KickOffPointsQueue = std::priority_queue< + KickOffPoint, + std::vector, + KickOffPointPriorityOrder + >; + /// Name of well being explored. std::string_view well_{}; @@ -1019,12 +1062,18 @@ namespace { /// Callback routine for visiting a new branch. NewBranchCallback runNewBranchCallback_{}; - /// Segments from which to kick off searching new branches. - std::queue kickOffSegments_{}; + /// Points from which to kick off new branches. + KickOffPointsQueue kickOffPoints_{}; /// One-based segment number of currently visited segment. int currentSegment_{}; + /// Create a new branch from "top" kick-off point. + /// + /// Invokes new branch callback and makes the "current" segment + /// be the kick-off segment. + void createNewBranch(); + /// Find all segments on current branch, in order from heel to /// toe. /// @@ -1042,20 +1091,6 @@ namespace { void discoverNewBranches(const int outletSegment, const std::vector& children); - /// Enqueue new branch. - /// - /// Will be visited later. Invokes new branch callback. - /// - /// \param[in] branchId Branch number for new branch. - /// - /// \param[in] kickOffSegment First segment on new branch. - /// - /// \param[in] outletSegment Segment on branch from which the - /// new branch kicks off. - void discoverNewBranch(const int branchId, - const int kickOffSegment, - const int outletSegment); - /// Split child segments of current segment into groups /// based on their associate branch number. /// @@ -1083,16 +1118,33 @@ namespace { void Topology::traverseStructure() { - this->kickOffSegments_.push(1); - - while (! this->kickOffSegments_.empty()) { - this->currentSegment_ = this->kickOffSegments_.front(); - this->kickOffSegments_.pop(); + // Search starts at top segment: branch 1, segment 1. + this->kickOffPoints_.push(KickOffPoint { 1, 0, 1 }); + while (! this->kickOffPoints_.empty()) { + this->createNewBranch(); this->buildCurrentBranch(); } } + void Topology::createNewBranch() + { + const auto kickOff = this->kickOffPoints_.top(); + this->kickOffPoints_.pop(); + + if (kickOff.branch != 1) { + // All branches other than the main branch (branch ID 1) + // create new branch objects. Run new branch callback for + // these. + this->runNewBranchCallback_(this->well_, + kickOff.branch, + kickOff.segment, + kickOff.outlet); + } + + this->currentSegment_ = kickOff.segment; + } + void Topology::buildCurrentBranch() { while (true) { @@ -1121,23 +1173,14 @@ namespace { const std::vector& children) { for (const auto& child : children) { - const auto branch = this->segment(child).branchNumber(); - this->discoverNewBranch(branch, child, outletSegment); + this->kickOffPoints_ + .push(KickOffPoint { + child, outletSegment, + this->segment(child).branchNumber() + }); } } - void Topology::discoverNewBranch(const int branchId, - const int kickOffSegment, - const int outletSegment) - { - this->runNewBranchCallback_(this->well_, - branchId, - kickOffSegment, - outletSegment); - - this->kickOffSegments_.push(kickOffSegment); - } - std::pair, std::optional> Topology::characteriseChildSegments() const { diff --git a/tests/test_AggregateMSWData.cpp b/tests/test_AggregateMSWData.cpp index 0fcc69853..aa167a76f 100644 --- a/tests/test_AggregateMSWData.cpp +++ b/tests/test_AggregateMSWData.cpp @@ -243,6 +243,113 @@ END )"); } +//------------------------------------------------------------------+ +// Multi-lateral well with several ICDs and/or valves in a segment | +// structure of the form | +//------------------------------------------------------------------+ +// | +// 8 9 10 11 12 | +// o-----o-----o------o------o------o (2) | +// 7 / \ | +// / \ | +// / 19 \ 20 26 | +// 1 2 3 / 4 5 6 o----o----o (4) | +// ---o---o---o-----o---o---o (1) | +// \ | +// 13 \ (7) (6) | +// \ o o | +// o 25 / 24 / | +// \ / / | +// 14 \ 15 / 16 17 / 18 | +// o-----o--------o-------o--------o (3) | +// \ | +// 21 \ 22 23 | +// o-------o--------o (5) | +// \ \ \ | +// 27 \ 28 \ 29 \ | +// o o o | +// (8) (9) (10) | +// | +//------------------------------------------------------------------+ +// Branch ( 1): 1, 2, 3, 4, 5, 6 | +// Branch ( 2): 7, 8, 9, 10, 11, 12 | +// Branch ( 3): 13, 14, 15, 16, 17, 18 | +// Branch ( 4): 19, 20, 26 | +// Branch ( 5): 21, 22, 23 | +// Branch ( 6): 24 | +// Branch ( 7): 25 | +// Branch ( 8): 27 | +// Branch ( 9): 24 | +// Branch (10): 29 | +//------------------------------------------------------------------+ +Opm::Deck multilaterals_with_icd_valve() +{ + return Opm::Parser{}.parseString(R"(RUNSPEC +START +23 'NOV' 2023 / +DIMENS +10 10 3 / +OIL +GAS +WATER +DISGAS +VAPOIL +GRID +DXV +10*100.0 / +DYV +10*100.0 / +DZV +3*5.0 / +PERMX +300*100.0 / +COPY +PERMX PERMY / +PERMX PERMZ / +/ +MULTIPLY +PERMZ 0.1 / +/ +PORO +300*0.3 / +DEPTHZ +121*2000.0 / +SCHEDULE +WELSPECS + 'MLP' 'G' 10 10 2002.5 'OIL' / +/ +COMPDAT + 'MLP' 10 10 3 3 'OPEN' 1* 123.4 / +/ +WELSEGS + 'MLP' 2002.5 0.0 1* 'INC' 'H--' / +-- + 2 6 1 1 0.1 0.1 0.2 0.01 / + 7 12 2 3 0.1 0.1 0.2 0.01 / + 13 18 3 2 0.1 0.1 0.2 0.01 / + 19 20 4 10 0.1 0.1 0.2 0.01 / + 26 26 4 20 0.1 0.1 0.2 0.01 / + 21 23 5 15 0.1 0.1 0.2 0.01 / + 24 24 6 17 0.1 0.1 0.2 0.01 / + 25 25 7 15 0.1 0.1 0.2 0.01 / + 27 27 8 21 0.1 0.1 0.2 0.01 / + 28 28 9 22 0.1 0.1 0.2 0.01 / + 29 29 10 23 0.1 0.1 0.2 0.01 / +/ +COMPSEGS + 'MLP' / +-- + 10 10 3 5 0.0 1.0 'Z' / +/ +WCONPROD + 'MLP' 'OPEN' 'ORAT' 321.0 4* 10.0 / +/ +TSTEP +5*30 / +END +)"); +} + } // Anonymous namespace struct SimulationCase @@ -503,9 +610,9 @@ BOOST_AUTO_TEST_CASE (Declared_MSW_Data) // 1, 2, 3, 4, 5, 6 -- Branch (1) // 11, 12, 13, 14, 15, 16 -- Branch (2) // 7, 8, 9, 10 -- Branch (3) +// 17, 18, 19 -- Branch (4) // 20, 22, 23, 24 -- Branch (5) // 21, -- Branch (6) -// 17, 18, 19 -- Branch (4) // BOOST_AUTO_TEST_CASE(Multilateral_Branches) { @@ -540,7 +647,7 @@ BOOST_AUTO_TEST_CASE(Multilateral_Branches) BOOST_CHECK_EQUAL(ilbs.size(), std::vector::size_type{6}); const auto expect = std::vector { - 11, 7, 20, 21, 17, 0, + 11, 7, 17, 20, 21, 0, }; BOOST_CHECK_EQUAL_COLLECTIONS(ilbs .begin(), ilbs .end(), @@ -593,7 +700,7 @@ BOOST_AUTO_TEST_CASE(Multilateral_Branches) BOOST_CHECK_EQUAL(ilbr[VI::ILbr::NumBranchSegments], 3); BOOST_CHECK_EQUAL(ilbr[VI::ILbr::FirstSegment], 17); BOOST_CHECK_EQUAL(ilbr[VI::ILbr::LastSegment], 19); - BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 5); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 3); } // ILBR, branch 5 @@ -604,7 +711,7 @@ BOOST_AUTO_TEST_CASE(Multilateral_Branches) BOOST_CHECK_EQUAL(ilbr[VI::ILbr::NumBranchSegments], 4); BOOST_CHECK_EQUAL(ilbr[VI::ILbr::FirstSegment], 20); BOOST_CHECK_EQUAL(ilbr[VI::ILbr::LastSegment], 24); - BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 3); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 4); } // ILBR, branch 6 @@ -615,7 +722,7 @@ BOOST_AUTO_TEST_CASE(Multilateral_Branches) BOOST_CHECK_EQUAL(ilbr[VI::ILbr::NumBranchSegments], 1); BOOST_CHECK_EQUAL(ilbr[VI::ILbr::FirstSegment], 21); BOOST_CHECK_EQUAL(ilbr[VI::ILbr::LastSegment], 21); - BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 4); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 5); } } @@ -680,6 +787,247 @@ BOOST_AUTO_TEST_CASE(Multilateral_Segments_ISEG_0) } } +// The segments and branches must appear in the following order in the +// ILBS/ILBR output arrays. +// +// 1, 2, 3, 4, 5, 6 -- Branch ( 1) +// 13, 14, 15, 16, 17, 18 -- Branch ( 3) +// 7, 8, 9, 10, 11, 12 -- Branch ( 2) +// 19, 20, 26 -- Branch ( 4) +// 21, 22, 23 -- Branch ( 5) +// 25, -- Branch ( 7) +// 24, -- Branch ( 6) +// 27, -- Branch ( 8) +// 28, -- Branch ( 9) +// 29, -- Branch (10) +// +BOOST_AUTO_TEST_CASE(Multilateral_Branches_ICD_Valve) +{ + const auto cse = SimulationCase { multilaterals_with_icd_valve() }; + + const auto& es = cse.es; + const auto& grid = cse.grid; + const auto& sched = cse.sched; + + // Report Step 1: 2023-11-23 --> 2023-12-23 + const auto rptStep = std::size_t {1}; + + const auto secs_elapsed = 30 * 86'400.0; + const auto ih = Opm::RestartIO::Helpers:: + createInteHead(es, grid, sched, secs_elapsed, + rptStep, rptStep + 1, rptStep); + + // No dynamic data needed for this test. We're checking the structure + // only. + const auto smry = Opm::SummaryState { Opm::TimeService::now() }; + const auto xw = Opm::data::Wells {}; + + auto amswd = Opm::RestartIO::Helpers::AggregateMSWData {ih}; + amswd.captureDeclaredMSWData(sched, rptStep, es.getUnits(), + ih, grid, smry, xw); + + // ILBS--First segment on each branch other than branch 1. Ordered by + // discovery. + { + const auto& ilbs = amswd.getILBs(); + + // No WSEGDIMS => size = maximum branch number + BOOST_CHECK_EQUAL(ilbs.size(), std::vector::size_type{10}); + + const auto expect = std::vector { + 13, 7, 19, 21, 25, 24, 27, 28, 29, 0, + }; + + BOOST_CHECK_EQUAL_COLLECTIONS(ilbs .begin(), ilbs .end(), + expect.begin(), expect.end()); + } + + auto ilbrOffset = [&ih](const int branch) + { + return ih[VI::intehead::NILBRZ] * (branch - 1); + }; + + // ILBR, branch 1 + { + const auto* ilbr = &amswd.getILBr()[ilbrOffset(1)]; + + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::OutletSegment], 0); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::NumBranchSegments], 6); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::FirstSegment], 1); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::LastSegment], 6); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 0); + } + + // ILBR, branch 2 + { + const auto* ilbr = &amswd.getILBr()[ilbrOffset(2)]; + + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::OutletSegment], 3); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::NumBranchSegments], 6); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::FirstSegment], 7); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::LastSegment], 12); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 2); + } + + // ILBR, branch 3 + { + const auto* ilbr = &amswd.getILBr()[ilbrOffset(3)]; + + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::OutletSegment], 2); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::NumBranchSegments], 6); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::FirstSegment], 13); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::LastSegment], 18); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 1); + } + + // ILBR, branch 4 + { + const auto* ilbr = &amswd.getILBr()[ilbrOffset(4)]; + + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::OutletSegment], 10); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::NumBranchSegments], 3); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::FirstSegment], 19); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::LastSegment], 26); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 3); + } + + // ILBR, branch 5 + { + const auto* ilbr = &amswd.getILBr()[ilbrOffset(5)]; + + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::OutletSegment], 15); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::NumBranchSegments], 3); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::FirstSegment], 21); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::LastSegment], 23); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 4); + } + + // ILBR, branch 6 + { + const auto* ilbr = &amswd.getILBr()[ilbrOffset(6)]; + + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::OutletSegment], 17); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::NumBranchSegments], 1); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::FirstSegment], 24); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::LastSegment], 24); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 6); + } + + // ILBR, branch 7 + { + const auto* ilbr = &amswd.getILBr()[ilbrOffset(7)]; + + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::OutletSegment], 15); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::NumBranchSegments], 1); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::FirstSegment], 25); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::LastSegment], 25); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 5); + } + + // ILBR, branch 8 + { + const auto* ilbr = &amswd.getILBr()[ilbrOffset(8)]; + + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::OutletSegment], 21); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::NumBranchSegments], 1); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::FirstSegment], 27); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::LastSegment], 27); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 7); + } + + // ILBR, branch 9 + { + const auto* ilbr = &amswd.getILBr()[ilbrOffset(9)]; + + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::OutletSegment], 22); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::NumBranchSegments], 1); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::FirstSegment], 28); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::LastSegment], 28); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 8); + } + + // ILBR, branch 10 + { + const auto* ilbr = &amswd.getILBr()[ilbrOffset(10)]; + + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::OutletSegment], 23); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::NumBranchSegments], 1); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::FirstSegment], 29); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::LastSegment], 29); + BOOST_CHECK_EQUAL(ilbr[VI::ILbr::KickOffDiscoveryOffset], 9); + } +} + +// The segments must appear in the following depth first search toe-to-heel +// order in ISEG[0]. We furthermore, go along kick-off branches before +// searching the main branch. Note that this order is *different* from +// ILBS/ILBR. +// +// 27 -- Branch ( 8) +// 28 -- Branch ( 9) +// 29 -- Branch (10) +// 23, 22, 21 -- Branch ( 5) +// 25 -- Branch ( 7) +// 24 -- Branch ( 6) +// 18, 17, 16, 15, 14, 13 -- Branch ( 3) +// 26, 20, 19 -- Branch ( 4) +// 12, 11, 10, 9, 8, 7 -- Branch ( 2) +// 6, 5, 4, 3, 2, 1 -- Branch ( 1) +// +BOOST_AUTO_TEST_CASE(Multilateral_ICD_Valve_ISEG_0) +{ + const auto cse = SimulationCase { multilaterals_with_icd_valve() }; + + const auto& es = cse.es; + const auto& grid = cse.grid; + const auto& sched = cse.sched; + + // Report Step 1: 2023-11-23 --> 2023-12-23 + const auto rptStep = std::size_t {1}; + + const auto secs_elapsed = 30 * 86'400.0; + const auto ih = Opm::RestartIO::Helpers:: + createInteHead(es, grid, sched, secs_elapsed, + rptStep, rptStep + 1, rptStep); + + // No dynamic data needed for this test. We're checking the structure + // only. + const auto smry = Opm::SummaryState { Opm::TimeService::now() }; + const auto xw = Opm::data::Wells {}; + + auto amswd = Opm::RestartIO::Helpers::AggregateMSWData {ih}; + amswd.captureDeclaredMSWData(sched, rptStep, es.getUnits(), + ih, grid, smry, xw); + + auto isegOffset = [&ih](const int ix) + { + return ih[VI::intehead::NISEGZ] * ix; + }; + + const auto expect = std::vector { + 27, // Branch ( 8) + 28, // Branch ( 9) + 29, // Branch (10) + 23, 22, 21, // Branch ( 5) + 25, // Branch ( 7) + 24, // Branch ( 6) + 18, 17, 16, 15, 14, 13, // Branch ( 3) + 26, 20, 19, // Branch ( 4) + 12, 11, 10, 9, 8, 7, // Branch ( 2) + 6, 5, 4, 3, 2, 1, // Branch ( 1) + }; + + const auto& iseg = amswd.getISeg(); + + for (auto i = 0*expect.size(); i < expect.size(); ++i) { + BOOST_CHECK_MESSAGE(iseg[isegOffset(i)] == expect[i], + "ISEG[0](" << i << ") == " + << iseg[isegOffset(i)] + << " differs from expected value " + << expect[i]); + } +} + BOOST_AUTO_TEST_CASE(MSW_AICD) { const auto simCase = SimulationCase {first_sim("TEST_AGGREGATE_MSW.DATA")};