Merge pull request #5744 from akva2/add_ratio_calc

Add dedicated class for some ratio calculation in StandardWell

CMakeLists_files.cmake

@@ -186,6 +186,7 @@ list (APPEND MAIN_SOURCE_FILES
   opm/simulators/wells/ParallelWellInfo.cpp
   opm/simulators/wells/PerfData.cpp
   opm/simulators/wells/RateConverter.cpp
+  opm/simulators/wells/RatioCalculator.cpp
   opm/simulators/wells/SegmentState.cpp
   opm/simulators/wells/SingleWellState.cpp
   opm/simulators/wells/StandardWellAssemble.cpp
@@ -988,6 +989,7 @@ list (APPEND PUBLIC_HEADER_FILES
   opm/simulators/wells/PerfData.hpp
   opm/simulators/wells/PerforationData.hpp
   opm/simulators/wells/RateConverter.hpp
+  opm/simulators/wells/RatioCalculator.hpp
   opm/simulators/wells/RegionAttributeHelpers.hpp
   opm/simulators/wells/RegionAverageCalculator.hpp
   opm/simulators/wells/SingleWellState.hpp

opm/simulators/wells/RatioCalculator.cpp

@@ -0,0 +1,268 @@
+/*
+  Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
+  Copyright 2017 Statoil ASA.
+  Copyright 2016 - 2017 IRIS AS.
+
+  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/RatioCalculator.hpp>
+
+#include <opm/material/densead/Evaluation.hpp>
+#include <opm/material/densead/EvaluationFormat.hpp>
+
+#include <opm/simulators/utils/DeferredLogger.hpp>
+#include <opm/simulators/wells/PerforationData.hpp>
+
+#include <fmt/format.h>
+
+namespace {
+    template<class dValue, class Value>
+    auto dValueError(const dValue& d,
+                     const std::string& name,
+                     const std::string& methodName,
+                     const Value& Rs,
+                     const Value& Rv,
+                     const Value& pressure)
+    {
+        return fmt::format("Problematic d value {} obtained for well {}"
+                           " during {} calculations with rs {}"
+                           ", rv {} and pressure {}."
+                           " Continue as if no dissolution (rs = 0) and vaporization (rv = 0) "
+                           " for this connection.", d, name, methodName, Rs, Rv, pressure);
+    }
+}
+
+namespace Opm {
+
+template<class Value>
+RatioCalculator<Value>::
+RatioCalculator(unsigned gasCompIdx,
+                unsigned oilCompIdx,
+                unsigned waterCompIdx,
+                std::string_view name)
+    : gasComp_{gasCompIdx}
+    , oilComp_(oilCompIdx)
+    , waterComp_{waterCompIdx}
+    , name_(name)
+{
+}
+
+template<class Value>
+void
+RatioCalculator<Value>::
+disOilVapWatVolumeRatio(Value& volumeRatio,
+                        const Value& rvw,
+                        const Value& rsw,
+                        const Value& pressure,
+                        const std::vector<Value>& cmix_s,
+                        const std::vector<Value>& b_perfcells_dense,
+                        DeferredLogger& deferred_logger) const
+{
+    // Incorporate RSW/RVW factors if both water and gas active
+    const Value d = 1.0 - rvw * rsw;
+
+    if (d <= 0.0) {
+        deferred_logger.debug(dValueError(d, name_,
+                                          "disOilVapWatVolumeRatio",
+                                          rsw, rvw, pressure));
+    }
+    const Value tmp_wat = d > 0.0 ? (cmix_s[waterComp_] - rvw * cmix_s[gasComp_]) / d
+                                  : cmix_s[waterComp_];
+    volumeRatio += tmp_wat / b_perfcells_dense[waterComp_];
+
+    const Value tmp_gas =  d > 0.0 ? (cmix_s[gasComp_] - rsw * cmix_s[waterComp_]) / d
+                                   : cmix_s[gasComp_];
+    volumeRatio += tmp_gas / b_perfcells_dense[gasComp_];
+}
+
+template<class Value>
+void
+RatioCalculator<Value>::
+gasOilPerfRateInj(const std::vector<Value>& cq_s,
+                  PerforationRates<Scalar>& perf_rates,
+                  const Value& rv,
+                  const Value& rs,
+                  const Value& pressure,
+                  const Value& rvw,
+                  const bool waterActive,
+                  DeferredLogger& deferred_logger) const
+{
+    // TODO: the formulations here remain to be tested with cases with strong crossflow through production wells
+    // s means standard condition, r means reservoir condition
+    // q_os = q_or * b_o + rv * q_gr * b_g
+    // q_gs = q_gr * b_g + rs * q_or * b_o
+    // d = 1.0 - rs * rv
+    // q_or = 1 / (b_o * d) * (q_os - rv * q_gs)
+    // q_gr = 1 / (b_g * d) * (q_gs - rs * q_os)
+
+    const Scalar d = 1.0 - getValue(rv) * getValue(rs);
+
+    if (d <= 0.0) {
+        deferred_logger.debug(dValueError(d, name_,
+                                          "gasOilPerfRateInj",
+                                          rs, rv, pressure));
+    } else {
+        // vaporized oil into gas
+        // rv * q_gr * b_g = rv * (q_gs - rs * q_os) / d
+        perf_rates.vap_oil = getValue(rv) * (getValue(cq_s[gasComp_]) -
+                                             getValue(rs) * getValue(cq_s[oilComp_])) / d;
+        // dissolved of gas in oil
+        // rs * q_or * b_o = rs * (q_os - rv * q_gs) / d
+        perf_rates.dis_gas = getValue(rs) * (getValue(cq_s[oilComp_]) -
+                                             getValue(rv) * getValue(cq_s[gasComp_])) / d;
+    }
+
+    if (waterActive) {
+        // q_ws = q_wr * b_w + rvw * q_gr * b_g
+        // q_wr = 1 / b_w * (q_ws - rvw * q_gr * b_g) = 1 / b_w * (q_ws - rvw * 1 / d  (q_gs - rs * q_os))
+        // vaporized water in gas
+        // rvw * q_gr * b_g = q_ws -q_wr *b_w = rvw * (q_gs -rs *q_os) / d
+        perf_rates.vap_wat = getValue(rvw) * (getValue(cq_s[gasComp_]) -
+                                              getValue(rs) * getValue(cq_s[oilComp_])) / d;
+    }
+}
+
+template<class Value>
+void
+RatioCalculator<Value>::
+gasOilPerfRateProd(std::vector<Value>& cq_s,
+                   PerforationRates<Scalar>& perf_rates,
+                   const Value& rv,
+                   const Value& rs,
+                   const Value& rvw,
+                   const bool waterActive,
+                   const bool isProducer) const
+{
+    const Value cq_sOil = cq_s[oilComp_];
+    const Value cq_sGas = cq_s[gasComp_];
+    const Value dis_gas = rs * cq_sOil;
+    const Value vap_oil = rv * cq_sGas;
+
+    cq_s[gasComp_] += dis_gas;
+    cq_s[oilComp_] += vap_oil;
+
+    // recording the perforation solution gas rate and solution oil rates
+    if (isProducer) {
+        perf_rates.dis_gas = getValue(dis_gas);
+        perf_rates.vap_oil = getValue(vap_oil);
+    }
+
+    if (waterActive) {
+        const Value vap_wat = rvw * cq_sGas;
+        cq_s[waterComp_] += vap_wat;
+        if (isProducer) {
+            perf_rates.vap_wat = getValue(vap_wat);
+        }
+    }
+}
+
+template<class Value>
+void
+RatioCalculator<Value>::
+gasOilVolumeRatio(Value& volumeRatio,
+                  const Value& rv,
+                  const Value& rs,
+                  const Value& pressure,
+                  const std::vector<Value>& cmix_s,
+                  const std::vector<Value>& b_perfcells_dense,
+                  DeferredLogger& deferred_logger) const
+{
+    // Incorporate RS/RV factors if both oil and gas active
+    const Value d = 1.0 - rv * rs;
+
+    if (d <= 0.0) {
+        deferred_logger.debug(dValueError(d, name_,
+                                          "gasOilVolumeRatio",
+                                          rs, rv, pressure));
+    }
+    const Value tmp_oil = d > 0.0 ? (cmix_s[oilComp_] - rv * cmix_s[gasComp_]) / d
+                                  : cmix_s[oilComp_];
+    volumeRatio += tmp_oil / b_perfcells_dense[oilComp_];
+
+    const Value tmp_gas =  d > 0.0 ? (cmix_s[gasComp_] - rs * cmix_s[oilComp_]) / d
+                                   : cmix_s[gasComp_];
+    volumeRatio += tmp_gas / b_perfcells_dense[gasComp_];
+}
+
+template<class Value>
+void
+RatioCalculator<Value>::
+gasWaterPerfRateInj(const std::vector<Value>& cq_s,
+                     PerforationRates<Scalar>& perf_rates,
+                     const Value& rvw,
+                     const Value& rsw,
+                     const Value& pressure,
+                     DeferredLogger& deferred_logger) const
+{
+    const Scalar dw = 1.0 - getValue(rvw) * getValue(rsw);
+
+    if (dw <= 0.0) {
+        deferred_logger.debug(dValueError(dw, name_,
+                                          "gasWaterPerfRateInj",
+                                          rsw, rvw, pressure));
+    } else {
+        // vaporized water into gas
+        // rvw * q_gr * b_g = rvw * (q_gs - rsw * q_ws) / dw
+        perf_rates.vap_wat = getValue(rvw) * (getValue(cq_s[gasComp_]) -
+                                              getValue(rsw) * getValue(cq_s[waterComp_])) / dw;
+        // dissolved gas in water
+        // rsw * q_wr * b_w = rsw * (q_ws - rvw * q_gs) / dw
+        perf_rates.dis_gas_in_water = getValue(rsw) * (getValue(cq_s[waterComp_]) -
+                                                       getValue(rvw) * getValue(cq_s[gasComp_])) / dw;
+    }
+}
+
+template<class Value>
+void
+RatioCalculator<Value>::
+gasWaterPerfRateProd(std::vector<Value>& cq_s,
+                     PerforationRates<Scalar>& perf_rates,
+                     const Value& rvw,
+                     const Value& rsw,
+                     const bool isProducer) const
+{
+    const Value cq_sWat = cq_s[waterComp_];
+    const Value cq_sGas = cq_s[gasComp_];
+    const Value vap_wat = rvw * cq_sGas;
+    const Value dis_gas_wat = rsw * cq_sWat;
+    cq_s[waterComp_] += vap_wat;
+    cq_s[gasComp_]   += dis_gas_wat;
+    if (isProducer) {
+        perf_rates.vap_wat = getValue(vap_wat);
+        perf_rates.dis_gas_in_water = getValue(dis_gas_wat);
+    }
+}
+
+#define INSTANTIATE_TYPE(T)                                          \
+    template class RatioCalculator<T>;                               \
+    template class RatioCalculator<DenseAd::Evaluation<T, -1, 4u>>;  \
+    template class RatioCalculator<DenseAd::Evaluation<T, -1, 5u>>;  \
+    template class RatioCalculator<DenseAd::Evaluation<T, -1, 6u>>;  \
+    template class RatioCalculator<DenseAd::Evaluation<T, -1, 7u>>;  \
+    template class RatioCalculator<DenseAd::Evaluation<T, -1, 8u>>;  \
+    template class RatioCalculator<DenseAd::Evaluation<T, -1, 9u>>;  \
+    template class RatioCalculator<DenseAd::Evaluation<T, -1, 10u>>; \
+    template class RatioCalculator<DenseAd::Evaluation<T, -1, 11u>>;
+
+INSTANTIATE_TYPE(double)
+
+#if FLOW_INSTANTIATE_FLOAT
+INSTANTIATE_TYPE(float)
+#endif
+
+}

opm/simulators/wells/RatioCalculator.hpp

@@ -0,0 +1,102 @@
+/*
+  Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
+  Copyright 2017 Statoil ASA.
+  Copyright 2016 - 2017 IRIS AS.
+
+  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 RATIO_CALCULATOR_HPP
+#define RATIO_CALCULATOR_HPP
+
+#include <opm/material/densead/Math.hpp>
+
+#include <string>
+#include <string_view>
+#include <vector>
+
+namespace Opm {
+
+class DeferredLogger;
+template<class Scalar> struct PerforationRates;
+
+template<class Value>
+class RatioCalculator
+{
+public:
+    using Scalar = decltype(getValue(Value{}));
+
+    RatioCalculator(unsigned gasCompIdx,
+                    unsigned oilCompIdx,
+                    unsigned waterCompIdx,
+                    std::string_view name);
+
+    void disOilVapWatVolumeRatio(Value& volumeRatio,
+                                 const Value& rvw,
+                                 const Value& rsw,
+                                 const Value& pressure,
+                                 const std::vector<Value>& cmix_s,
+                                 const std::vector<Value>& b_perfcells_dense,
+                                 DeferredLogger& deferred_logger) const;
+
+    void gasOilPerfRateInj(const std::vector<Value>& cq_s,
+                           PerforationRates<Scalar>& perf_rates,
+                           const Value& rv,
+                           const Value& rs,
+                           const Value& pressure,
+                           const Value& rvw,
+                           const bool waterActive,
+                           DeferredLogger& deferred_logger) const;
+
+    void gasOilPerfRateProd(std::vector<Value>& cq_s,
+                            PerforationRates<Scalar>& perf_rates,
+                            const Value& rv,
+                            const Value& rs,
+                            const Value& rvw,
+                            const bool waterActive,
+                            const bool isProducer) const;
+
+    void gasOilVolumeRatio(Value& volumeRatio,
+                           const Value& rv,
+                           const Value& rs,
+                           const Value& pressure,
+                           const std::vector<Value>& cmix_s,
+                           const std::vector<Value>& b_perfcells_dense,
+                           DeferredLogger& deferred_logger) const;
+
+    void gasWaterPerfRateInj(const std::vector<Value>& cq_s,
+                             PerforationRates<Scalar>& perf_rates,
+                             const Value& rvw,
+                             const Value& rsw,
+                             const Value& pressure,
+                             DeferredLogger& deferred_logger) const;
+
+    void gasWaterPerfRateProd(std::vector<Value>& cq_s,
+                              PerforationRates<Scalar>& perf_rates,
+                              const Value& rvw,
+                              const Value& rsw,
+                              const bool isProducer) const;
+
+private:
+    unsigned gasComp_;
+    unsigned oilComp_;
+    unsigned waterComp_;
+    std::string name_;
+};
+
+} // namespace Opm
+
+#endif // RATIO_CALCULATOR_HPP

opm/simulators/wells/StandardWell.hpp

@@ -25,6 +25,7 @@
 
 #include <opm/simulators/timestepping/ConvergenceReport.hpp>
 #include <opm/simulators/wells/RateConverter.hpp>
+#include <opm/simulators/wells/RatioCalculator.hpp>
 #include <opm/simulators/wells/VFPInjProperties.hpp>
 #include <opm/simulators/wells/VFPProdProperties.hpp>
 #include <opm/simulators/wells/WellInterface.hpp>
@@ -466,54 +467,6 @@ namespace Opm
                                   const std::vector<EvalWell>& cq_s,
                                   const IntensiveQuantities& intQuants,
                                   DeferredLogger& deferred_logger) const;
-
-        template<class Value>
-        void gasOilPerfRateInj(const std::vector<Value>& cq_s,
-                               PerforationRates<Scalar>& perf_rates,
-                               const Value& rv,
-                               const Value& rs,
-                               const Value& pressure,
-                               const Value& rvw,
-                               DeferredLogger& deferred_logger) const;
-
-        template<class Value>
-        void gasOilPerfRateProd(std::vector<Value>& cq_s,
-                                PerforationRates<Scalar>& perf_rates,
-                                const Value& rv,
-                                const Value& rs,
-                                const Value& rvw) const;
-
-        template<class Value>
-        void gasWaterPerfRateProd(std::vector<Value>& cq_s,
-                                  PerforationRates<Scalar>& perf_rates,
-                                  const Value& rvw,
-                                  const Value& rsw) const;
-
-        template<class Value>
-        void gasWaterPerfRateInj(const std::vector<Value>& cq_s,
-                                 PerforationRates<Scalar>& perf_rates,
-                                 const Value& rvw,
-                                 const Value& rsw,
-                                 const Value& pressure,
-                                 DeferredLogger& deferred_logger) const;
-
-        template<class Value>
-        void disOilVapWatVolumeRatio(Value& volumeRatio,
-                                     const Value& rvw,
-                                     const Value& rsw,
-                                     const Value& pressure,
-                                     const std::vector<Value>& cmix_s,
-                                     const std::vector<Value>& b_perfcells_dense,
-                                     DeferredLogger& deferred_logger) const;
-
-        template<class Value>
-        void gasOilVolumeRatio(Value& volumeRatio,
-                               const Value& rv,
-                               const Value& rs,
-                               const Value& pressure,
-                               const std::vector<Value>& cmix_s,
-                               const std::vector<Value>& b_perfcells_dense,
-                               DeferredLogger& deferred_logger) const;
     };
 
 }

opm/simulators/wells/StandardWellEval.cpp

@@ -22,6 +22,7 @@
 #include <config.h>
 #include <opm/simulators/wells/StandardWellEval.hpp>
 
+#include <opm/material/densead/EvaluationFormat.hpp>
 #include <opm/material/fluidsystems/BlackOilFluidSystem.hpp>
 
 #include <opm/models/blackoil/blackoilindices.hh>
@@ -39,8 +40,9 @@
 #include <cmath>
 #include <cstddef>
 
-namespace Opm
+#include <fmt/format.h>
-{
+
+namespace Opm {
 
 template<class FluidSystem, class Indices>
 StandardWellEval<FluidSystem,Indices>::

opm/simulators/wells/StandardWellEval.hpp

@@ -19,7 +19,6 @@
   along with OPM.  If not, see <http://www.gnu.org/licenses/>.
 */
 
-
 #ifndef OPM_STANDARDWELL_EVAL_HEADER_INCLUDED
 #define OPM_STANDARDWELL_EVAL_HEADER_INCLUDED
 

opm/simulators/wells/StandardWell_impl.hpp

@@ -30,39 +30,17 @@
 
 #include <opm/input/eclipse/Units/Units.hpp>
 
-#include <opm/material/densead/EvaluationFormat.hpp>
-
 #include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
 #include <opm/simulators/wells/StandardWellAssemble.hpp>
 #include <opm/simulators/wells/VFPHelpers.hpp>
 #include <opm/simulators/wells/WellBhpThpCalculator.hpp>
 #include <opm/simulators/wells/WellConvergence.hpp>
 
-#include <fmt/format.h>
-
 #include <algorithm>
 #include <cstddef>
 #include <functional>
-#include <numeric>
 
-namespace {
+#include <fmt/format.h>
-
-template<class dValue, class Value>
-auto dValueError(const dValue& d,
-                 const std::string& name,
-                 const std::string& methodName,
-                 const Value& Rs,
-                 const Value& Rv,
-                 const Value& pressure)
-{
-    return fmt::format("Problematic d value {} obtained for well {}"
-                       " during {} calculations with rs {}"
-                       ", rv {} and pressure {}."
-                       " Continue as if no dissolution (rs = 0) and vaporization (rv = 0) "
-                       " for this connection.", d, name, methodName, Rs, Rv, pressure);
-}
-
-}
 
 namespace Opm
 {
@@ -246,6 +224,19 @@ namespace Opm
             drawdown += skin_pressure;
         }
 
+        RatioCalculator<Value> ratioCalc{
+            FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)
+                ? Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx)
+                : -1,
+            FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)
+                ? Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx)
+                : -1,
+            FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)
+                ? Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx)
+                : -1,
+            this->name()
+        };
+
         // producing perforations
         if (drawdown > 0)  {
             // Do nothing if crossflow is not allowed
@@ -259,10 +250,16 @@ namespace Opm
                 cq_s[componentIdx] = b_perfcells_dense[componentIdx] * cq_p;
             }
 
