add WellInterfaceIndices

while not a lot of class in this layer, having it helps
in downstream well code.

CMakeLists_files.cmake

@@ -68,6 +68,7 @@ list (APPEND MAIN_SOURCE_FILES
   opm/simulators/wells/WellGroupHelpers.cpp
   opm/simulators/wells/WellInterfaceFluidSystem.cpp
   opm/simulators/wells/WellInterfaceGeneric.cpp
+  opm/simulators/wells/WellInterfaceIndices.cpp
   opm/simulators/wells/WellProdIndexCalculator.cpp
   opm/simulators/wells/WellState.cpp
   opm/simulators/wells/WGState.cpp

opm/simulators/wells/WellInterface.hpp

@@ -59,7 +59,7 @@ namespace Opm {
 #include <opm/material/densead/Math.hpp>
 #include <opm/material/densead/Evaluation.hpp>
 
-#include <opm/simulators/wells/WellInterfaceFluidSystem.hpp>
+#include <opm/simulators/wells/WellInterfaceIndices.hpp>
 
 #include <array>
 #include <cassert>
@@ -71,7 +71,9 @@ namespace Opm
 {
 
 template<typename TypeTag>
-class WellInterface : public WellInterfaceFluidSystem<GetPropType<TypeTag, Properties::FluidSystem>>
+class WellInterface : public WellInterfaceIndices<GetPropType<TypeTag, Properties::FluidSystem>,
+                                                  GetPropType<TypeTag, Properties::Indices>,
+                                                  GetPropType<TypeTag, Properties::Scalar>>
 {
 public:
 
@@ -292,10 +294,6 @@ protected:
 
     bool changed_to_stopped_this_step_ = false;
 
-    int flowPhaseToEbosCompIdx( const int phaseIdx ) const;
-
-    int ebosCompIdxToFlowCompIdx( const unsigned compIdx ) const;
-
     double wpolymer() const;
 
     double wfoam() const;
@@ -310,8 +308,6 @@ protected:
     // Component fractions for each phase for the well
     const std::vector<double>& compFrac() const;
 
-    double scalingFactor(const int comp_idx) const;
-
     std::vector<double> initialWellRateFractions(const Simulator& ebosSimulator, const WellState& well_state) const;
 
     // check whether the well is operable under BHP limit with current reservoir condition

opm/simulators/wells/WellInterfaceFluidSystem.hpp

@@ -52,14 +52,14 @@ public:
 
     int flowPhaseToEbosPhaseIdx(const int phaseIdx) const;
 
-protected:
-    using RateConverterType = RateConverter::
-    SurfaceToReservoirVoidage<FluidSystem, std::vector<int>>;
-
     static constexpr int Water = BlackoilPhases::Aqua;
     static constexpr int Oil = BlackoilPhases::Liquid;
     static constexpr int Gas = BlackoilPhases::Vapour;
 
+protected:
+    using RateConverterType = RateConverter::
+    SurfaceToReservoirVoidage<FluidSystem, std::vector<int>>;
+
     // to indicate a invalid completion
     static constexpr int INVALIDCOMPLETION = std::numeric_limits<int>::max();
 

opm/simulators/wells/WellInterfaceIndices.cpp

@@ -0,0 +1,145 @@
+/*
+  Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
+  Copyright 2017 Statoil ASA.
+  Copyright 2018 IRIS
+
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 3 of the License, or
+  (at your option) any later version.
+
+  OPM is distributed in the hope that it will be useful,
+  but WITHOUT ANY WARRANTY; without even the implied warranty of
+  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+  GNU General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+*/
+
+#include <config.h>
+#include <opm/simulators/wells/WellInterfaceIndices.hpp>
+
+#include <ebos/eclalternativeblackoilindices.hh>
+
+#include <opm/material/fluidsystems/BlackOilFluidSystem.hpp>
+#include <opm/models/blackoil/blackoilindices.hh>
+#include <opm/models/blackoil/blackoilonephaseindices.hh>
+#include <opm/models/blackoil/blackoiltwophaseindices.hh>
+
+#include <cassert>
+
+namespace Opm
+{
+
+template<class FluidSystem, class Indices, class Scalar>
+WellInterfaceIndices<FluidSystem,Indices,Scalar>::
+WellInterfaceIndices(const Well& well,
+                     const ParallelWellInfo& parallel_well_info,
+                     const int time_step,
+                     const typename WellInterfaceFluidSystem<FluidSystem>::RateConverterType& rate_converter,
+                     const int pvtRegionIdx,
+                     const int num_components,
+                     const int num_phases,
+                     const int index_of_well,
+                     const int first_perf_index,
+                     const std::vector<PerforationData>& perf_data)
+    : WellInterfaceFluidSystem<FluidSystem>(well,
+                                            parallel_well_info,
+                                            time_step,
+                                            rate_converter,
+                                            pvtRegionIdx,
+                                            num_components,
+                                            num_phases,
+                                            index_of_well,
+                                            first_perf_index,
+                                            perf_data)
+{
+}
+
+template<class FluidSystem, class Indices, class Scalar>
+int
+WellInterfaceIndices<FluidSystem,Indices,Scalar>::
+flowPhaseToEbosCompIdx(const int phaseIdx) const
+{
+    const auto& pu = this->phaseUsage();
+    if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) && pu.phase_pos[Water] == phaseIdx)
+        return Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx);
+    if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && pu.phase_pos[Oil] == phaseIdx)
+        return Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx);
+    if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) && pu.phase_pos[Gas] == phaseIdx)
+        return Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx);
+
+    // for other phases return the index
+    return phaseIdx;
+}
+
+template<class FluidSystem, class Indices, class Scalar>
+int
+WellInterfaceIndices<FluidSystem,Indices,Scalar>::
+ebosCompIdxToFlowCompIdx(const unsigned compIdx) const
+{
+    const auto& pu = this->phaseUsage();
+    if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) && Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx) == compIdx)
+        return pu.phase_pos[Water];
+    if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx) == compIdx)
+        return pu.phase_pos[Oil];
+    if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) && Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx) == compIdx)
+        return pu.phase_pos[Gas];
+
+    // for other phases return the index
+    return compIdx;
+}
+
+template<class FluidSystem, class Indices, class Scalar>
+double
+WellInterfaceIndices<FluidSystem,Indices,Scalar>::
+scalingFactor(const int phaseIdx) const
+{
+    const auto& pu = this->phaseUsage();
+    if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) && pu.phase_pos[Water] == phaseIdx)
+        return 1.0;
+    if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && pu.phase_pos[Oil] == phaseIdx)
+        return 1.0;
+    if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) && pu.phase_pos[Gas] == phaseIdx)
+        return 0.01;
+    if (Indices::enableSolvent && phaseIdx == Indices::contiSolventEqIdx )
+        return 0.01;
+
+    // we should not come this far
+    assert(false);
+    return 1.0;
+}
+
+#define INSTANCE(A, ...) \
+template class WellInterfaceIndices<BlackOilFluidSystem<double, A>, \
+                                    __VA_ARGS__, \
+                                    double>;
+
+// One phase
+INSTANCE(BlackOilDefaultIndexTraits,BlackOilOnePhaseIndices<0u,0u,0u,0u,false,false,0u,1u>)
+INSTANCE(BlackOilDefaultIndexTraits,BlackOilOnePhaseIndices<0u,0u,0u,1u,false,false,0u,1u>)
+
+// Two phase
+INSTANCE(BlackOilDefaultIndexTraits,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,0u>)
+INSTANCE(BlackOilDefaultIndexTraits,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,1u>)
+INSTANCE(BlackOilDefaultIndexTraits,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,2u>)
+INSTANCE(BlackOilDefaultIndexTraits,BlackOilTwoPhaseIndices<0u,0u,1u,0u,false,false,0u,2u>)
+INSTANCE(BlackOilDefaultIndexTraits,BlackOilTwoPhaseIndices<0u,0u,2u,0u,false,false,0u,2u>)
+INSTANCE(BlackOilDefaultIndexTraits,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,true,0u,2u>)
+
+// Blackoil
+INSTANCE(BlackOilDefaultIndexTraits,BlackOilIndices<0u,0u,0u,0u,false,false,0u>)
+INSTANCE(BlackOilDefaultIndexTraits,BlackOilIndices<0u,0u,0u,0u,true,false,0u>)
+INSTANCE(BlackOilDefaultIndexTraits,BlackOilIndices<0u,0u,0u,0u,false,true,0u>)
+INSTANCE(BlackOilDefaultIndexTraits,BlackOilIndices<1u,0u,0u,0u,false,false,0u>)
+INSTANCE(BlackOilDefaultIndexTraits,BlackOilIndices<0u,1u,0u,0u,false,false,0u>)
+INSTANCE(BlackOilDefaultIndexTraits,BlackOilIndices<0u,0u,1u,0u,false,false,0u>)
+INSTANCE(BlackOilDefaultIndexTraits,BlackOilIndices<0u,0u,0u,1u,false,false,0u>)
+
+// Alternative indices
+INSTANCE(EclAlternativeBlackOilIndexTraits,BlackOilIndices<0u,0u,0u,0u,false,false,0u>)
+
+} // namespace Opm

