diff --git a/opm/core/fluid/IncompPropertiesBasic.cpp b/opm/core/fluid/IncompPropertiesBasic.cpp
index 004aed94..1b78481a 100644
--- a/opm/core/fluid/IncompPropertiesBasic.cpp
+++ b/opm/core/fluid/IncompPropertiesBasic.cpp
@@ -136,6 +136,20 @@ namespace Opm
     }
 
 
+    /// Obtain the range of allowable saturation values.
+    /// In cell cells[i], saturation of phase p is allowed to be
+    /// in the interval [smin[i*P + p], smax[i*P + p]].
+    /// \param[in]  n      Number of data points.
+    /// \param[in]  cells  Array of n cell indices.
+    /// \param[out] smin   Array of nP minimum s values, array must be valid before calling.
+    /// \param[out] smax   Array of nP maximum s values, array must be valid before calling.
+    void IncompPropertiesBasic::satRange(const int n,
+					 const int* /*cells*/,
+					 double* smin,
+					 double* smax) const
+    {
+	satprops_.satRange(n, smin, smax);
+    }
 
 } // namespace Opm
 
diff --git a/opm/core/fluid/IncompPropertiesBasic.hpp b/opm/core/fluid/IncompPropertiesBasic.hpp
index db8ab535..e52044af 100644
--- a/opm/core/fluid/IncompPropertiesBasic.hpp
+++ b/opm/core/fluid/IncompPropertiesBasic.hpp
@@ -114,6 +114,19 @@ namespace Opm
                               const int* cells,
                               double* pc,
                               double* dpcds) const;
+
+
+	/// Obtain the range of allowable saturation values.
+	/// In cell cells[i], saturation of phase p is allowed to be
+	/// in the interval [smin[i*P + p], smax[i*P + p]].
+        /// \param[in]  n      Number of data points.
+        /// \param[in]  cells  Array of n cell indices.
+        /// \param[out] smin   Array of nP minimum s values, array must be valid before calling.
+        /// \param[out] smax   Array of nP maximum s values, array must be valid before calling.
+        virtual void satRange(const int n,
+                              const int* cells,
+                              double* smin,
+                              double* smax) const;
     private:
         RockBasic rock_;
 	PvtPropertiesBasic pvt_;
diff --git a/opm/core/fluid/IncompPropertiesFromDeck.cpp b/opm/core/fluid/IncompPropertiesFromDeck.cpp
index 4f800568..519615ea 100644
--- a/opm/core/fluid/IncompPropertiesFromDeck.cpp
+++ b/opm/core/fluid/IncompPropertiesFromDeck.cpp
@@ -128,6 +128,20 @@ namespace Opm
     }
 
 
+    /// Obtain the range of allowable saturation values.
+    /// In cell cells[i], saturation of phase p is allowed to be
+    /// in the interval [smin[i*P + p], smax[i*P + p]].
+    /// \param[in]  n      Number of data points.
+    /// \param[in]  cells  Array of n cell indices.
+    /// \param[out] smin   Array of nP minimum s values, array must be valid before calling.
+    /// \param[out] smax   Array of nP maximum s values, array must be valid before calling.
+    void IncompPropertiesFromDeck::satRange(const int n,
+					    const int* /*cells*/,
+					    double* smin,
+					    double* smax) const
+    {
+	satprops_.satRange(n, smin, smax);
+    }
 
 } // namespace Opm
 
diff --git a/opm/core/fluid/IncompPropertiesFromDeck.hpp b/opm/core/fluid/IncompPropertiesFromDeck.hpp
index de3ba043..1472f043 100644
--- a/opm/core/fluid/IncompPropertiesFromDeck.hpp
+++ b/opm/core/fluid/IncompPropertiesFromDeck.hpp
@@ -111,6 +111,19 @@ namespace Opm
                               const int* cells,
                               double* pc,
                               double* dpcds) const;
