[Poincare/IEEE754] Changed methods to use std function

Implemented std::nextafter to replace hand made one.
Methods next and previous are no longer making the difference between -0
and +0

Change-Id: I42e1a073623b70656d9df954694803840cf3088c

liba/Makefile

@@ -70,6 +70,8 @@ liba_src += $(addprefix liba/src/external/openbsd/, \
   s_logbf.c \
   s_modf.c \
   s_modff.c \
+  s_nextafter.c \
+  s_nextafterf.c \
   s_rint.c \
   s_roundf.c \
   s_scalbnf.c \

liba/include/math.h

@@ -70,6 +70,7 @@ float log10f(float x);
 float logf(float x);
 float modff(float x, float *iptr);
 float nearbyintf(float x);
+float nextafterf(float from, float to);
 float powf(float x, float y);
 float roundf(float x);
 float scalbnf(float x, int exp);
@@ -110,6 +111,7 @@ double log2(double x);
 double logb(double x);
 double modf(double x, double *iptr);
 double nearbyint(double x);
+double nextafter(double from, double to);
 double pow(double x, double y);
 double rint(double x);
 double round(double x);
@@ -156,6 +158,7 @@ extern int signgam;
 #define modff(x, iptr) __builtin_modff(x, iptr)
 #define nanf(tagp) __builtin_nanf(tagp)
 #define nearbyintf(x) __builtin_nearbyintf(x)
+#define nextafterf(from, to) __builtin_nextafterf(from, to)
 #define powf(x, y) __builtin_powf(x, y)
 #define roundf(x) __builtin_roundf(x)
 #define scalbnf(x, exp) __builtin_scalbnf(x, exp)
@@ -196,6 +199,7 @@ extern int signgam;
 #define modf(x, iptr) __builtin_modf(x, iptr)
 #define nan(tagp) __builtin_nan(tagp)
 #define nearbyint(x) __builtin_nearbyint(x)
+#define nextafter(from, to) __builtin_nextafter(from, to)
 #define pow(x, y) __builtin_pow(x, y)
 #define rint(x) __builtin_rint(x)
 #define round(x) __builtin_round(x)

liba/src/external/openbsd/s_nextafter.c

@@ -0,0 +1,71 @@
+/* @(#)s_nextafter.c 5.1 93/09/24 */
+/*
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice 
+ * is preserved.
+ * ====================================================
+ */
+
+/* IEEE functions
+ *	nextafter(x,y)
+ *	return the next machine floating-point number of x in the
+ *	direction toward y.
+ *   Special cases:
+ */
+
+#include "math.h"
+#include "math_private.h"
+
+double
+nextafter(double x, double y)
+{
+	int32_t hx,hy,ix,iy;
+	u_int32_t lx,ly;
+
+	EXTRACT_WORDS(hx,lx,x);
+	EXTRACT_WORDS(hy,ly,y);
+	ix = hx&0x7fffffff;		/* |x| */
+	iy = hy&0x7fffffff;		/* |y| */
+
+	if(((ix>=0x7ff00000)&&((ix-0x7ff00000)|lx)!=0) ||   /* x is nan */
+	   ((iy>=0x7ff00000)&&((iy-0x7ff00000)|ly)!=0))     /* y is nan */
+	   return x+y;
+	if(x==y) return x;		/* x=y, return x */
+	if((ix|lx)==0) {			/* x == 0 */
+	    INSERT_WORDS(x,hy&0x80000000,1);	/* return +-minsubnormal */
+	    y = x*x;
+	    if(y==x) return y; else return x;	/* raise underflow flag */
+	}
+	if(hx>=0) {				/* x > 0 */
+	    if(hx>hy||((hx==hy)&&(lx>ly))) {	/* x > y, x -= ulp */
+		if(lx==0) hx -= 1;
+		lx -= 1;
+	    } else {				/* x < y, x += ulp */
+		lx += 1;
+		if(lx==0) hx += 1;
+	    }
+	} else {				/* x < 0 */
+	    if(hy>=0||hx>hy||((hx==hy)&&(lx>ly))){/* x < y, x -= ulp */
+		if(lx==0) hx -= 1;
+		lx -= 1;
+	    } else {				/* x > y, x += ulp */
+		lx += 1;
+		if(lx==0) hx += 1;
+	    }
+	}
+	hy = hx&0x7ff00000;
+	if(hy>=0x7ff00000) return x+x;	/* overflow  */
+	if(hy<0x00100000) {		/* underflow */
+	    y = x*x;
+	    if(y!=x) {		/* raise underflow flag */
+	        INSERT_WORDS(y,hx,lx);
+		return y;
+	    }
+	}
+	INSERT_WORDS(x,hx,lx);
+	return x;
+}

liba/src/external/openbsd/s_nextafterf.c

@@ -0,0 +1,62 @@
+/* s_nextafterf.c -- float version of s_nextafter.c.
+ * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
+ */
+
+/*
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice 
+ * is preserved.
+ * ====================================================
+ */
+
+#include "math.h"
+#include "math_private.h"
+
+float
+nextafterf(float x, float y)
+{
+	int32_t hx,hy,ix,iy;
+
+	GET_FLOAT_WORD(hx,x);
+	GET_FLOAT_WORD(hy,y);
+	ix = hx&0x7fffffff;		/* |x| */
+	iy = hy&0x7fffffff;		/* |y| */
+
+	if((ix>0x7f800000) ||   /* x is nan */
+	   (iy>0x7f800000))     /* y is nan */
+	   return x+y;
+	if(x==y) return x;		/* x=y, return x */
+	if(ix==0) {				/* x == 0 */
+	    SET_FLOAT_WORD(x,(hy&0x80000000)|1);/* return +-minsubnormal */
+	    y = x*x;
+	    if(y==x) return y; else return x;	/* raise underflow flag */
+	}
+	if(hx>=0) {				/* x > 0 */
+	    if(hx>hy) {				/* x > y, x -= ulp */
+		hx -= 1;
+	    } else {				/* x < y, x += ulp */
+		hx += 1;
+	    }
+	} else {				/* x < 0 */
+	    if(hy>=0||hx>hy){			/* x < y, x -= ulp */
+		hx -= 1;
+	    } else {				/* x > y, x += ulp */
+		hx += 1;
+	    }
+	}
+	hy = hx&0x7f800000;
+	if(hy>=0x7f800000) return x+x;	/* overflow  */
+	if(hy<0x00800000) {		/* underflow */
+	    y = x*x;
+	    if(y!=x) {		/* raise underflow flag */
+	        SET_FLOAT_WORD(y,hx);
+		return y;
+	    }
+	}
+	SET_FLOAT_WORD(x,hx);
+	return x;
+}

libaxx/include/cmath

@@ -39,6 +39,7 @@
 #undef modff
 #undef nanf
 #undef nearbyintf
+#undef nextafter
 #undef powf
 #undef roundf
 #undef scalbnf
@@ -81,6 +82,7 @@
 #undef modf
 #undef nan
 #undef nearbyint
+#undef nextafter
 #undef pow
 #undef rint
 #undef round
@@ -134,6 +136,7 @@ inline constexpr float logb(float x) { return __builtin_logbf(x); }
 inline constexpr float modf(float x, float *iptr) { return __builtin_modff(x, iptr); }
 inline constexpr float nanf(const char *tagp) { return __builtin_nanf(tagp); }
 inline constexpr float nearbyint(float x) { return __builtin_nearbyintf(x); }
+inline constexpr float nextafter(float from, float to) { return __builtin_nextafterf(from, to); }
 inline constexpr float pow(float x, float y) { return __builtin_powf(x, y); }
 inline constexpr float round(float x) { return __builtin_roundf(x); }
 inline constexpr float scalbn(float x, int exp) { return __builtin_scalbnf(x, exp); }
@@ -183,6 +186,7 @@ inline constexpr double logb(double x) { return __builtin_logb(x); }
 inline constexpr double modf(double x, double *iptr) { return __builtin_modf(x, iptr); }
 inline constexpr double nan(const char *tagp) { return __builtin_nan(tagp); }
 inline constexpr double nearbyint(double x) { return __builtin_nearbyint(x); }
+inline constexpr double nextafter(double from, double to) { return __builtin_nextafter(from, to); }
 inline constexpr double pow(double x, double y) { return __builtin_pow(x, y); }
 inline constexpr double rint(double x) { return __builtin_rint(x); }
 inline constexpr double round(double x) { return __builtin_round(x); }

poincare/Makefile

@@ -178,7 +178,6 @@ tests_src += $(addprefix poincare/test/,\
   function_solver.cpp\
   helper.cpp\
   helpers.cpp\
-  ieee754.cpp\
   integer.cpp\
   layout.cpp\
   layout_cursor.cpp\

poincare/include/poincare/ieee754.h

@@ -73,40 +73,12 @@ public:
     }
     return exponentBase10;
   }
