adding num_phases_ to MultisegmentWells

and removing all the np in the argument of the public methods.

opm/autodiff/BlackoilMultiSegmentModel_impl.hpp

@@ -456,7 +456,7 @@ namespace Opm {
             // Compute initial accumulation contributions
             // and well connection pressures.
             asImpl().computeAccum(state0, 0);
-            msWells().computeSegmentFluidProperties(state0, phaseCondition(), active_, fluid_, numPhases());
+            msWells().computeSegmentFluidProperties(state0, phaseCondition(), active_, fluid_);
             const int np = numPhases();
             assert(np == int(msWells().segmentCompSurfVolumeInitial().size()));
             for (int phase = 0; phase < np; ++phase) {
@@ -484,7 +484,7 @@ namespace Opm {
         }
 
         // asImpl().computeSegmentFluidProperties(state);
-        msWells().computeSegmentFluidProperties(state, phaseCondition(), active_, fluid_, numPhases());
+        msWells().computeSegmentFluidProperties(state, phaseCondition(), active_, fluid_);
 
         // asImpl().computeSegmentPressuresDelta(state);
         const double gravity = detail::getGravity(geo_.gravity(), UgGridHelpers::dimensions(grid_));
@@ -508,12 +508,12 @@ namespace Opm {
         const V perf_press_diffs = stdWells().wellPerforationPressureDiffs();
         msWells().computeWellFlux(state, fluid_.phaseUsage(), active_,
                                   perf_press_diffs, compi,
-                                  mob_perfcells, b_perfcells, np, aliveWells, cq_s);
+                                  mob_perfcells, b_perfcells, aliveWells, cq_s);
         asImpl().updatePerfPhaseRatesAndPressures(cq_s, state, well_state);
-        msWells().addWellFluxEq(cq_s, state, np, residual_);
+        msWells().addWellFluxEq(cq_s, state, residual_);
         asImpl().addWellContributionToMassBalanceEq(cq_s, state, well_state);
         // asImpl().addWellControlEq(state, well_state, aliveWells);
-        msWells().addWellControlEq(state, well_state, aliveWells, np, active_, residual_);
+        msWells().addWellControlEq(state, well_state, aliveWells, active_, residual_);
     }
 
 
@@ -616,7 +616,7 @@ namespace Opm {
     BlackoilMultiSegmentModel<Grid>::updateWellState(const V& dwells,
                                                      WellState& well_state)
     {
-        msWells().updateWellState(dwells, fluid_.numPhases(), dpMaxRel(), well_state);
+        msWells().updateWellState(dwells, dpMaxRel(), well_state);
     }
 
 
@@ -667,12 +667,12 @@ namespace Opm {
             const V perf_press_diffs = stdWells().wellPerforationPressureDiffs();
             msWells().computeWellFlux(wellSolutionState, fluid_.phaseUsage(), active_,
                                       perf_press_diffs, compi,
-                                      mob_perfcells_const, b_perfcells_const, np, aliveWells, cq_s);
+                                      mob_perfcells_const, b_perfcells_const, aliveWells, cq_s);
 
             updatePerfPhaseRatesAndPressures(cq_s, wellSolutionState, well_state);
-            msWells().addWellFluxEq(cq_s, wellSolutionState, np, residual_);
+            msWells().addWellFluxEq(cq_s, wellSolutionState, residual_);
             // addWellControlEq(wellSolutionState, well_state, aliveWells);
-            msWells().addWellControlEq(wellSolutionState, well_state, aliveWells, np, active_, residual_);
+            msWells().addWellControlEq(wellSolutionState, well_state, aliveWells, active_, residual_);
             converged = Base::getWellConvergence(it);
 
             if (converged) {

opm/autodiff/MultisegmentWells.cpp

@@ -143,11 +143,12 @@ namespace Opm {
     MultisegmentWells(const std::vector<WellMultiSegmentConstPtr>& wells_ms, const int np)
       : wells_multisegment_(wells_ms)
       , wops_ms_(wells_ms)
+      , num_phases_(np)
       , well_segment_perforation_pressure_diffs_(ADB::null())
       , well_segment_densities_(ADB::null())
       , well_segment_pressures_delta_(ADB::null())
-      , segment_comp_surf_volume_initial_(np)
+      , segment_comp_surf_volume_initial_(num_phases_)
-      , segment_comp_surf_volume_current_(np, ADB::null())
+      , segment_comp_surf_volume_current_(num_phases_, ADB::null())
       , segment_mass_flow_rates_(ADB::null())
       , segment_viscosities_(ADB::null())
     {

opm/autodiff/MultisegmentWells.hpp

@@ -85,7 +85,10 @@ namespace Opm {
             const std::vector<WellMultiSegmentConstPtr>& wells() const;
             const MultisegmentWellOps& wellOps() const;
 
+            int numPhases() const { return num_phases_; };
+
             int numSegment() const { return nseg_total_; };
+
             int numPerf() const { return nperf_total_; };
 
             const Vector& wellPerforationCellPressureDiffs() const { return well_perforation_cell_pressure_diffs_; };
@@ -129,7 +132,6 @@ namespace Opm {
             template <class WellState>
             void
             updateWellState(const Vector& dwells,
-                            const int np,
                             const double dpmaxrel,
                             WellState& well_state) const;
 
@@ -144,7 +146,6 @@ namespace Opm {
                             const DataBlock& compi,
                             const std::vector<ADB>& mob_perfcells,
                             const std::vector<ADB>& b_perfcells,
-                            const int np,
                             Vector& aliveWells,
                             std::vector<ADB>& cq_s) const;
 
@@ -156,8 +157,7 @@ namespace Opm {
             computeSegmentFluidProperties(const SolutionState& state,
                                           const std::vector<PhasePresence>& pc,
                                           const std::vector<bool>& active,
-                                          const BlackoilPropsAdInterface& fluid,
+                                          const BlackoilPropsAdInterface& fluid);
-                                          const int np);
 
             void
             computeSegmentPressuresDelta(const double grav);
@@ -166,7 +166,6 @@ namespace Opm {
             void
             addWellFluxEq(const std::vector<ADB>& cq_s,
                           const SolutionState& state,
-                          const int np,
                           LinearisedBlackoilResidual& residual);
 
             template <class SolutionState, class WellState>
@@ -174,7 +173,6 @@ namespace Opm {
             addWellControlEq(const SolutionState& state,
                              const WellState& xw,
                              const Vector& aliveWells,
-                             const int np,
                              const std::vector<bool>& active,
                              LinearisedBlackoilResidual& residual);
 
@@ -187,6 +185,7 @@ namespace Opm {
         // TODO: probably a wells_active_ will be required here.
         const std::vector<WellMultiSegmentConstPtr> wells_multisegment_;
         const MultisegmentWellOps wops_ms_;
+        const int num_phases_;
         int nseg_total_;
         int nperf_total_;
 

opm/autodiff/MultisegmentWells_impl.hpp

@@ -51,7 +51,6 @@ namespace Opm
     void
     MultisegmentWells::
     updateWellState(const Vector& dwells,
-                    const int np,
                     const double dpmaxrel,
                     WellState& well_state) const
     {
@@ -59,6 +58,7 @@ namespace Opm
         {
             const int nw = wells().size();
             const int nseg_total = nseg_total_;
+            const int np = numPhases();
 
             // Extract parts of dwells corresponding to each part.
             int varstart = 0;
@@ -128,12 +128,12 @@ namespace Opm
                     const DataBlock& compi,
                     const std::vector<ADB>& mob_perfcells,
                     const std::vector<ADB>& b_perfcells,
-                    const int np,
                     Vector& aliveWells,
                     std::vector<ADB>& cq_s) const
     {
         if (wells().size() == 0) return;
 
+        const int np = numPhases();
         const int nw = wells().size();
 
         aliveWells = Vector::Constant(nw, 1.0);
@@ -316,9 +316,9 @@ namespace Opm
     computeSegmentFluidProperties(const SolutionState& state,
                                   const std::vector<PhasePresence>& pc,
                                   const std::vector<bool>& active,
-                                  const BlackoilPropsAdInterface& fluid,
+                                  const BlackoilPropsAdInterface& fluid)
-                                  const int np)
     {
+        const int np = numPhases();
         const int nw = wells().size();
         const int nseg_total = nseg_total_;
 
@@ -524,7 +524,6 @@ namespace Opm
     MultisegmentWells::
     addWellFluxEq(const std::vector<ADB>& cq_s,
                   const SolutionState& state,
-                  const int np,
                   LinearisedBlackoilResidual& residual)
     {
         // the well flux equations are for each segment and each phase.
@@ -538,6 +537,7 @@ namespace Opm
         // 3. for the third term, it is the inflow through the perforations.
         // 4. for the last term, it is the outlet rates and also the segment rates,
         //    which are the primary variable.
+        const int np = numPhases();
         const int nseg_total = nseg_total_;
 
         ADB segqs = state.segqs;
@@ -580,7 +580,6 @@ namespace Opm
     addWellControlEq(const SolutionState& state,
                      const WellState& xw,
                      const Vector& aliveWells,
-                     const int np,
                      const std::vector<bool>& active,
                      LinearisedBlackoilResidual& residual)
     {
@@ -589,6 +588,7 @@ namespace Opm
         // And also, it work as the control equation when it is the segment
         if( wells().empty() ) return;
 
+        const int np = numPhases();
         const int nw = wells().size();
         const int nseg_total = nseg_total_;
 
@@ -753,7 +753,7 @@ namespace Opm
         std::string modestring[4] = { "BHP", "THP", "RESERVOIR_RATE", "SURFACE_RATE" };
         // Find, for each well, if any constraints are broken. If so,
         // switch control to first broken constraint.
-        const int np = wells()[0]->numberOfPhases();
+        const int np = numPhases();
         const int nw = wells().size();
         for (int w = 0; w < nw; ++w) {
             const WellControls* wc = wells()[w]->wellControls();