First set of licensing changes
[mono.git] / mono / metadata / sysmath.c
index 6627009eebb98644ad62238f196cb4bbb7d41e1b..8d55a552f84dbc86d35b01deea7793039e7b007d 100644 (file)
-/* math.c - these are based on bob smith's csharp routines */
+/*
+ * sysmath.c: these are based on bob smith's csharp routines 
+ *
+ * Author:
+ *     Mono Project (http://www.mono-project.com)
+ *     Ludovic Henry (ludovic@xamarin.com)
+ *
+ * Copyright 2001-2003 Ximian, Inc (http://www.ximian.com)
+ * Copyright 2004-2009 Novell, Inc (http://www.novell.com)
+ * Copyright 2015 Xamarin, Inc (https://www.xamarin.com)
+ * Licensed under the MIT license. See LICENSE file in the project root for full license information.
+ */
+
+//
+// Copyright (c) Microsoft. All rights reserved.
+// Licensed under the MIT license. See LICENSE file in the project root for full license information.
+//
+// Files:
+//  - src/classlibnative/float/floatnative.cpp
+//  - src/pal/src/cruntime/floatnative.cpp
+//
+// Ported from C++ to C and adjusted to Mono runtime
 
 #define __USE_ISOC99
+
 #include <math.h>
 #include <mono/metadata/sysmath.h>
-#include <mono/metadata/exception.h>
-
-#ifndef NAN
-# if G_BYTE_ORDER == G_BIG_ENDIAN
-#  define __nan_bytes           { 0x7f, 0xc0, 0, 0 }
-# endif
-# if G_BYTE_ORDER == G_LITTLE_ENDIAN
-#  define __nan_bytes           { 0, 0, 0xc0, 0x7f }
-# endif
-
-static union { unsigned char __c[4]; float __d; } __nan_union = { __nan_bytes };
-# define NAN    (__nan_union.__d)
-#endif
 
-#ifndef HUGE_VAL
-#define __huge_val_t   union { unsigned char __c[8]; double __d; }
-# if G_BYTE_ORDER == G_BIG_ENDIAN
-#  define __HUGE_VAL_bytes       { 0x7f, 0xf0, 0, 0, 0, 0, 0, 0 }
-# endif
-# if G_BYTE_ORDER == G_LITTLE_ENDIAN
-#  define __HUGE_VAL_bytes       { 0, 0, 0, 0, 0, 0, 0xf0, 0x7f }
-# endif
-static __huge_val_t __huge_val = { __HUGE_VAL_bytes };
-#  define HUGE_VAL      (__huge_val.__d)
-#endif
+#include "number-ms.h"
+#include "utils/mono-compiler.h"
+
+static const MonoDouble_double NaN = { .s = { .sign = 0x0, .exp = 0x7FF, .mantHi = 0x80000, .mantLo = 0x0 } };
+
+/* +Infinity */
+static const MonoDouble_double PInfinity = { .s = { .sign = 0x0, .exp = 0x7FF, .mantHi = 0x0, .mantLo = 0x0 } };
+
+/* -Infinity */
+static const MonoDouble_double MInfinity = { .s = { .sign = 0x1, .exp = 0x7FF, .mantHi = 0x0, .mantLo = 0x0 } };
+
+/* +1 */
+static const MonoDouble_double POne = { .s = { .sign = 0x0, .exp = 0x3FF, .mantHi = 0x0, .mantLo = 0x0 } };
+
+/* -1 */
+static const MonoDouble_double MOne = { .s = { .sign = 0x1, .exp = 0x3FF, .mantHi = 0x0, .mantLo = 0x0 } };
+
+static MONO_ALWAYS_INLINE gboolean
+isplusinfinity (gdouble d)
+{
+       return d == PInfinity.d;
+}
+
+static MONO_ALWAYS_INLINE gboolean
+isminusinfinity (gdouble d)
+{
+       return d == MInfinity.d;
+}
+
+static MONO_ALWAYS_INLINE gboolean
+isinfinity (gdouble d)
+{
+       return isplusinfinity (d) || isminusinfinity (d);
+}
+
+static MONO_ALWAYS_INLINE gboolean
+isplusone (gdouble d)
+{
+       return d == POne.d;
+}
 
+static MONO_ALWAYS_INLINE gboolean
+isminusone (gdouble d)
+{
+       return d == MOne.d;
+}
 
