src/libm/s_scalbn.c
author Sam Lantinga <slouken@libsdl.org>
Sat, 04 Nov 2017 15:53:19 -0700
changeset 11683 48bcba563d9c
parent 7678 286c42d7c5ed
child 11684 eccdf37f8996
permissions -rw-r--r--
Updated math code from the uClibc 0.9.33 release
     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  * scalbln(double x, long n)
    14  * scalbln(x,n) returns x * 2**n computed by exponent
    15  * manipulation rather than by actually performing an
    16  * exponentiation or a multiplication.
    17  */
    18 
    19 #include "math_libm.h"
    20 #include "math_private.h"
    21 #include <limits.h>
    22 
    23 static const double
    24 two54  = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
    25 twom54 = 5.55111512312578270212e-17, /* 0x3C900000, 0x00000000 */
    26 huge   = 1.0e+300,
    27 tiny   = 1.0e-300;
    28 
    29 double scalbln(double x, long n)
    30 {
    31 	int32_t k, hx, lx;
    32 
    33 	EXTRACT_WORDS(hx, lx, x);
    34 	k = (hx & 0x7ff00000) >> 20; /* extract exponent */
    35 	if (k == 0) { /* 0 or subnormal x */
    36 		if ((lx | (hx & 0x7fffffff)) == 0)
    37 			return x; /* +-0 */
    38 		x *= two54;
    39 		GET_HIGH_WORD(hx, x);
    40 		k = ((hx & 0x7ff00000) >> 20) - 54;
    41 	}
    42 	if (k == 0x7ff)
    43 		return x + x; /* NaN or Inf */
    44 	k = k + n;
    45 	if (k > 0x7fe)
    46 		return huge * copysign(huge, x); /* overflow */
    47 	if (n < -50000)
    48 		return tiny * copysign(tiny, x); /* underflow */
    49 	if (k > 0) { /* normal result */
    50 		SET_HIGH_WORD(x, (hx & 0x800fffff) | (k << 20));
    51 		return x;
    52 	}
    53 	if (k <= -54) {
    54 		if (n > 50000) /* in case integer overflow in n+k */
    55 			return huge * copysign(huge, x); /* overflow */
    56 		return tiny * copysign(tiny, x); /* underflow */
    57 	}
    58 	k += 54; /* subnormal result */
    59 	SET_HIGH_WORD(x, (hx & 0x800fffff) | (k << 20));
    60 	return x * twom54;
    61 }
    62 libm_hidden_def(scalbln)
    63 
    64 #if LONG_MAX == INT_MAX
    65 /* strong_alias(scalbln, scalbn) - "error: conflicting types for 'scalbn'"
    66  * because it tries to declare "typeof(scalbln) scalbn;"
    67  * which tries to give "long" parameter to scalbn.
    68  * Doing it by hand:
    69  */
    70 __typeof(scalbn) scalbn __attribute__((alias("scalbln")));
    71 #else
    72 double scalbn(double x, int n)
    73 {
    74 	return scalbln(x, n);
    75 }
    76 #endif
    77 libm_hidden_def(scalbn)