+
+
+	/// Obtain the range of allowable saturation values.
+	/// In cell cells[i], saturation of phase p is allowed to be
+	/// in the interval [smin[i*P + p], smax[i*P + p]].
+        /// \param[in]  n      Number of data points.
+        /// \param[in]  cells  Array of n cell indices.
+        /// \param[out] smin   Array of nP minimum s values, array must be valid before calling.
+        /// \param[out] smax   Array of nP maximum s values, array must be valid before calling.
+        virtual void satRange(const int n,
+                              const int* cells,
+                              double* smin,
+                              double* smax) const;
     private:
         RockFromDeck rock_;
 	PvtPropertiesIncompFromDeck pvt_;
diff --git a/opm/core/fluid/IncompPropertiesInterface.hpp b/opm/core/fluid/IncompPropertiesInterface.hpp
index 5ac84fbb..f3c77caa 100644
--- a/opm/core/fluid/IncompPropertiesInterface.hpp
+++ b/opm/core/fluid/IncompPropertiesInterface.hpp
@@ -95,6 +95,18 @@ namespace Opm
                               const int* cells,
                               double* pc,
                               double* dpcds) const = 0;
+
+	/// Obtain the range of allowable saturation values.
+	/// In cell cells[i], saturation of phase p is allowed to be
+	/// in the interval [smin[i*P + p], smax[i*P + p]].
+        /// \param[in]  n      Number of data points.
+        /// \param[in]  cells  Array of n cell indices.
+        /// \param[out] smin   Array of nP minimum s values, array must be valid before calling.
+        /// \param[out] smax   Array of nP maximum s values, array must be valid before calling.
+        virtual void satRange(const int n,
+                              const int* cells,
+                              double* smin,
+                              double* smax) const = 0;
     };
 
 
diff --git a/opm/core/fluid/SaturationPropsBasic.cpp b/opm/core/fluid/SaturationPropsBasic.cpp
index f804af15..5bbdd939 100644
--- a/opm/core/fluid/SaturationPropsBasic.cpp
+++ b/opm/core/fluid/SaturationPropsBasic.cpp
@@ -197,6 +197,17 @@ namespace Opm
 
 
 
+    /// Obtain the range of allowable saturation values.
+    /// \param[in]  n      Number of data points.
+    /// \param[out] smin   Array of nP minimum s values, array must be valid before calling.
+    /// \param[out] smax   Array of nP maximum s values, array must be valid before calling.
+    void SaturationPropsBasic::satRange(const int n,
+					double* smin,
+					double* smax) const
+    {
+	std::fill(smin, smin + num_phases_*n, 0.0);
+	std::fill(smax, smax + num_phases_*n, 1.0);
+    }
 
 
 
diff --git a/opm/core/fluid/SaturationPropsBasic.hpp b/opm/core/fluid/SaturationPropsBasic.hpp
index 58ee3278..a00db3f2 100644
--- a/opm/core/fluid/SaturationPropsBasic.hpp
+++ b/opm/core/fluid/SaturationPropsBasic.hpp
@@ -76,6 +76,13 @@ namespace Opm
                       double* pc,
                       double* dpcds) const;
 
+	/// Obtain the range of allowable saturation values.
+        /// \param[in]  n      Number of data points.
+        /// \param[out] smin   Array of nP minimum s values, array must be valid before calling.
+        /// \param[out] smax   Array of nP maximum s values, array must be valid before calling.
+	void satRange(const int n,
+		      double* smin,
+		      double* smax) const;
     private:
 	enum RelPermFunc { Constant, Linear, Quadratic };
 	int num_phases_;
diff --git a/opm/core/fluid/SaturationPropsFromDeck.cpp b/opm/core/fluid/SaturationPropsFromDeck.cpp
index 2733feae..7882f926 100644
--- a/opm/core/fluid/SaturationPropsFromDeck.cpp
+++ b/opm/core/fluid/SaturationPropsFromDeck.cpp
@@ -42,6 +42,7 @@ namespace Opm
             THROW("SaturationPropsFromDeck::init()   --  oil phase must be active.");
         }
         const int samples = 200;
