diff --git a/opm/models/blackoil/blackoilnewtonmethod.hh b/opm/models/blackoil/blackoilnewtonmethod.hh
index 583d1c7eb..05eed406a 100644
--- a/opm/models/blackoil/blackoilnewtonmethod.hh
+++ b/opm/models/blackoil/blackoilnewtonmethod.hh
@@ -55,6 +55,10 @@ template<class TypeTag, class MyTypeTag>
 struct TemperatureMax { using type = UndefinedProperty; };
 template<class TypeTag, class MyTypeTag>
 struct TemperatureMin { using type = UndefinedProperty; };
+template<class TypeTag, class MyTypeTag>
+struct MaximumWaterSaturation { using type = UndefinedProperty; };
+template<class TypeTag, class MyTypeTag>
+struct WaterOnlyThreshold { using type = UndefinedProperty; };
 template<class TypeTag>
 struct DpMaxRel<TypeTag, TTag::NewtonMethod>
 {
@@ -93,6 +97,18 @@ struct TemperatureMin<TypeTag, TTag::NewtonMethod>
     using type = GetPropType<TypeTag, Scalar>;
     static constexpr type value = 0.0; //Kelvin
 };
+template<class TypeTag>
+struct MaximumWaterSaturation<TypeTag, TTag::NewtonMethod>
+{
+    using type = GetPropType<TypeTag, Scalar>;
+    static constexpr type value = 1.0;
+};
+template<class TypeTag>
+struct WaterOnlyThreshold<TypeTag, TTag::NewtonMethod>
+{
+    using type = GetPropType<TypeTag, Scalar>;
+    static constexpr type value = 1.0;
+};
 } // namespace Opm::Properties
 
 namespace Opm {
@@ -130,7 +146,8 @@ public:
         maxTempChange_ = EWOMS_GET_PARAM(TypeTag, Scalar, MaxTemperatureChange);
         tempMax_ = EWOMS_GET_PARAM(TypeTag, Scalar, TemperatureMax);
         tempMin_ = EWOMS_GET_PARAM(TypeTag, Scalar, TemperatureMin);
-
+        waterSaturationMax_ = EWOMS_GET_PARAM(TypeTag, Scalar, MaximumWaterSaturation);
+        waterOnlyThreshold_ = EWOMS_GET_PARAM(TypeTag, Scalar, WaterOnlyThreshold);
     }
 
     /*!
@@ -159,6 +176,8 @@ public:
         EWOMS_REGISTER_PARAM(TypeTag, Scalar, MaxTemperatureChange, "Maximum absolute change of temperature in a single iteration");
         EWOMS_REGISTER_PARAM(TypeTag, Scalar, TemperatureMax, "Maximum absolute temperature");
         EWOMS_REGISTER_PARAM(TypeTag, Scalar, TemperatureMin, "Minimum absolute temperature");
+        EWOMS_REGISTER_PARAM(TypeTag, Scalar, MaximumWaterSaturation, "Maximum water saturation");
+        EWOMS_REGISTER_PARAM(TypeTag, Scalar, WaterOnlyThreshold, "Cells with water saturation above or equal is considered one-phase water only");
     }
 
     /*!
@@ -433,9 +452,9 @@ protected:
         // use a threshold value after a switch to make it harder to switch back
         // immediately.
         if (wasSwitched_[globalDofIdx])
-            wasSwitched_[globalDofIdx] = nextValue.adaptPrimaryVariables(this->problem(), globalDofIdx, priVarOscilationThreshold_);
+            wasSwitched_[globalDofIdx] = nextValue.adaptPrimaryVariables(this->problem(), globalDofIdx, waterSaturationMax_, waterOnlyThreshold_, priVarOscilationThreshold_);
         else
-            wasSwitched_[globalDofIdx] = nextValue.adaptPrimaryVariables(this->problem(), globalDofIdx);
+            wasSwitched_[globalDofIdx] = nextValue.adaptPrimaryVariables(this->problem(), globalDofIdx, waterSaturationMax_, waterOnlyThreshold_);
 
         if (wasSwitched_[globalDofIdx])
             ++ numPriVarsSwitched_;
@@ -450,6 +469,9 @@ private:
     int numPriVarsSwitched_;
 
     Scalar priVarOscilationThreshold_;
+    Scalar waterSaturationMax_;
+    Scalar waterOnlyThreshold_;
+
     Scalar dpMaxRel_;
     Scalar dsMax_;
     bool projectSaturations_;
diff --git a/opm/models/blackoil/blackoilprimaryvariables.hh b/opm/models/blackoil/blackoilprimaryvariables.hh
index 88aba097c..b7925e2fb 100644
--- a/opm/models/blackoil/blackoilprimaryvariables.hh
+++ b/opm/models/blackoil/blackoilprimaryvariables.hh
@@ -493,10 +493,8 @@ public:
      *
      * \return true Iff the interpretation of one of the switching variables was changed
      */
-    bool adaptPrimaryVariables(const Problem& problem, unsigned globalDofIdx, Scalar eps = 0.0)
+    bool adaptPrimaryVariables(const Problem& problem, unsigned globalDofIdx, Scalar swMaximum, Scalar thresholdWaterFilledCell, Scalar eps = 0.0)
     {
-        static const Scalar thresholdWaterFilledCell = 1.0 - eps;
-
         // this function accesses quite a few black-oil specific low-level functions
         // directly for better performance (instead of going the canonical way through
         // the IntensiveQuantities). The reason is that most intensive quantities are not
@@ -546,16 +544,17 @@ public:
         // keep track if any meaning has changed
         bool changed = false;
 
-        // special case cells with almost only water
-        // use both saturations (if the phase is enabled)
-        // if dissolved gas in water is enabled we shouldn't enter
+        // Special case for cells with almost only water
+        // for these cells both saturations (if the phase is enabled) is used
+        // to avoid singular systems.
+        // If dissolved gas in water is enabled we shouldn't enter
         // here but instead switch to Rsw as primary variable
         // as sw >= 1.0 -> gas <= 0 (i.e. gas phase disappears)
         if (sw >= thresholdWaterFilledCell && !FluidSystem::enableDissolvedGasInWater()) {
 
-            // make sure water saturations does not exceed 1.0
+            // make sure water saturations does not exceed sw_maximum. Default to 1.0
             if constexpr (waterEnabled) {
-                (*this)[Indices::waterSwitchIdx] = 1.0;
+                (*this)[Indices::waterSwitchIdx] = std::min(swMaximum, sw);
                 assert(primaryVarsMeaningWater() == WaterMeaning::Sw);
             }
             // the hydrocarbon gas saturation is set to 0.0