math: improve POSIX compliance

include/math.h

@@ -25,6 +25,10 @@
 extern "C" {
 #endif
 
+#define MATH_ERRNO       1
+#define MATH_ERREXCEPT   2
+#define math_errhandling (MATH_ERRNO)
+
 #define FP_NAN       0
 #define FP_INFINITE  1
 #define FP_ZERO      2

math/common.c

@@ -14,32 +14,36 @@
  */
 
 #include <stdint.h>
+#include <errno.h>
+#include <math.h>
+#include <float.h>
 #include "common.h"
 
 
 void normalizeSub(double *x, int *exp)
 {
 	conv_t *conv = (conv_t *)x;
 
-	if (conv->i.mantisa == 0)
+	if (conv->i.mantisa == 0) {
 		return;
+	}
 
-	while (!(conv->i.mantisa & 0xffffff0000000LL)) {
+	while ((conv->i.mantisa & 0xffffff0000000LL) == 0) {
 		conv->i.mantisa <<= 24;
 		*exp -= 24;
 	}
 
-	while (!(conv->i.mantisa & 0xff00000000000LL)) {
+	while ((conv->i.mantisa & 0xff00000000000LL) == 0) {
 		conv->i.mantisa <<= 8;
 		*exp -= 8;
 	}
 
-	while (!(conv->i.mantisa & 0xf000000000000LL)) {
+	while ((conv->i.mantisa & 0xf000000000000LL) == 0) {
 		conv->i.mantisa <<= 4;
 		*exp -= 4;
 	}
 
-	while (!(conv->i.mantisa & 0x8000000000000LL)) {
+	while ((conv->i.mantisa & 0x8000000000000LL) == 0) {
 		conv->i.mantisa <<= 1;
 		*exp -= 1;
 	}
@@ -89,8 +93,9 @@ void createSub(double *x, int exp)
 		++exp;
 	}
 
-	if (conv->i.mantisa == 0)
+	if (conv->i.mantisa == 0) {
 		*x = 0.0;
+	}
 }
 
 
@@ -99,34 +104,43 @@ double quickPow(double x, int e)
 	double res;
 	unsigned int eabs;
 
-	if (e == 0)
+	if (e == 0) {
 		return 1.0;
-	else if (e == 1)
+	}
+	else if (e == 1) {
 		return x;
+	}
 
 	res = 1.0;
 
 	if (e < 0) {
-		eabs = -e;
+		eabs = (unsigned int)(-e);
 
-		while (eabs != 0) {
-			if (eabs & 1)
+		while (eabs != 0u) {
+			if (eabs & 1u) {
 				res /= x;
+			}
 			x *= x;
 			eabs >>= 1;
 		}
 	}
 	else {
-		eabs = e;
+		eabs = (unsigned int)e;
 
-		while (eabs != 0) {
-			if (eabs & 1)
+		while (eabs != 0u) {
+			if (eabs & 1u) {
 				res *= x;
+			}
 			x *= x;
 			eabs >>= 1;
 		}
 	}
 
+	/* ia32 use 80-bit extended precision registers so isinf() may not be true */
+	if ((isinf(res) != 0) || (fabs(res) > DBL_MAX)) {
+		errno = ERANGE;
+	}
+
 	return res;
 }
 
@@ -135,15 +149,17 @@ int isInteger(double x)
 {
 	conv_t *conv = (conv_t *)&x;
 	int exp = conv->i.exponent - 1023;
-	uint64_t mask = 0xfffffffffffffLL;
+	uint64_t mask = 0xfffffffffffffULL;
 	uint64_t m;
 
-	if (exp > 52)
+	if (exp > 52) {
 		return 1;
-	else if (exp < 0)
+	}
+	else if (exp < 0) {
 		return 0;
+	}
 
 	m = conv->i.mantisa & ~(mask >> exp);
 
-	return m == conv->i.mantisa;
+	return (m == conv->i.mantisa);
 }

math/common.h

@@ -20,20 +20,20 @@
 
 
 typedef union {
-		struct {
+	struct {
 #if __FLOAT_WORD_ORDER__ == __ORDER_LITTLE_ENDIAN__
-			uint64_t mantisa:52;
-			uint16_t exponent:11;
-			uint8_t sign:1;
+		uint64_t mantisa : 52;
+		uint16_t exponent : 11;
+		uint8_t sign : 1;
 #elif __FLOAT_WORD_ORDER__ == __ORDER_BIG_ENDIAN__
-			uint8_t sign:1;
-			uint16_t exponent:11;
-			uint64_t mantisa:52;
+		uint8_t sign : 1;
+		uint16_t exponent : 11;
+		uint64_t mantisa : 52;
 #else
 #error "Unsupported byte order"
 #endif
-		} i;
-		double d;
+	} i;
+	double d;
 } conv_t;