-            if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
+            if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) &&
-                gasOilPerfRateProd(cq_s, perf_rates, rv, rs, rvw);
+                FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx))
-            } else if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
+            {
-                gasWaterPerfRateProd(cq_s, perf_rates, rvw, rsw);
+                ratioCalc.gasOilPerfRateProd(cq_s, perf_rates, rv, rs, rvw,
+                                             FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx),
+                                             this->isProducer());
+            } else if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) &&
+                       FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx))
+            {
+                ratioCalc.gasWaterPerfRateProd(cq_s, perf_rates, rvw, rsw, this->isProducer());
             }
         } else {
             // Do nothing if crossflow is not allowed
@@ -280,7 +277,7 @@ namespace Opm
             Value volumeRatio = bhp * 0.0; // initialize it with the correct type
 
             if (FluidSystem::enableVaporizedWater() && FluidSystem::enableDissolvedGasInWater()) {
-                disOilVapWatVolumeRatio(volumeRatio, rvw, rsw, pressure,
+                ratioCalc.disOilVapWatVolumeRatio(volumeRatio, rvw, rsw, pressure,
                                                   cmix_s, b_perfcells_dense, deferred_logger);
                 // DISGASW only supported for gas-water CO2STORE/H2STORE case
                 // and the simulator will throw long before it reach to this point in the code
@@ -299,9 +296,12 @@ namespace Opm
                     volumeRatio += cmix_s[Indices::contiSolventEqIdx] / b_perfcells_dense[Indices::contiSolventEqIdx];
                 }
 
-                if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
+                if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) &&
-                    gasOilVolumeRatio(volumeRatio, rv, rs, pressure,
+                    FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx))
-                                      cmix_s, b_perfcells_dense, deferred_logger);
+                {
+                    ratioCalc.gasOilVolumeRatio(volumeRatio, rv, rs, pressure,
+                                                cmix_s, b_perfcells_dense,
+                                                deferred_logger);
                 } else {
                     if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
                         const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx);
