diff --git a/Cartesian_kernel/include/CGAL/Cartesian/solve_3.h b/Cartesian_kernel/include/CGAL/Cartesian/solve_3.h
index 2c9137f1596..81d1f422a2d 100644
--- a/Cartesian_kernel/include/CGAL/Cartesian/solve_3.h
+++ b/Cartesian_kernel/include/CGAL/Cartesian/solve_3.h
@@ -29,13 +29,13 @@ void solve (const VectorC3<R> &v0,
             const VectorC3<R> &v1,
             const VectorC3<R> &v2,
             const VectorC3<R> &d,
-            typename R::FT &alpha, typename R::FT &beta, typename R::FT &gamma)
+            typename R::FT &alpha, typename R::FT &beta, typename R::FT &gamma, typename R::FT &denom)
 {
   CGAL::solve(v0.x(), v0.y(), v0.z(),
-        v1.x(), v1.y(), v1.z(),
-        v2.x(), v2.y(), v2.z(),
-        d.x(),  d.y(),  d.z(),
-        alpha, beta, gamma);
+              v1.x(), v1.y(), v1.z(),
+              v2.x(), v2.y(), v2.z(),
+              d.x(),  d.y(),  d.z(),
+              alpha, beta, gamma, denom);
 }
 
 } // namespace Cartesian_internal
diff --git a/Kernel_23/include/CGAL/Kernel/solve.h b/Kernel_23/include/CGAL/Kernel/solve.h
index f4866259898..7b915836523 100644
--- a/Kernel_23/include/CGAL/Kernel/solve.h
+++ b/Kernel_23/include/CGAL/Kernel/solve.h
@@ -58,6 +58,37 @@ void solve (const FT &a1, const FT &a2, const FT &a3,
 }
 
 
+template <class FT>
+void solve (const FT &a1, const FT &a2, const FT &a3,
+            const FT &b1, const FT &b2, const FT &b3,
+            const FT &c1, const FT &c2, const FT &c3,
+            const FT &d1, const FT &d2, const FT &d3,
+            FT &x, FT &y, FT &z, FT& denom)
+{
+  FT ab23 = a3*b2 - a2*b3;
+  FT ab13 = a3*b1 - a1*b3;
+  FT ab12 = a2*b1 - a1*b2;
+
+  denom = ab23*c1 - ab13*c2 + ab12*c3;
+
+  FT cd23 = c3*d2 - c2*d3;
+  FT cd13 = c3*d1 - c1*d3;
+  FT cd12 = c2*d1 - c1*d2;
+
+  x = (b3*cd12 - b2*cd13 + b1*cd23);
+
+  y = (a2*cd13 - cd12*a3 - cd23*a1);
+
+  z = (ab23*d1 + ab12*d3 - ab13*d2);
+  if(denom < 0){
+    denom = -denom;
+    x = -x;
+    y = -y;
+    z = -z;
+  }
+}
+
+
 // this is for a parabola c1, c2, c3 are equal to 1
 template <class FT>
 void solve_quadratic (const FT &a1, const FT &a2, const FT &a3,