[469] | 1 | /* |
---|

| 2 | * ==================================================== |
---|

| 3 | * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. |
---|

| 4 | * |
---|

| 5 | * Developed at SunPro, a Sun Microsystems, Inc. business. |
---|

| 6 | * Permission to use, copy, modify, and distribute this |
---|

| 7 | * software is freely granted, provided that this notice |
---|

| 8 | * is preserved. |
---|

| 9 | * ==================================================== |
---|

| 10 | */ |
---|

| 11 | |
---|

| 12 | /* |
---|

| 13 | * Modified for ALMOS-MKH OS at UPMC, France, August 2018. (Alain Greiner) |
---|

| 14 | */ |
---|

| 15 | |
---|

| 16 | /* |
---|

| 17 | * scalbln(double x, long n) |
---|

| 18 | * scalbln(x,n) returns x * 2**n computed by exponent |
---|

| 19 | * manipulation rather than by actually performing an |
---|

| 20 | * exponentiation or a multiplication. |
---|

| 21 | */ |
---|

| 22 | |
---|

| 23 | #include "math.h" |
---|

| 24 | #include "math_private.h" |
---|

| 25 | #include <limits.h> |
---|

| 26 | |
---|

| 27 | static const double |
---|

| 28 | two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */ |
---|

| 29 | twom54 = 5.55111512312578270212e-17, /* 0x3C900000, 0x00000000 */ |
---|

| 30 | huge = 1.0e+300, |
---|

| 31 | tiny = 1.0e-300; |
---|

| 32 | |
---|

| 33 | double scalbln(double x, long n) |
---|

| 34 | { |
---|

| 35 | int32_t k, hx, lx; |
---|

| 36 | |
---|

| 37 | EXTRACT_WORDS(hx, lx, x); |
---|

| 38 | k = (hx & 0x7ff00000) >> 20; /* extract exponent */ |
---|

| 39 | if (k == 0) { /* 0 or subnormal x */ |
---|

| 40 | if ((lx | (hx & 0x7fffffff)) == 0) |
---|

| 41 | return x; /* +-0 */ |
---|

| 42 | x *= two54; |
---|

| 43 | GET_HIGH_WORD(hx, x); |
---|

| 44 | k = ((hx & 0x7ff00000) >> 20) - 54; |
---|

| 45 | } |
---|

| 46 | if (k == 0x7ff) |
---|

| 47 | return x + x; /* NaN or Inf */ |
---|

| 48 | k = k + n; |
---|

| 49 | if (k > 0x7fe) |
---|

| 50 | return huge * copysign(huge, x); /* overflow */ |
---|

| 51 | if (n < -50000) |
---|

| 52 | return tiny * copysign(tiny, x); /* underflow */ |
---|

| 53 | if (k > 0) { /* normal result */ |
---|

| 54 | SET_HIGH_WORD(x, (hx & 0x800fffff) | (k << 20)); |
---|

| 55 | return x; |
---|

| 56 | } |
---|

| 57 | if (k <= -54) { |
---|

| 58 | if (n > 50000) /* in case integer overflow in n+k */ |
---|

| 59 | return huge * copysign(huge, x); /* overflow */ |
---|

| 60 | return tiny * copysign(tiny, x); /* underflow */ |
---|

| 61 | } |
---|

| 62 | k += 54; /* subnormal result */ |
---|

| 63 | SET_HIGH_WORD(x, (hx & 0x800fffff) | (k << 20)); |
---|

| 64 | return x * twom54; |
---|

| 65 | } |
---|

| 66 | |
---|

| 67 | double scalbn(double x, int n) |
---|

| 68 | { |
---|

| 69 | return scalbln(x, n); |
---|

| 70 | } |
---|

| 71 | |
---|