From 72548ecb2d2e15b2fb7dd19a970f92d3e9e28b68 Mon Sep 17 00:00:00 2001
From: Marc Glisse <marc.glisse@inria.fr>
Date: Wed, 30 Oct 2019 17:27:00 +0100
Subject: [PATCH] More NaN fixes for the non-SSE path.

---
 Number_types/include/CGAL/Interval_nt.h       |  7 +++++--
 .../test/Number_types/Interval_nt_new.cpp     | 19 +++++++++++++++++++
 2 files changed, 24 insertions(+), 2 deletions(-)

diff --git a/Number_types/include/CGAL/Interval_nt.h b/Number_types/include/CGAL/Interval_nt.h
index b8908f2b05b..e19e1668d04 100644
--- a/Number_types/include/CGAL/Interval_nt.h
+++ b/Number_types/include/CGAL/Interval_nt.h
@@ -639,6 +639,7 @@ private:
     return IA (IA_opacify128(r));
 # endif
 #else
+    // TODO: try to move some NaN tests out of the hot path (test a.inf()>0 instead of >=0?).
     if (a.inf() >= 0.0)					// a>=0
     {
       // b>=0     [a.inf()*b.inf(); a.sup()*b.sup()]
@@ -652,7 +653,8 @@ private:
         if (b.sup() < 0.0)
           bb = a.inf();
       }
-      return IA(-CGAL_IA_MUL(aa, -b.inf()), CGAL_IA_MUL(bb, b.sup()));
+      double r = (b.sup() == 0) ? 0. : CGAL_IA_MUL(bb, b.sup()); // In case bb is infinite, avoid NaN.
+      return IA(-CGAL_IA_MUL(aa, -b.inf()), r);
     }
     else if (a.sup()<=0.0)				// a<=0
     {
@@ -663,9 +665,10 @@ private:
       if (b.inf() < 0.0)
       {
         aa=bb;
-        if (b.sup() < 0.0)
+        if (b.sup() <= 0.0)
           bb=a.sup();
       }
+      else if (b.sup() <= 0) return 0.; // In case a has an infinite bound, avoid NaN.
       return IA(-CGAL_IA_MUL(-bb, b.sup()), CGAL_IA_MUL(-aa, -b.inf()));
     }
     else						// 0 \in a
diff --git a/Number_types/test/Number_types/Interval_nt_new.cpp b/Number_types/test/Number_types/Interval_nt_new.cpp
index 291e4e41d23..95d4f5521fe 100644
--- a/Number_types/test/Number_types/Interval_nt_new.cpp
+++ b/Number_types/test/Number_types/Interval_nt_new.cpp
@@ -170,12 +170,31 @@ int main() {
         Interval d(0,inf);
         // For CGAL's purposes, [0,0]*anything is 0, but we tolerate larger intervals if they can be produced faster.
         for (Interval I : { all*zero, zero*all, all*0., 0.*all, b*zero, zero*b, -b*zero, zero*-b, c*zero, zero*c, -c*zero, zero*-c })
+        {
+          //std::cout << I << '\n';
           assert(I.inf()<=0 && I.sup()>=0);
+        }
         // This should be [0,inf], but again we tolerate more.
         for (Interval I : { a*b, b*a, -a*-b, -b*-a, d*d, -d*-d })
+        {
+          //std::cout << I << '\n';
           assert(I.inf()<=0 && I.sup()==inf);
+        }
+        for (Interval I : { -a*b, -b*a, a*-b, b*-a, -d*d, d*-d })
+        {
+          //std::cout << I << '\n';
+          assert(I.inf()==-inf && I.sup()>=0);
+        }
         for (Interval I : { all*a, a*all, all*-a, -a*all })
+        {
+          //std::cout << I << '\n';
           assert(I.inf()==-inf && I.sup()==inf);
+        }
+        for (Interval I : { all*d, d*all, all*-d, -d*all })
+        {
+          //std::cout << I << '\n';
+          assert(I.inf()==-inf && I.sup()==inf);
+        }
     }
     {// external functions on Intervals
      // functions (abs, square, sqrt, pow)