opm/simulators/wells/WellInterfaceIndices.hpp

@@ -0,0 +1,61 @@
+/*
+  Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
+  Copyright 2017 Statoil ASA.
+  Copyright 2017 IRIS
+  Copyright 2019 Norce
+
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 3 of the License, or
+  (at your option) any later version.
+
+  OPM is distributed in the hope that it will be useful,
+  but WITHOUT ANY WARRANTY; without even the implied warranty of
+  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+  GNU General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+*/
+
+#ifndef OPM_WELLINTERFACE_INDICES_HEADER_INCLUDED
+#define OPM_WELLINTERFACE_INDICES_HEADER_INCLUDED
+
+#include <opm/simulators/wells/WellInterfaceFluidSystem.hpp>
+
+namespace Opm
+{
+
+template<class FluidSystem, class Indices, class Scalar>
+class WellInterfaceIndices : public WellInterfaceFluidSystem<FluidSystem> {
+public:
+    using WellInterfaceFluidSystem<FluidSystem>::Gas;
+    using WellInterfaceFluidSystem<FluidSystem>::Oil;
+    using WellInterfaceFluidSystem<FluidSystem>::Water;
+
+protected:
+    WellInterfaceIndices(const Well& well,
+                         const ParallelWellInfo& parallel_well_info,
+                         const int time_step,
+                         const typename WellInterfaceFluidSystem<FluidSystem>::RateConverterType& rate_converter,
+                         const int pvtRegionIdx,
+                         const int num_components,
+                         const int num_phases,
+                         const int index_of_well,
+                         const int first_perf_index,
+                         const std::vector<PerforationData>& perf_data);
+
+    int flowPhaseToEbosCompIdx( const int phaseIdx ) const;
+    int ebosCompIdxToFlowCompIdx( const unsigned compIdx ) const;
+
+    double scalingFactor(const int phaseIdx) const;
+
+    Scalar volumetricSurfaceRateForConnection(int cellIdx,
+                                              int phaseIdx) const;
+};
+
+}
+
+#endif // OPM_WELLINTERFACE_INDICES_HEADER_INCLUDED

