src/libm/s_scalbn.c
author Sam Lantinga <slouken@libsdl.org>
Sat, 04 Nov 2017 17:35:03 -0700
changeset 11684 eccdf37f8996
parent 11683 48bcba563d9c
child 12420 4a6c91d9cc33
permissions -rw-r--r--
Added SDL_fmod() and SDL_fmodf()
     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 
    65 double scalbn(double x, int n)
    66 {
    67 	return scalbln(x, n);
    68 }
    69 libm_hidden_def(scalbn)