-gdouble ves_icall_System_Math_Floor (gdouble x) {
-       MONO_ARCH_SAVE_REGS;
+gdouble
+ves_icall_System_Math_Floor (gdouble x)
+{
        return floor(x);
 }
 
-gdouble ves_icall_System_Math_Round (gdouble x) {
-       double int_part, dec_part;
-       MONO_ARCH_SAVE_REGS;
-       int_part = floor(x);
-       dec_part = x - int_part;
-       if (((dec_part == 0.5) &&
-               ((2.0 * ((int_part / 2.0) - floor(int_part / 2.0))) != 0.0)) ||
-               (dec_part > 0.5)) {
-               int_part++;
+gdouble
+ves_icall_System_Math_Round (gdouble x)
+{
+       gdouble tmp, floor_tmp;
+
+       /* If the number has no fractional part do nothing This shortcut is necessary
+        * to workaround precision loss in borderline cases on some platforms */
+       if (x == (gdouble)(gint64) x)
+               return x;
+
+       tmp = x + 0.5;
+       floor_tmp = floor (tmp);
+
+       if (floor_tmp == tmp) {
+               if (fmod (tmp, 2.0) != 0)
+                       floor_tmp -= 1.0;
        }
-       return int_part;
-}
-
-gdouble ves_icall_System_Math_Round2 (gdouble value, gint32 digits) {
-       double p, int_part, dec_part;
-       MONO_ARCH_SAVE_REGS;
-       if (value == HUGE_VAL)
-               return HUGE_VAL;
-       if (value == -HUGE_VAL)
-               return -HUGE_VAL;
-       if (digits == 0)
-               return ves_icall_System_Math_Round(value);
-       p = pow(10, digits);
-       dec_part = modf (value, &int_part);
-       dec_part *= 1000000000000000ULL;
-       dec_part = floor(dec_part);
-       dec_part /= (1000000000000000ULL / p);
-       dec_part = ves_icall_System_Math_Round(dec_part);
-       dec_part /= p;
-       return int_part + dec_part;
+
+       return copysign (floor_tmp, x);
 }
 
 gdouble 
 ves_icall_System_Math_Sin (gdouble x)
 {
-       MONO_ARCH_SAVE_REGS;
-
        return sin (x);
 }
 
 gdouble 
 ves_icall_System_Math_Cos (gdouble x)
 {
-       MONO_ARCH_SAVE_REGS;
-
        return cos (x);
 }
 
 gdouble 
 ves_icall_System_Math_Tan (gdouble x)
 {
-       MONO_ARCH_SAVE_REGS;
-
        return tan (x);
 }
 
 gdouble 
 ves_icall_System_Math_Sinh (gdouble x)
 {
-       MONO_ARCH_SAVE_REGS;
-
        return sinh (x);
 }
 
 gdouble 
 ves_icall_System_Math_Cosh (gdouble x)
 {
-       MONO_ARCH_SAVE_REGS;
-
        return cosh (x);
 }
 
 gdouble 
 ves_icall_System_Math_Tanh (gdouble x)
 {
-       MONO_ARCH_SAVE_REGS;
-
        return tanh (x);
 }
 
 gdouble 
 ves_icall_System_Math_Acos (gdouble x)
 {
-       MONO_ARCH_SAVE_REGS;
-
        if (x < -1 || x > 1)
-               return NAN;
+               return NaN.d;
 
        return acos (x);
 }
@@ -129,10 +148,8 @@ ves_icall_System_Math_Acos (gdouble x)
 gdouble 
 ves_icall_System_Math_Asin (gdouble x)
 {
-       MONO_ARCH_SAVE_REGS;
-
        if (x < -1 || x > 1)
-               return NAN;
+               return NaN.d;
 
        return asin (x);
 }
@@ -140,31 +157,26 @@ ves_icall_System_Math_Asin (gdouble x)
 gdouble 
 ves_icall_System_Math_Atan (gdouble x)
 {
-       MONO_ARCH_SAVE_REGS;
-
        return atan (x);
 }
 
 gdouble 
 ves_icall_System_Math_Atan2 (gdouble y, gdouble x)
 {
-       double result;
-       MONO_ARCH_SAVE_REGS;
+       gdouble result;
+
+       if (isinfinity (x) && isinfinity (y))
+               return NaN.d;
 
-       if ((y == HUGE_VAL && x == HUGE_VAL) ||
-               (y == HUGE_VAL && x == -HUGE_VAL) ||
-               (y == -HUGE_VAL && x == HUGE_VAL) ||
-               (y == -HUGE_VAL && x == -HUGE_VAL)) {
-               return NAN;
-       }
        result = atan2 (y, x);
-       return (result == -0)? 0: result;
+       return result == -0.0 ? 0.0: result;
 }
 
 gdouble 
 ves_icall_System_Math_Exp (gdouble x)
 {
-       MONO_ARCH_SAVE_REGS;
+       if (isinfinity (x))
+               return x < 0 ? 0.0 : x;
 
        return exp (x);
 }
@@ -172,12 +184,10 @@ ves_icall_System_Math_Exp (gdouble x)
 gdouble 
 ves_icall_System_Math_Log (gdouble x)
 {
-       MONO_ARCH_SAVE_REGS;
-
        if (x == 0)
-               return -HUGE_VAL;
+               return MInfinity.d;
        else if (x < 0)
-               return NAN;
+               return NaN.d;
 
        return log (x);
 }
@@ -185,12 +195,10 @@ ves_icall_System_Math_Log (gdouble x)
 gdouble 
 ves_icall_System_Math_Log10 (gdouble x)
 {
-       MONO_ARCH_SAVE_REGS;
-
        if (x == 0)
-               return -HUGE_VAL;
+               return MInfinity.d;
        else if (x < 0)
-               return NAN;
+               return NaN.d;
 
        return log10 (x);
 }
@@ -198,46 +206,112 @@ ves_icall_System_Math_Log10 (gdouble x)
 gdouble 
 ves_icall_System_Math_Pow (gdouble x, gdouble y)
 {
-       double result;
-       MONO_ARCH_SAVE_REGS;
-
-       if (isnan(x) || isnan(y)) {
-               return NAN;
+       gdouble result;
+
+       if (isnan (y))
+               return y;
+       if (isnan (x))
+               return x;
+
+       if (isinfinity (y)) {
+               if (isplusone (x))
+                       return x;
+               if (isminusone (x))
+                       return NaN.d;
        }
 
-       if ((x == 1 || x == -1) && (y == HUGE_VAL || y == -HUGE_VAL)) {
-               return NAN;
+       /* following are cases from PAL_pow which abstract the implementation of pow for posix and win32 platforms
+        * (https://github.com/dotnet/coreclr/blob/master/src/pal/src/cruntime/finite.cpp#L331) */
+
+       if (isplusinfinity (y) && !isnan (x)) {
+               if (isplusone (x) || isminusone (x))
+                       result = NaN.d;
+               else if (x > MOne.d && x < POne.d)
+                       result = 0.0;
+               else
+                       result = PInfinity.d;
+       } else if (isminusinfinity (y) && !isnan (x)) {
+               if (isplusone (x) || isminusone (x))
+                       result = NaN.d;
+               if (x > MOne.d && x < POne.d)
+                       result = PInfinity.d;
+               else
+                       result = 0.0;
+       } else if (x == 0.0 && y < 0.0) {
+               result = PInfinity.d;
+       } else if (y == 0.0 && isnan (x)) {
+               /* Windows returns NaN for pow(NaN, 0), but POSIX specifies
+                * a return value of 1 for that case.  We need to return
+                * the same result as Windows. */
+               result = NaN.d;
+       } else {
+               result = pow (x, y);
        }
 
-       /* This code is for return the same results as MS.NET for certain
-        * limit values */
-       if (x < -9007199254740991.0) {
-               if (y > 9007199254740991.0)
-                       return HUGE_VAL;
-               if (y < -9007199254740991.0)
-                       return 0;
-       }
-
-       result = pow (x, y);
-
-       /* This code is for return the same results as MS.NET for certain
-        * limit values */
-       if (isnan(result) &&
-               (x == -1.0) &&
-               ((y > 9007199254740991.0) || (y < -9007199254740991.0))) {
-               return 1;
+       if (result == PInfinity.d && x < 0.0 && isfinite (x) && ceil (y / 2) != floor (y / 2))
+               result = MInfinity.d;
+
+       /*
+        * The even/odd test in the if (this one and the one above) used to be ((long long) y % 2 == 0)
+        * on SPARC (long long) y for large y (>2**63) is always 0x7fffffff7fffffff, which
+        * is an odd number, so the test ((long long) y % 2 == 0) will always fail for
+        * large y. Since large double numbers are always even (e.g., the representation of
+        * 1E20+1 is the same as that of 1E20, the last .+1. is too insignificant to be part
+        * of the representation), this test will always return the wrong result for large y.
+        *
+        * The (ceil(y/2) == floor(y/2)) test is slower, but more robust.
+        */
+       if (result == MInfinity.d && x < 0.0 && isfinite (x) && ceil (y / 2) == floor (y / 2))
+               result = PInfinity.d;
+
+#if defined (__linux__) && SIZEOF_VOID_P == 4
+       /* On Linux 32bits, some tests erroneously return NaN */
+       if (isnan (result)) {
+               if (isminusone (x) && (y > 9007199254740991.0 || y < -9007199254740991.0)) {
+                       /* Math.Pow (-1, Double.MaxValue) and Math.Pow (-1, Double.MinValue) should return 1 */
+                       result = POne.d;
+               } else if (x < -9007199254740991.0 && y < -9007199254740991.0) {
+                       /* Math.Pow (Double.MinValue, Double.MinValue) should return 0 */
+                       result = 0.0;
+               } else if (x < -9007199254740991.0 && y > 9007199254740991.0) {
+                       /* Math.Pow (Double.MinValue, Double.MaxValue) should return Double.PositiveInfinity */
+                       result = PInfinity.d;
+               }
        }
+#endif
 
-       return (result == -0)? 0: result;
+       return result == -0.0 ? 0 : result;
 }
 
 gdouble 
 ves_icall_System_Math_Sqrt (gdouble x)
 {
-       MONO_ARCH_SAVE_REGS;
-
        if (x < 0)
-               return NAN;
+               return NaN.d;
 
        return sqrt (x);
 }
+
+gdouble
+ves_icall_System_Math_Abs_double (gdouble v)
+{
+       return fabs (v);
+}
+
+gfloat
+ves_icall_System_Math_Abs_single (gfloat v)
+{
+       return fabsf (v);
+}
+
+gdouble
+ves_icall_System_Math_Ceiling (gdouble v)
+{
+       return ceil (v);
+}
+
+gdouble
+ves_icall_System_Math_SplitFractionDouble (gdouble *v)
+{
+       return modf (*v, v);
+}