1 files changed, 36 insertions, 13 deletions

diff --git a/boost/geometry/strategy/cartesian/side_non_robust.hpp b/boost/geometry/strategy/cartesian/side_non_robust.hpp
index 2ef109cc1b..25074cc894 100644
--- a/boost/geometry/strategy/cartesian/side_non_robust.hpp
+++ b/boost/geometry/strategy/cartesian/side_non_robust.hpp
@@ -1,6 +1,6 @@
 // Boost.Geometry
 
-// Copyright (c) 2020, Oracle and/or its affiliates.
+// Copyright (c) 2020-2021, Oracle and/or its affiliates.
 
 // Contributed and/or modified by Vissarion Fysikopoulos, on behalf of Oracle
 
@@ -14,6 +14,8 @@
 #include <boost/geometry/util/select_calculation_type.hpp>
 #include <boost/geometry/util/precise_math.hpp>
 
+#include <boost/geometry/arithmetic/determinant.hpp>
+
 namespace boost { namespace geometry
 {
 
@@ -21,7 +23,7 @@ namespace strategy { namespace side
 {
 
 /*!
-\brief Adaptive precision predicate to check at which side of a segment a point lies:
+\brief Predicate to check at which side of a segment a point lies:
     left of segment (>0), right of segment (< 0), on segment (0).
 \ingroup strategies
 \tparam CalculationType \tparam_calculation
@@ -35,8 +37,6 @@ struct side_non_robust
 {
 public:
     //! \brief Computes double the signed area of the CCW triangle p1, p2, p
-
-#ifndef DOXYGEN_SHOULD_SKIP_THIS
     template
     <
         typename P1,
@@ -51,21 +51,44 @@ public:
                 P1,
                 P2,
                 P
-            >::type coordinate_type;
+            >::type CoordinateType;
         typedef typename select_most_precise
             <
-                coordinate_type,
+                CoordinateType,
                 double
-            >::type promoted_type;
+            >::type PromotedType;
+
+        CoordinateType const x = get<0>(p);
+        CoordinateType const y = get<1>(p);
+
+        CoordinateType const sx1 = get<0>(p1);
+        CoordinateType const sy1 = get<1>(p1);
+        CoordinateType const sx2 = get<0>(p2);
+        CoordinateType const sy2 = get<1>(p2);
+
+        //non-robust 1
+        //the following is 2x slower in some generic cases when compiled with g++
+        //(tested versions 9 and 10)
+        //
+        //auto detleft = (sx1 - x) * (sy2 - y);
+        //auto detright = (sy1 - y) * (sx2 - x);
+        //return detleft > detright ? 1 : (detleft < detright ? -1 : 0 );
+
+        //non-robust 2
+        PromotedType const dx = sx2 - sx1;
+        PromotedType const dy = sy2 - sy1;
+        PromotedType const dpx = x - sx1;
+        PromotedType const dpy = y - sy1;
 
-        auto detleft = (promoted_type(get<0>(p1)) - promoted_type(get<0>(p)))
-                * (promoted_type(get<1>(p2)) - promoted_type(get<1>(p)));
-        auto detright = (promoted_type(get<1>(p1)) - promoted_type(get<1>(p)))
-                * (promoted_type(get<0>(p2)) - promoted_type(get<0>(p)));
-        return detleft > detright ? 1 : (detleft < detright ? -1 : 0 );
+        PromotedType sv = geometry::detail::determinant<PromotedType>
+                (
+                    dx, dy,
+                    dpx, dpy
+                );
+        PromotedType const zero = PromotedType();
 
+        return sv == zero ? 0 : sv > zero ? 1 : -1;
     }
-#endif
 
 };