opm/simulators/wells/WellInterface_impl.hpp

@@ -41,9 +41,16 @@ namespace Opm
                   const int index_of_well,
                   const int first_perf_index,
                   const std::vector<PerforationData>& perf_data)
-      : WellInterfaceFluidSystem<FluidSystem>(well, pw_info, time_step, rate_converter,
-                                              pvtRegionIdx, num_components, num_phases,
-                                              index_of_well, first_perf_index, perf_data)
+      : WellInterfaceIndices<FluidSystem,Indices,Scalar>(well,
+                                                         pw_info,
+                                                         time_step,
+                                                         rate_converter,
+                                                         pvtRegionIdx,
+                                                         num_components,
+                                                         num_phases,
+                                                         index_of_well,
+                                                         first_perf_index,
+                                                         perf_data)
       , param_(param)
     {
         connectionRates_.resize(this->number_of_perforations_);
@@ -74,41 +81,6 @@ namespace Opm
     }
 
 
-    template<typename TypeTag>
-    int
-    WellInterface<TypeTag>::
-    flowPhaseToEbosCompIdx( const int phaseIdx ) const
-    {
-        const auto& pu = this->phaseUsage();
-        if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) && pu.phase_pos[Water] == phaseIdx)
-            return Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx);
-        if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && pu.phase_pos[Oil] == phaseIdx)
-            return Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx);
-        if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) && pu.phase_pos[Gas] == phaseIdx)
-            return Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx);
-
-        // for other phases return the index
-        return phaseIdx;
-    }
-
-    template<typename TypeTag>
-    int
-    WellInterface<TypeTag>::
-    ebosCompIdxToFlowCompIdx( const unsigned compIdx ) const
-    {
-        const auto& pu = this->phaseUsage();
-        if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) && Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx) == compIdx)
-            return pu.phase_pos[Water];
-        if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx) == compIdx)
-            return pu.phase_pos[Oil];
-        if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) && Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx) == compIdx)
-            return pu.phase_pos[Gas];
-
-        // for other phases return the index
-        return compIdx;
-    }
-
-
 
 
     template<typename TypeTag>
@@ -368,28 +340,6 @@ namespace Opm
 
 
 
-    template<typename TypeTag>
-    double
-    WellInterface<TypeTag>::scalingFactor(const int phaseIdx) const
-    {
-        const auto& pu = this->phaseUsage();
-        if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) && pu.phase_pos[Water] == phaseIdx)
-            return 1.0;
-        if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && pu.phase_pos[Oil] == phaseIdx)
-            return 1.0;
-        if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) && pu.phase_pos[Gas] == phaseIdx)
-            return 0.01;
-        if (has_solvent && phaseIdx == contiSolventEqIdx )
-            return 0.01;
-
-        // we should not come this far
-        assert(false);
-        return 1.0;
-    }
-
-
-
-
     template<typename TypeTag>
     bool
     WellInterface<TypeTag>::
@@ -1448,10 +1398,10 @@ namespace Opm
 
         // Set the currently-zero phase flows to be nonzero in proportion to well_q_s.
         const double initial_nonzero_rate = well_state.wellRates(this->index_of_well_)[nonzero_rate_index];
-        const int comp_idx_nz = flowPhaseToEbosCompIdx(nonzero_rate_index);
+        const int comp_idx_nz = this->flowPhaseToEbosCompIdx(nonzero_rate_index);
         for (int p = 0; p < this->number_of_phases_; ++p) {
             if (p != nonzero_rate_index) {
-                const int comp_idx = flowPhaseToEbosCompIdx(p);
+                const int comp_idx = this->flowPhaseToEbosCompIdx(p);
                 double& rate = well_state.wellRates(this->index_of_well_)[p];
                 rate = (initial_nonzero_rate/well_q_s[comp_idx_nz]) * (well_q_s[comp_idx]);
             }