Merge pull request #3296 from joakim-hove/segment-use-well-index

Segment use well index

opm/simulators/wells/BlackoilWellModel_impl.hpp

@@ -2011,14 +2011,13 @@ namespace Opm {
 
                 const WellSegments& segment_set = well_ecl.getSegments();
 
-                const int top_segment_index = well_state.topSegmentIndex(well_index);
                 const auto& rst_segments = rst_well.segments;
 
                 // \Note: eventually we need to hanlde the situations that some segments are shut
                 assert(0u + segment_set.size() == rst_segments.size());
 
-                auto * segment_pressure = &well_state.segPress()[top_segment_index];
+                auto segment_pressure = well_state.segPress(well_index);
-                auto * segment_rates  = &well_state.segRates()[top_segment_index*np];
+                auto segment_rates  = well_state.segRates(well_index);
                 for (const auto& rst_segment : rst_segments) {
                     const int segment_index = segment_set.segmentNumberToIndex(rst_segment.first);
 

opm/simulators/wells/MultisegmentWell_impl.hpp

@@ -278,8 +278,7 @@ namespace Opm
     MultisegmentWell<TypeTag>::
     scaleSegmentRatesWithWellRates(WellState& well_state) const
     {
-        const int top_segment_index = well_state.topSegmentIndex(index_of_well_);
+        auto segment_rates = well_state.segRates(index_of_well_);
-        auto * segment_rates = &well_state.segRates()[top_segment_index*this->number_of_phases_];
         for (int phase = 0; phase < number_of_phases_; ++phase) {
             const double unscaled_top_seg_rate = segment_rates[phase];
             const double well_phase_rate = well_state.wellRates(index_of_well_)[phase];
@@ -316,9 +315,8 @@ namespace Opm
     scaleSegmentPressuresWithBhp(WellState& well_state) const
     {
         //scale segment pressures
-        const int top_segment_index = well_state.topSegmentIndex(index_of_well_);
         const double bhp = well_state.bhp(index_of_well_);
-        auto * segment_pressure = &well_state.segPress()[top_segment_index];
+        auto segment_pressure = well_state.segPress(index_of_well_);
         const double unscaled_top_seg_pressure = segment_pressure[0];
         for (int seg = 0; seg < numberOfSegments(); ++seg) {
             segment_pressure[seg] *= bhp/unscaled_top_seg_pressure;
@@ -722,9 +720,8 @@ namespace Opm
         const Well& well = Base::wellEcl();
 
         // the index of the top segment in the WellState
-        const int top_segment_index = well_state.topSegmentIndex(index_of_well_);
+        const auto segment_rates = well_state.segRates(index_of_well_);
-        const auto * segment_rates = &well_state.segRates()[top_segment_index * this->number_of_phases_];
+        const auto segment_pressure = well_state.segPress(index_of_well_);
-        const auto * segment_pressure = &well_state.segPress()[top_segment_index];
         const PhaseUsage& pu = phaseUsage();
 
         for (int seg = 0; seg < numberOfSegments(); ++seg) {
@@ -2340,8 +2337,8 @@ namespace Opm
         assert( FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) );
         const int oil_pos = pu.phase_pos[Oil];
 
-        auto * segment_rates = &well_state.segRates()[well_state.topSegmentIndex(this->index_of_well_) * this->number_of_phases_];
+        auto segment_rates = well_state.segRates(this->index_of_well_);
-        auto * segment_pressure = &well_state.segPress()[well_state.topSegmentIndex(this->index_of_well_)];
+        auto segment_pressure = well_state.segPress(this->index_of_well_);
         for (int seg = 0; seg < numberOfSegments(); ++seg) {
             std::vector<double> fractions(number_of_phases_, 0.0);
             fractions[oil_pos] = 1.0;

opm/simulators/wells/WellContainer.hpp

@@ -163,7 +163,7 @@ public:
         if (index_iter != this->index_map.end())
             return index_iter->second;
 
-        return {};
+        return std::nullopt;
     }
 
 

opm/simulators/wells/WellState.cpp

@@ -965,7 +965,6 @@ void WellState::initWellStateMSWell(const std::vector<Well>& wells_ecl,
                     continue;
                 }
 
-                const int old_top_segment_index = prev_well_state->topSegmentIndex(old_index_well);
                 const int new_top_segment_index = topSegmentIndex(new_index_well);
                 int number_of_segment = 0;
                 // if it is the last well in list
@@ -975,11 +974,11 @@ void WellState::initWellStateMSWell(const std::vector<Well>& wells_ecl,
                     number_of_segment = topSegmentIndex(new_index_well + 1) - new_top_segment_index;
                 }
 
-                auto * segment_rates = &this->seg_rates_[new_top_segment_index*np];
+                auto segment_rates = this->segRates(w);
-                auto * segment_pressure = &this->seg_press_[new_top_segment_index];
+                auto segment_pressure = this->segPress(w);
 
-                const auto * prev_segment_rates = &prev_well_state->segRates()[old_top_segment_index*np];
+                const auto prev_segment_rates = prev_well_state->segRates(old_index_well);
-                const auto * prev_segment_pressure = &prev_well_state->segPress()[new_top_segment_index];
+                const auto prev_segment_pressure = prev_well_state->segPress(old_index_well);
 
                 for (int seg=0; seg < number_of_segment; ++seg) {
                     for (int p = 0; p < np; ++p)
@@ -1152,23 +1151,20 @@ WellState::reportSegmentResults(const PhaseUsage& pu,
                                                      const int         seg_no) const
 {
     auto seg_res = data::Segment{};
-
     const auto seg_dof =
             this->topSegmentIndex(well_id) + seg_ix;
-
-    const auto* rate =
-            &this->segRates()[seg_dof * this->numPhases()];
-
     {
         using Value = data::SegmentPressures::Value;
         auto& segpress = seg_res.pressures;
-        segpress[Value::Pressure] = this->segPress()[seg_dof];
+        segpress[Value::Pressure] = this->segPress(well_id)[seg_ix];
         segpress[Value::PDrop] = this->segPressDrop()[seg_dof];
         segpress[Value::PDropHydrostatic] = this->segPressDropHydroStatic()[seg_dof];
         segpress[Value::PDropFriction] = this->segPressDropFriction()[seg_dof];
         segpress[Value::PDropAccel] = this->segPressDropAcceleration()[seg_dof];
     }
 
+    const auto segment_rates = this->segRates(well_id);
+    const auto rate = &segment_rates[seg_ix * this->numPhases()];
     if (pu.phase_used[Water]) {
         seg_res.rates.set(data::Rates::opt::wat,
                           rate[pu.phase_pos[Water]]);

opm/simulators/wells/WellState.hpp

@@ -192,19 +192,29 @@ public:
         return well_vaporized_oil_rates_[well_index];
     }
 
-    const std::vector<double>& segRates() const
+
+    const double * segRates(std::size_t well_index) const
     {
-        return seg_rates_;
+        const int top_segment_index = this->top_segment_index_[well_index];
+        return &this->seg_rates_[top_segment_index * this->phase_usage_.num_phases];
     }
 
-    std::vector<double>& segRates()
+    double * segRates(std::size_t well_index)
     {
-        return seg_rates_;
+        const int top_segment_index = this->top_segment_index_[well_index];
+        return &this->seg_rates_[top_segment_index * this->phase_usage_.num_phases];
     }
 
-    const std::vector<double>& segPress() const
+    double * segPress(std::size_t well_index)
     {
-        return seg_press_;
+        const int top_segment_index = this->top_segment_index_[well_index];
+        return &seg_press_[top_segment_index];
+    }
+
+    const double * segPress(std::size_t well_index) const
+    {
+        const int top_segment_index = this->top_segment_index_[well_index];
+        return &seg_press_[top_segment_index];
     }
 
     std::vector<double>& segPressDrop()
@@ -247,11 +257,6 @@ public:
         return seg_pressdrop_acceleration_;
     }
 
-    std::vector<double>& segPress()
-    {
-        return seg_press_;
-    }
-
     int numSegment() const
     {
         return nseg_;

tests/test_wellstate.cpp

@@ -164,14 +164,12 @@ namespace {
     {
         const auto nWell = wells.size();
 
-        auto& segPress = wstate.segPress();
 
         for (auto wellID = 0*nWell; wellID < nWell; ++wellID) {
             const auto& well     = wells[wellID];
-            const auto  topSegIx = wstate.topSegmentIndex(wellID);
             const auto  pressTop = 100.0 * wellID;
 
-            auto* press = &segPress[topSegIx];
+            auto* press = wstate.segPress(wellID);
 
             press[0] = pressTop;
 
@@ -208,16 +206,15 @@ namespace {
 
         const auto nWell = wells.size();
 
-        auto& segRates = wstate.segRates();
 
         for (auto wellID = 0*nWell; wellID < nWell; ++wellID) {
             const auto& well     = wells[wellID];
-            const auto  topSegIx = wstate.topSegmentIndex(wellID);
             const auto  rateTop  = 1000.0 * wellID;
+            auto segRates = wstate.segRates(wellID);
 
-            if (wat) { segRates[np*topSegIx + iw] = rateTop; }
+            if (wat) { segRates[iw] = rateTop; }
-            if (oil) { segRates[np*topSegIx + io] = rateTop; }
+            if (oil) { segRates[io] = rateTop; }
-            if (gas) { segRates[np*topSegIx + ig] = rateTop; }
+            if (gas) { segRates[ig] = rateTop; }
 
             if (! well.isMultiSegment()) {
                 continue;
@@ -230,7 +227,7 @@ namespace {
                 // One-based numbering scheme for segments.
                 const auto segNo = segSet[segID].segmentNumber();
 
-                auto* rates = &segRates[(topSegIx + segNo - 1) * np];
+                auto* rates = &segRates[(segNo - 1) * np];
 
                 if (wat) { rates[iw] = rateTop + 100.0*(segNo - 1); }
                 if (oil) { rates[io] = rateTop + 200.0*(segNo - 1); }