-  static T next(T f) {
-    return nextOrPrevious(f, true);
-  }
-
-  static T previous(T f) {
-    return nextOrPrevious(f, false);
-  }
 
 private:
   union uint_float {
     uint64_t ui;
     T f;
   };
-  static T nextOrPrevious(T f, bool isNext) {
-    if (std::isinf(f) || std::isnan(f)) {
-      return f;
-    }
-    uint_float u;
-    u.ui = 0;
-    u.f = f;
-    uint64_t oneBitOnSignBit = (uint64_t)1 << (k_exponentNbBits + k_mantissaNbBits);
-    if ((isNext && (u.ui & oneBitOnSignBit) > 0) // next: Negative float
-     || (!isNext && (u.ui & oneBitOnSignBit) == 0)) { // previous: Positive float
-      if ((isNext && u.ui == oneBitOnSignBit) // next: -0.0
-       || (!isNext && u.ui == 0.0)) { // previous: 0.0
-        u.ui = isNext ? 0 : oneBitOnSignBit;
-      } else {
-        u.ui -= 1;
-      }
-    } else { // next: Positive float, previous: Negative float
-      u.ui += 1;
-    }
-    return u.f;
-  }
 
   constexpr static size_t k_signNbBits = 1;
   constexpr static size_t k_exponentNbBits = sizeof(T) == sizeof(float) ? 8 : 11;