+	double swco = 0.0;
         if (phase_usage_.phase_used[Aqua]) {
             const SWOF::table_t& swof_table = deck.getSWOF().swof_;
             if (swof_table.size() != 1) {
@@ -55,6 +56,9 @@ namespace Opm
             buildUniformMonotoneTable(sw, krow, samples, krow_);
             buildUniformMonotoneTable(sw, pcow, samples, pcow_);
             krocw_ = krow[0]; // At connate water -> ecl. SWOF
+	    swco = sw[0];
+	    smin_[phase_usage_.phase_pos[Aqua]] = sw[0];
+	    smax_[phase_usage_.phase_pos[Aqua]] = sw.back();
         }
         if (phase_usage_.phase_used[Vapour]) {
             const SGOF::table_t& sgof_table = deck.getSGOF().sgof_;
@@ -68,7 +72,15 @@ namespace Opm
             buildUniformMonotoneTable(sg, krg,  samples, krg_);
             buildUniformMonotoneTable(sg, krog, samples, krog_);
             buildUniformMonotoneTable(sg, pcog, samples, pcog_);
+	    smin_[phase_usage_.phase_pos[Vapour]] = sg[0];
+	    if (std::fabs(sg.back() + swco - 1.0) > 1e-2) {
+		THROW("Gas maximum saturation in SGOF table = " << sg.back() <<
+		      ", should equal (1.0 - connate water sat) = " << (1.0 - swco));
+	    }
+	    smax_[phase_usage_.phase_pos[Vapour]] = sg.back();
         }
+	smin_[phase_usage_.phase_pos[Liquid]] = 0.0;
+	smax_[phase_usage_.phase_pos[Liquid]] = 1.0 - swco;
     }
 
 
@@ -142,6 +154,32 @@ namespace Opm
         }
     }
 
+
+
+
+    /// Obtain the range of allowable saturation values.
+    /// \param[in]  n      Number of data points.
+    /// \param[out] smin   Array of nP minimum s values, array must be valid before calling.
+    /// \param[out] smax   Array of nP maximum s values, array must be valid before calling.
+    void SaturationPropsFromDeck::satRange(const int n,
+					   double* smin,
+					   double* smax) const
+    {
+	const int np = phase_usage_.num_phases;
+	for (int i = 0; i < n; ++i) {
+	    for (int p = 0; p < np; ++p) {
+		smin[np*i + p] = smin_[p];
+		smax[np*i + p] = smax_[p];
+	    }
+	}
+    }
+
+
+
+
+    // Private methods below.
+
+
     void SaturationPropsFromDeck::evalKr(const double* s, double* kr) const
     {
         if (phase_usage_.num_phases == 3) {
diff --git a/opm/core/fluid/SaturationPropsFromDeck.hpp b/opm/core/fluid/SaturationPropsFromDeck.hpp
index b991a760..0e5d3b2d 100644
--- a/opm/core/fluid/SaturationPropsFromDeck.hpp
+++ b/opm/core/fluid/SaturationPropsFromDeck.hpp
@@ -67,6 +67,13 @@ namespace Opm
                       double* pc,
                       double* dpcds) const;
 
+	/// Obtain the range of allowable saturation values.
+        /// \param[in]  n      Number of data points.
+        /// \param[out] smin   Array of nP minimum s values, array must be valid before calling.
+        /// \param[out] smax   Array of nP maximum s values, array must be valid before calling.
+	void satRange(const int n,
+		      double* smin,
+		      double* smax) const;
     private:
         void evalKr(const double* s, double* kr) const;
         void evalKrDeriv(const double* s, double* kr, double* dkrds) const;
@@ -82,6 +89,8 @@ namespace Opm
         utils::UniformTableLinear<double> krog_;
         utils::UniformTableLinear<double> pcog_;
         double krocw_; // = krow_(s_wc)
+	double smin_[PhaseUsage::MaxNumPhases];
+	double smax_[PhaseUsage::MaxNumPhases];
     };