3 * System.Double, System.Single and System.Number runtime support
6 * Ludovic Henry (ludovic@xamarin.com)
8 * Copyright 2015 Xamarin, Inc (http://www.xamarin.com)
9 * Licensed under the MIT license. See LICENSE file in the project root for full license information.
13 // Copyright (c) Microsoft. All rights reserved.
14 // Licensed under the MIT license. See LICENSE file in the project root for full license information.
17 // - src/classlibnative/bcltype/number.cpp
19 // Ported from C++ to C and adjusted to Mono runtime
23 #include "number-ms.h"
25 static const guint64 rgval64Power10[] = {
61 static const gint8 rgexp64Power10[] = {
62 /* exponents for both powers of 10 and 0.1 */
80 static const guint64 rgval64Power10By16[] = {
100 0xe950df20247c83f8LL,
101 0x81842f29f2cce373LL,
102 0x8fcac257558ee4e2LL,
104 /* powers of 0.1^16 */
105 0xe69594bec44de160LL,
106 0xcfb11ead453994c3LL,
107 0xbb127c53b17ec165LL,
108 0xa87fea27a539e9b3LL,
109 0x97c560ba6b0919b5LL,
110 0x88b402f7fd7553abLL,
111 0xf64335bcf065d3a0LL,
112 0xddd0467c64bce4c4LL,
113 0xc7caba6e7c5382edLL,
114 0xb3f4e093db73a0b7LL,
115 0xa21727db38cb0053LL,
116 0x91ff83775423cc29LL,
117 0x8380dea93da4bc82LL,
118 0xece53cec4a314f00LL,
119 0xd5605fcdcf32e217LL,
120 0xc0314325637a1978LL,
121 0xad1c8eab5ee43ba2LL,
122 0x9becce62836ac5b0LL,
123 0x8c71dcd9ba0b495cLL,
124 0xfd00b89747823938LL,
125 0xe3e27a444d8d991aLL,
128 static const gint16 rgexp64Power10By16[] = {
129 /* exponents for both powers of 10^16 and 0.1^16 */
153 static inline guint64
154 digits_to_int (guint16 *p, int count)
156 g_assert (1 <= count && count <= 9);
160 case 9: res += 100000000 * (p [i++] - '0');
161 case 8: res += 10000000 * (p [i++] - '0');
162 case 7: res += 1000000 * (p [i++] - '0');
163 case 6: res += 100000 * (p [i++] - '0');
164 case 5: res += 10000 * (p [i++] - '0');
165 case 4: res += 1000 * (p [i++] - '0');
166 case 3: res += 100 * (p [i++] - '0');
167 case 2: res += 10 * (p [i++] - '0');
168 case 1: res += 1 * (p [i++] - '0');
173 static inline guint64
174 mul_64_lossy (guint64 a, guint64 b, gint *pexp)
176 /* it's ok to losse some precision here - it will be called
177 * at most twice during the conversion, so the error won't
178 * propagate to any of the 53 significant bits of the result */
180 ((((guint64) (guint32) (a >> 32)) * ((guint64) (guint32) (b >> 32))) )
181 + ((((guint64) (guint32) (a >> 32)) * ((guint64) (guint32) (b ))) >> 32)
182 + ((((guint64) (guint32) (a )) * ((guint64) (guint32) (b >> 32))) >> 32);
185 if ((val & 0x8000000000000000LL) == 0) {
194 number_to_double (MonoNumber *number, gdouble *value)
198 gint exp, remaining, total, count, scale, absscale, index;
201 src = number->digits;
202 while (*src++) total ++;
206 src = number->digits;
207 while (*src == '0') {
212 if (remaining == 0) {
217 count = MIN (remaining, 9);
219 val = digits_to_int (src, count);
222 count = MIN (remaining, 9);
225 /* get the denormalized power of 10 */
226 guint32 mult = (guint32) (rgval64Power10 [count - 1] >> (64 - rgexp64Power10 [count - 1]));
227 val = ((guint64) (guint32) val) * ((guint64) mult) + digits_to_int (src + 9, count);
230 scale = number->scale - (total - remaining);
231 absscale = abs (scale);
233 if (absscale >= 22 * 16) {
234 /* overflow / underflow */
235 *(guint64*) value = (scale > 0) ? 0x7FF0000000000000LL : 0;
241 /* normalize the mantiss */
242 if ((val & 0xFFFFFFFF00000000LL) == 0) { val <<= 32; exp -= 32; }
243 if ((val & 0xFFFF000000000000LL) == 0) { val <<= 16; exp -= 16; }
244 if ((val & 0xFF00000000000000LL) == 0) { val <<= 8; exp -= 8; }
245 if ((val & 0xF000000000000000LL) == 0) { val <<= 4; exp -= 4; }
246 if ((val & 0xC000000000000000LL) == 0) { val <<= 2; exp -= 2; }
247 if ((val & 0x8000000000000000LL) == 0) { val <<= 1; exp -= 1; }
249 index = absscale & 15;
251 gint multexp = rgexp64Power10 [index - 1];
252 /* the exponents are shared between the inverted and regular table */
253 exp += (scale < 0) ? (-multexp + 1) : multexp;
255 guint64 multval = rgval64Power10 [index + ((scale < 0) ? 15 : 0) - 1];
256 val = mul_64_lossy (val, multval, &exp);
259 index = absscale >> 4;
261 gint multexp = rgexp64Power10By16 [index - 1];
262 /* the exponents are shared between the inverted and regular table */
263 exp += (scale < 0) ? (-multexp + 1) : multexp;
265 guint64 multval = rgval64Power10By16 [index + ((scale < 0) ? 21 : 0) - 1];
266 val = mul_64_lossy (val, multval, &exp);
269 if ((guint32) val & (1 << 10)) {
270 /* IEEE round to even */
271 guint64 tmp = val + ((1 << 10) - 1) + (((guint32) val >> 11) & 1);
274 tmp = (tmp >> 1) | 0x8000000000000000LL;
280 /* return the exponent to a biased state */
283 /* handle overflow, underflow, "Epsilon - 1/2 Epsilon", denormalized, and the normal case */
285 if (exp == -52 && (val >= 0x8000000000000058LL)) {
286 /* round X where {Epsilon > X >= 2.470328229206232730000000E-324} up to Epsilon (instead of down to zero) */
287 val = 0x0000000000000001LL;
288 } else if (exp <= -52) {
293 val >>= (-exp + 11 + 1);
295 } else if (exp >= 0x7FF) {
297 val = 0x7FF0000000000000LL;
299 /* normal postive exponent case */
300 val = ((guint64) exp << 52) + ((val >> 11) & 0x000FFFFFFFFFFFFFLL);
303 *(guint64*) value = val;
307 *(guint64*) value |= 0x8000000000000000LL;
311 mono_double_from_number (gpointer from, MonoDouble *target)
313 MonoDouble_double res;
314 guint e, mant_lo, mant_hi;
318 number_to_double ((MonoNumber*) from, &res.d);
320 mant_lo = res.s.mantLo;
321 mant_hi = res.s.mantHi;
326 if (e == 0 && mant_lo == 0 && mant_hi == 0)