@@ -323,13 +323,19 @@ namespace Opm
 
             // calculating the perforation solution gas rate and solution oil rates
             if (this->isProducer()) {
-                if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
+                if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) &&
-                    gasOilPerfRateInj(cq_s, perf_rates,
+                    FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx))
-                                      rv, rs, pressure, rvw, deferred_logger);
+                {
+                    ratioCalc.gasOilPerfRateInj(cq_s, perf_rates,
+                                                rv, rs, pressure, rvw,
+                                                FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx),
+                                                deferred_logger);
                 }
-                if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
+                if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) &&
+                    FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx))
+                {
                     //no oil
-                    gasWaterPerfRateInj(cq_s, perf_rates, rvw, rsw,
+                    ratioCalc.gasWaterPerfRateInj(cq_s, perf_rates, rvw, rsw,
                                                   pressure, deferred_logger);
                 }
             }
@@ -2639,207 +2645,4 @@ namespace Opm
 
         return result * this->well_efficiency_factor_;
     }
-
-
-    template <typename TypeTag>
-    template<class Value>
-    void
-    StandardWell<TypeTag>::
-    gasOilPerfRateInj(const std::vector<Value>& cq_s,
-                      PerforationRates<Scalar>& perf_rates,
-                      const Value& rv,
-                      const Value& rs,
-                      const Value& pressure,
-                      const Value& rvw,
-                      DeferredLogger& deferred_logger) const
-    {
-        const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx);
-        const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx);
-        // TODO: the formulations here remain to be tested with cases with strong crossflow through production wells
-        // s means standard condition, r means reservoir condition
-        // q_os = q_or * b_o + rv * q_gr * b_g
-        // q_gs = q_gr * b_g + rs * q_or * b_o
-        // d = 1.0 - rs * rv
-        // q_or = 1 / (b_o * d) * (q_os - rv * q_gs)
-        // q_gr = 1 / (b_g * d) * (q_gs - rs * q_os)
-
-        const Scalar d = 1.0 - getValue(rv) * getValue(rs);
-
-        if (d <= 0.0) {
-            deferred_logger.debug(dValueError(d, this->name(),
-                                              "gasOilPerfRateInj",
-                                              rs, rv, pressure));
-        } else {
-            // vaporized oil into gas
-            // rv * q_gr * b_g = rv * (q_gs - rs * q_os) / d
-            perf_rates.vap_oil = getValue(rv) * (getValue(cq_s[gasCompIdx]) - getValue(rs) * getValue(cq_s[oilCompIdx])) / d;
-            // dissolved of gas in oil
-            // rs * q_or * b_o = rs * (q_os - rv * q_gs) / d
-            perf_rates.dis_gas = getValue(rs) * (getValue(cq_s[oilCompIdx]) - getValue(rv) * getValue(cq_s[gasCompIdx])) / d;
-        }
-
-        if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
-            // q_ws = q_wr * b_w + rvw * q_gr * b_g
-            // q_wr = 1 / b_w * (q_ws - rvw * q_gr * b_g) = 1 / b_w * (q_ws - rvw * 1 / d  (q_gs - rs * q_os))
-            // vaporized water in gas
-            // rvw * q_gr * b_g = q_ws -q_wr *b_w = rvw * (q_gs -rs *q_os) / d
-            perf_rates.vap_wat = getValue(rvw) * (getValue(cq_s[gasCompIdx]) - getValue(rs) * getValue(cq_s[oilCompIdx])) / d;
-        }
-    }
-
-
-
-    template <typename TypeTag>
-    template<class Value>
-    void
-    StandardWell<TypeTag>::
-    gasOilPerfRateProd(std::vector<Value>& cq_s,
-                       PerforationRates<Scalar>& perf_rates,
-                       const Value& rv,
-                       const Value& rs,
-                       const Value& rvw) const
-    {
-        const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx);
-        const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx);
-        const Value cq_sOil = cq_s[oilCompIdx];
-        const Value cq_sGas = cq_s[gasCompIdx];
-        const Value dis_gas = rs * cq_sOil;
-        const Value vap_oil = rv * cq_sGas;
-
-        cq_s[gasCompIdx] += dis_gas;
-        cq_s[oilCompIdx] += vap_oil;
-
-        // recording the perforation solution gas rate and solution oil rates
-        if (this->isProducer()) {
-            perf_rates.dis_gas = getValue(dis_gas);
-            perf_rates.vap_oil = getValue(vap_oil);
-        }
-
-        if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
-            const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx);
-            const Value vap_wat = rvw * cq_sGas;
-            cq_s[waterCompIdx] += vap_wat;
-            if (this->isProducer())
-                perf_rates.vap_wat = getValue(vap_wat);
-        }
-    }
-
-
-    template <typename TypeTag>
-    template<class Value>
-    void
-    StandardWell<TypeTag>::
-    gasWaterPerfRateProd(std::vector<Value>& cq_s,
-                         PerforationRates<Scalar>& perf_rates,
-                         const Value& rvw,
-                         const Value& rsw) const
-    {
-        const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx);
-        const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx);
-        const Value cq_sWat = cq_s[waterCompIdx];
-        const Value cq_sGas = cq_s[gasCompIdx];
-        const Value vap_wat = rvw * cq_sGas;
-        const Value dis_gas_wat = rsw * cq_sWat;
-        cq_s[waterCompIdx] += vap_wat;
-        cq_s[gasCompIdx]   += dis_gas_wat;
-        if (this->isProducer()) {
-            perf_rates.vap_wat = getValue(vap_wat);
-            perf_rates.dis_gas_in_water = getValue(dis_gas_wat);
-        }
-    }
-
-
-    template <typename TypeTag>
-    template<class Value>
-    void
-    StandardWell<TypeTag>::
-    gasWaterPerfRateInj(const std::vector<Value>& cq_s,
-                        PerforationRates<Scalar>& perf_rates,
-                        const Value& rvw,
-                        const Value& rsw,
-                        const Value& pressure,
-                        DeferredLogger& deferred_logger) const
-
-    {
-        const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx);
-        const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx);
-
-        const Scalar dw = 1.0 - getValue(rvw) * getValue(rsw);
-
-        if (dw <= 0.0) {
-            deferred_logger.debug(dValueError(dw, this->name(),
-                                              "gasWaterPerfRateInj",
-                                              rsw, rvw, pressure));
-        } else {
-            // vaporized water into gas
-            // rvw * q_gr * b_g = rvw * (q_gs - rsw * q_ws) / dw
-            perf_rates.vap_wat = getValue(rvw) * (getValue(cq_s[gasCompIdx]) - getValue(rsw) * getValue(cq_s[waterCompIdx])) / dw;
-            // dissolved gas in water
-            // rsw * q_wr * b_w = rsw * (q_ws - rvw * q_gs) / dw
-            perf_rates.dis_gas_in_water = getValue(rsw) * (getValue(cq_s[waterCompIdx]) - getValue(rvw) * getValue(cq_s[gasCompIdx])) / dw;
-        }
-    }
-
-
-    template <typename TypeTag>
-    template<class Value>
-    void
-    StandardWell<TypeTag>::
-    disOilVapWatVolumeRatio(Value& volumeRatio,
-                            const Value& rvw,
-                            const Value& rsw,
-                            const Value& pressure,
-                            const std::vector<Value>& cmix_s,
-                            const std::vector<Value>& b_perfcells_dense,
-                            DeferredLogger& deferred_logger) const
-    {
-        const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx);
-        const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx);
-        // Incorporate RSW/RVW factors if both water and gas active
-        const Value d = 1.0 - rvw * rsw;
-
-        if (d <= 0.0) {
-            deferred_logger.debug(dValueError(d, this->name(),
-                                              "disOilVapWatVolumeRatio",
-                                              rsw, rvw, pressure));
-        }
-        const Value tmp_wat = d > 0.0 ? (cmix_s[waterCompIdx] - rvw * cmix_s[gasCompIdx]) / d
-                                      : cmix_s[waterCompIdx];
-        volumeRatio += tmp_wat / b_perfcells_dense[waterCompIdx];
-
-        const Value tmp_gas =  d > 0.0 ? (cmix_s[gasCompIdx] - rsw * cmix_s[waterCompIdx]) / d
-                                       : cmix_s[gasCompIdx];
-        volumeRatio += tmp_gas / b_perfcells_dense[gasCompIdx];
-    }
-
-
-    template <typename TypeTag>
-    template<class Value>
-    void
-    StandardWell<TypeTag>::
-    gasOilVolumeRatio(Value& volumeRatio,
-                      const Value& rv,
-                      const Value& rs,
-                      const Value& pressure,
-                      const std::vector<Value>& cmix_s,
-                      const std::vector<Value>& b_perfcells_dense,
-                      DeferredLogger& deferred_logger) const
-    {
-        const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx);
-        const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx);
-        // Incorporate RS/RV factors if both oil and gas active
-        const Value d = 1.0 - rv * rs;
-
-        if (d <= 0.0) {
-            deferred_logger.debug(dValueError(d, this->name(),
-                                              "gasOilVolumeRatio",
-                                              rs, rv, pressure));
-        }
-        const Value tmp_oil = d > 0.0? (cmix_s[oilCompIdx] - rv * cmix_s[gasCompIdx]) / d : cmix_s[oilCompIdx];
-        volumeRatio += tmp_oil / b_perfcells_dense[oilCompIdx];
-
-        const Value tmp_gas =  d > 0.0? (cmix_s[gasCompIdx] - rs * cmix_s[oilCompIdx]) / d : cmix_s[gasCompIdx];
-        volumeRatio += tmp_gas / b_perfcells_dense[gasCompIdx];
-    }
-
 } // namespace Opm