Add productivity_index to SingleWellState

opm/simulators/wells/BlackoilWellModelGeneric.cpp

@@ -1060,7 +1060,7 @@ BlackoilWellModelGeneric::
 wellPI(const int well_index) const
 {
     const auto& pu = this->phase_usage_;
-    const auto& pi = this->wellState().productivityIndex(well_index);
+    const auto& pi = this->wellState().well(well_index).productivity_index;
 
     const auto preferred = this->wells_ecl_[well_index].getPreferredPhase();
     switch (preferred) { // Should really have LIQUID = OIL + WATER here too...

opm/simulators/wells/MultisegmentWell_impl.hpp

@@ -597,9 +597,10 @@ namespace Opm
             std::transform(src, src + np, dest, dest, std::plus<>{});
         };
 
+        auto& ws = well_state.well(this->index_of_well_);
         auto& perf_data = well_state.perfData(this->index_of_well_);
-        auto* wellPI = well_state.productivityIndex(this->index_of_well_).data();
         auto* connPI = perf_data.prod_index.data();
+        auto* wellPI = ws.productivity_index.data();
 
         setToZero(wellPI);
 

opm/simulators/wells/SingleWellState.cpp

@@ -25,6 +25,7 @@ SingleWellState::SingleWellState(bool is_producer, std::size_t num_phases, doubl
     : producer(is_producer)
     , temperature(temp)
     , well_potentials(num_phases)
+    , productivity_index(num_phases)
 {}
 
 

opm/simulators/wells/SingleWellState.hpp

@@ -39,6 +39,7 @@ public:
     double dissolved_gas_rate{0};
     double vaporized_oil_rate{0};
     std::vector<double> well_potentials;
+    std::vector<double> productivity_index;
     Events events;
     Well::InjectorCMode injection_cmode{Well::InjectorCMode::CMODE_UNDEFINED};
     Well::ProducerCMode production_cmode{Well::ProducerCMode::CMODE_UNDEFINED};

opm/simulators/wells/StandardWell_impl.hpp

@@ -1374,8 +1374,9 @@ namespace Opm
             std::transform(src, src + np, dest, dest, std::plus<>{});
         };
 
+        auto& ws = well_state.well(this->index_of_well_);
         auto& perf_data = well_state.perfData(this->index_of_well_);
-        auto* wellPI = well_state.productivityIndex(this->index_of_well_).data();
+        auto* wellPI = ws.productivity_index.data();
         auto* connPI = perf_data.prod_index.data();
 
         setToZero(wellPI);

opm/simulators/wells/WellState.cpp

@@ -44,7 +44,6 @@ void WellState::base_init(const std::vector<double>& cellPressures,
     this->parallel_well_info_.clear();
     this->wellrates_.clear();
     this->segment_state.clear();
-    this->productivity_index_.clear();
     this->wells_.clear();
     {
         // const int nw = wells->number_of_wells;
@@ -95,7 +94,6 @@ void WellState::initSingleWell(const std::vector<double>& cellPressures,
     const int num_perf_this_well = well_info->communication().sum(well_perf_data.size());
     this->segment_state.add(well.name(), SegmentState{});
     this->perfdata.add(well.name(), PerfData{static_cast<std::size_t>(num_perf_this_well), well.isInjector(), this->phase_usage_});
-    this->productivity_index_.add(well.name(), std::vector<double>(np, 0));
 
     if ( num_perf_this_well == 0 )
         return;
@@ -391,7 +389,7 @@ void WellState::init(const std::vector<double>& cellPressures,
                 }
 
                 // Productivity index.
-                this->productivity_index_.copy_welldata( prevState->productivity_index_, wname );
+                new_well.productivity_index = prev_well.productivity_index;
             }
 
             // If in the new step, there is no THP related
@@ -481,7 +479,7 @@ WellState::report(const int* globalCellIdxMap,
 
         const auto& reservoir_rates = this->well_reservoir_rates_[well_index];
         const auto& well_potentials = ws.well_potentials;
-        const auto& wpi = this->productivity_index_[well_index];
+        const auto& wpi = ws.productivity_index;
         const auto& wv = this->wellRates(well_index);
 
         data::Well well;
@@ -822,7 +820,7 @@ void WellState::shutWell(int well_index)
     this->wellrates_[well_index].assign(np, 0);
 
     auto& resv = this->well_reservoir_rates_[well_index];
-    auto& wpi  = this->productivity_index_[well_index];
+    auto& wpi  = ws.productivity_index;
 
     for (int p = 0; p < np; ++p) {
         resv[p] = 0.0;

opm/simulators/wells/WellState.hpp

@@ -178,14 +178,6 @@ public:
         return this->segment_state[wname];
     }
 
-    std::vector<double>& productivityIndex(std::size_t well_index) {
-        return this->productivity_index_[well_index];
-    }
-
-    const std::vector<double>& productivityIndex(std::size_t well_index) const {
-        return this->productivity_index_[well_index];
-    }
-
     template<class Comm>
     void communicateGroupRates(const Comm& comm);
 
@@ -349,9 +341,6 @@ private:
 
     WellContainer<SegmentState> segment_state;
 
-    // Productivity Index
-    WellContainer<std::vector<double>> productivity_index_;
-
     data::Segment
     reportSegmentResults(const PhaseUsage& pu,
                          const int         well_id,