poincare/src/solver.cpp

@@ -196,10 +196,14 @@ Coordinate2D<double> Solver::IncreasingFunctionRoot(double ax, double bx, double
      * representable double between min and max strictly. If there is, we choose
      * it instead, otherwise, we reached the most precise result possible. */
     if (currentAbscissa == min) {
-      currentAbscissa = IEEE754<double>::next(min);
+      if (currentAbscissa != -INFINITY) {
+        currentAbscissa = std::nextafter(currentAbscissa, INFINITY);
+      }
     }
     if (currentAbscissa == max) {
-      currentAbscissa = IEEE754<double>::previous(max);
+      if (currentAbscissa != INFINITY) {
+        currentAbscissa = std::nextafter(currentAbscissa, -INFINITY);
+      }
     }
     if (currentAbscissa == min || currentAbscissa == max) {
       break;

poincare/test/ieee754.cpp

@@ -1,40 +0,0 @@
-#include <poincare/ieee754.h>
-#include <string.h>
-#include <stdlib.h>
-#include <cmath>
-#include "helper.h"
-
-using namespace Poincare;
-
-template<typename T>
-void assert_next_and_previous_IEEE754_is(T a, T b) {
-  T next = IEEE754<T>::next(a);
-  T previous = IEEE754<T>::previous(b);
-  quiz_assert((std::isnan(next) && std::isnan(b)) || next == b);
-  quiz_assert((std::isnan(previous) && std::isnan(a)) || previous == a);
-}
-
-QUIZ_CASE(ieee754_next_and_previous) {
-  assert_next_and_previous_IEEE754_is<float>(0.0f, 1.4E-45f);
-  assert_next_and_previous_IEEE754_is<float>(INFINITY, INFINITY);
-  assert_next_and_previous_IEEE754_is<float>(NAN, NAN);
-  assert_next_and_previous_IEEE754_is<float>(2.2566837E-10f, 2.2566839E-10f);
-  assert_next_and_previous_IEEE754_is<float>(-INFINITY, -INFINITY);
-  assert_next_and_previous_IEEE754_is<float>(-0.0f, 0.0f);
-  assert_next_and_previous_IEEE754_is<float>(-1.4E-45f, -0.0f);
-  assert_next_and_previous_IEEE754_is<float>(-3.4359738E10f, -3.43597363E10f);
-  quiz_assert(IEEE754<float>::next(3.4028235E38f) == INFINITY);
-  quiz_assert(IEEE754<float>::previous(-3.4028235E38f) == -INFINITY);
-
-  assert_next_and_previous_IEEE754_is<double>(0.0f, 4.94065645841246544176568792868E-324);
-  assert_next_and_previous_IEEE754_is<double>(INFINITY, INFINITY);
-  assert_next_and_previous_IEEE754_is<double>(NAN, NAN);
-  assert_next_and_previous_IEEE754_is<double>(1.936766735060658315512927142E-282, 1.93676673506065869777770528092E-282);
-  assert_next_and_previous_IEEE754_is<double>(-INFINITY, -INFINITY);
-  assert_next_and_previous_IEEE754_is<double>(-0.0, 0.0);
-  assert_next_and_previous_IEEE754_is<double>(-4.94065645841246544176568792868E-324, -0.0);
-  assert_next_and_previous_IEEE754_is<double>(-1.38737906372912431085182213247E201, -1.38737906372912403890916981028E201);
-  quiz_assert(IEEE754<double>::next(1.79769313486231570814527423732E308) == INFINITY);
-  quiz_assert(IEEE754<double>::previous(-1.79769313486231570814527423732E308) == -INFINITY);
-
-}