mono/metadata/sysmath.c

   1 /*
   2  * sysmath.c: these are based on bob smith's csharp routines
   3  *
   4  * Author:
   5  *      Mono Project (http://www.mono-project.com)
   6  *      Ludovic Henry (ludovic@xamarin.com)
   7  *
   8  * Copyright 2001-2003 Ximian, Inc (http://www.ximian.com)
   9  * Copyright 2004-2009 Novell, Inc (http://www.novell.com)
  10  * Copyright 2015 Xamarin, Inc (https://www.xamarin.com)
  11  * Licensed under the MIT license. See LICENSE file in the project root for full license information.
  12  */
  13
  14 //
  15 // Copyright (c) Microsoft. All rights reserved.
  16 // Licensed under the MIT license. See LICENSE file in the project root for full license information.
  17 //
  18 // Files:
  19 //  - src/classlibnative/float/floatnative.cpp
  20 //  - src/pal/src/cruntime/floatnative.cpp
  21 //
  22 // Ported from C++ to C and adjusted to Mono runtime
  23
  24 #define __USE_ISOC99
  25
  26 #include <math.h>
  27 #include <mono/metadata/sysmath.h>
  28
  29 #include "number-ms.h"
  30 #include "utils/mono-compiler.h"
  31
  32 static const MonoDouble_double NaN = { .s = { .sign = 0x0, .exp = 0x7FF, .mantHi = 0x80000, .mantLo = 0x0 } };
  33
  34 /* +Infinity */
  35 static const MonoDouble_double PInfinity = { .s = { .sign = 0x0, .exp = 0x7FF, .mantHi = 0x0, .mantLo = 0x0 } };
  36
  37 /* -Infinity */
  38 static const MonoDouble_double MInfinity = { .s = { .sign = 0x1, .exp = 0x7FF, .mantHi = 0x0, .mantLo = 0x0 } };
  39
  40 /* +1 */
  41 static const MonoDouble_double POne = { .s = { .sign = 0x0, .exp = 0x3FF, .mantHi = 0x0, .mantLo = 0x0 } };
  42
  43 /* -1 */
  44 static const MonoDouble_double MOne = { .s = { .sign = 0x1, .exp = 0x3FF, .mantHi = 0x0, .mantLo = 0x0 } };
  45
  46 static MONO_ALWAYS_INLINE gboolean
  47 isplusinfinity (gdouble d)
  48 {
  49         return d == PInfinity.d;
  50 }
  51
  52 static MONO_ALWAYS_INLINE gboolean
  53 isminusinfinity (gdouble d)
  54 {
  55         return d == MInfinity.d;
  56 }
  57
  58 static MONO_ALWAYS_INLINE gboolean
  59 isinfinity (gdouble d)
  60 {
  61         return isplusinfinity (d) || isminusinfinity (d);
  62 }
  63
  64 static MONO_ALWAYS_INLINE gboolean
  65 isplusone (gdouble d)
  66 {
  67         return d == POne.d;
  68 }
  69
  70 static MONO_ALWAYS_INLINE gboolean
  71 isminusone (gdouble d)
  72 {
  73         return d == MOne.d;
  74 }
  75
  76 gdouble
  77 ves_icall_System_Math_Floor (gdouble x)
  78 {
  79         return floor(x);
  80 }
  81
  82 gdouble
  83 ves_icall_System_Math_Round (gdouble x)
  84 {
  85         gdouble tmp, floor_tmp;
  86
  87         /* If the number has no fractional part do nothing This shortcut is necessary
  88          * to workaround precision loss in borderline cases on some platforms */
  89         if (x == (gdouble)(gint64) x)
  90                 return x;
  91
  92         tmp = x + 0.5;
  93         floor_tmp = floor (tmp);
  94
  95         if (floor_tmp == tmp) {
  96                 if (fmod (tmp, 2.0) != 0)
  97                         floor_tmp -= 1.0;
  98         }
  99
 100         return copysign (floor_tmp, x);
 101 }
 102
 103 gdouble
 104 ves_icall_System_Math_Sin (gdouble x)
 105 {
 106         return sin (x);
 107 }
 108
 109 gdouble
 110 ves_icall_System_Math_Cos (gdouble x)
 111 {
 112         return cos (x);
 113 }
 114
 115 gdouble
 116 ves_icall_System_Math_Tan (gdouble x)
 117 {
 118         return tan (x);
 119 }
 120
 121 gdouble
 122 ves_icall_System_Math_Sinh (gdouble x)
 123 {
 124         return sinh (x);
 125 }
 126
 127 gdouble
 128 ves_icall_System_Math_Cosh (gdouble x)
 129 {
 130         return cosh (x);
 131 }
 132
 133 gdouble
 134 ves_icall_System_Math_Tanh (gdouble x)
 135 {
 136         return tanh (x);
 137 }
 138
 139 gdouble
 140 ves_icall_System_Math_Acos (gdouble x)
 141 {
 142         if (x < -1 || x > 1)
 143                 return NaN.d;
 144
 145         return acos (x);
 146 }
 147
 148 gdouble
 149 ves_icall_System_Math_Asin (gdouble x)
 150 {
 151         if (x < -1 || x > 1)
 152                 return NaN.d;
 153
 154         return asin (x);
 155 }
 156
 157 gdouble
 158 ves_icall_System_Math_Atan (gdouble x)
 159 {
 160         return atan (x);
 161 }
 162
 163 gdouble
 164 ves_icall_System_Math_Atan2 (gdouble y, gdouble x)
 165 {
 166         gdouble result;
 167
 168         if (isinfinity (x) && isinfinity (y))
 169                 return NaN.d;
 170
 171         result = atan2 (y, x);
 172         return result == -0.0 ? 0.0: result;
 173 }
 174
 175 gdouble
 176 ves_icall_System_Math_Exp (gdouble x)
 177 {
 178         if (isinfinity (x))
 179                 return x < 0 ? 0.0 : x;
 180
 181         return exp (x);
 182 }
 183
 184 gdouble
 185 ves_icall_System_Math_Log (gdouble x)
 186 {
 187         if (x == 0)
 188                 return MInfinity.d;
 189         else if (x < 0)
 190                 return NaN.d;
 191
 192         return log (x);
 193 }
 194
 195 gdouble
 196 ves_icall_System_Math_Log10 (gdouble x)
 197 {
 198         if (x == 0)
 199                 return MInfinity.d;
 200         else if (x < 0)
 201                 return NaN.d;
 202
 203         return log10 (x);
 204 }
 205
 206 gdouble
 207 ves_icall_System_Math_Pow (gdouble x, gdouble y)
 208 {
 209         gdouble result;
 210
 211         if (isnan (y))
 212                 return y;
 213         if (isnan (x))
 214                 return x;
 215
 216         if (isinfinity (y)) {
 217                 if (isplusone (x))
 218                         return x;
 219                 if (isminusone (x))
 220                         return NaN.d;
 221         }
 222
 223         /* following are cases from PAL_pow which abstract the implementation of pow for posix and win32 platforms
 224          * (https://github.com/dotnet/coreclr/blob/master/src/pal/src/cruntime/finite.cpp#L331) */
 225
 226         if (isplusinfinity (y) && !isnan (x)) {
 227                 if (isplusone (x) || isminusone (x))
 228                         result = NaN.d;
 229                 else if (x > MOne.d && x < POne.d)
 230                         result = 0.0;
 231                 else
 232                         result = PInfinity.d;
 233         } else if (isminusinfinity (y) && !isnan (x)) {
 234                 if (isplusone (x) || isminusone (x))
 235                         result = NaN.d;
 236                 if (x > MOne.d && x < POne.d)
 237                         result = PInfinity.d;
 238                 else
 239                         result = 0.0;
 240         } else if (x == 0.0 && y < 0.0) {
 241                 result = PInfinity.d;
 242         } else if (y == 0.0 && isnan (x)) {
 243                 /* Windows returns NaN for pow(NaN, 0), but POSIX specifies
 244                  * a return value of 1 for that case.  We need to return
 245                  * the same result as Windows. */
 246                 result = NaN.d;
 247         } else {
 248                 result = pow (x, y);
 249         }
 250
 251         if (result == PInfinity.d && x < 0.0 && isfinite (x) && ceil (y / 2) != floor (y / 2))
 252                 result = MInfinity.d;
 253
 254         /*
 255          * The even/odd test in the if (this one and the one above) used to be ((long long) y % 2 == 0)
 256          * on SPARC (long long) y for large y (>2**63) is always 0x7fffffff7fffffff, which
 257          * is an odd number, so the test ((long long) y % 2 == 0) will always fail for
 258          * large y. Since large double numbers are always even (e.g., the representation of
 259          * 1E20+1 is the same as that of 1E20, the last .+1. is too insignificant to be part
 260          * of the representation), this test will always return the wrong result for large y.
 261          *
 262          * The (ceil(y/2) == floor(y/2)) test is slower, but more robust.
 263          */
 264         if (result == MInfinity.d && x < 0.0 && isfinite (x) && ceil (y / 2) == floor (y / 2))
 265                 result = PInfinity.d;
 266
 267 #if defined (__linux__) && SIZEOF_VOID_P == 4
 268         /* On Linux 32bits, some tests erroneously return NaN */
 269         if (isnan (result)) {
 270                 if (isminusone (x) && (y > 9007199254740991.0 || y < -9007199254740991.0)) {
 271                         /* Math.Pow (-1, Double.MaxValue) and Math.Pow (-1, Double.MinValue) should return 1 */
 272                         result = POne.d;
 273                 } else if (x < -9007199254740991.0 && y < -9007199254740991.0) {
 274                         /* Math.Pow (Double.MinValue, Double.MinValue) should return 0 */
 275                         result = 0.0;
 276                 } else if (x < -9007199254740991.0 && y > 9007199254740991.0) {
 277                         /* Math.Pow (Double.MinValue, Double.MaxValue) should return Double.PositiveInfinity */
 278                         result = PInfinity.d;
 279                 }
 280         }
 281 #endif
 282
 283         return result == -0.0 ? 0 : result;
 284 }
 285
 286 gdouble
 287 ves_icall_System_Math_Sqrt (gdouble x)
 288 {
 289         if (x < 0)
 290                 return NaN.d;
 291
 292         return sqrt (x);
 293 }
 294
 295 gdouble
 296 ves_icall_System_Math_Abs_double (gdouble v)
 297 {
 298         return fabs (v);
 299 }
 300
 301 gfloat
 302 ves_icall_System_Math_Abs_single (gfloat v)
 303 {
 304         return fabsf (v);
 305 }
 306
 307 gdouble
 308 ves_icall_System_Math_Ceiling (gdouble v)
 309 {
 310         return ceil (v);
 311 }
 312
 313 gdouble
 314 ves_icall_System_Math_SplitFractionDouble (gdouble *v)
 315 {
 316         return modf (*v, v);
 317 }