s_scalbn.c (2374B)
1/* @(#)s_scalbn.c 5.1 93/09/24 */ 2/* 3 * ==================================================== 4 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. 5 * 6 * Developed at SunPro, a Sun Microsystems, Inc. business. 7 * Permission to use, copy, modify, and distribute this 8 * software is freely granted, provided that this notice 9 * is preserved. 10 * ==================================================== 11 */ 12 13#if defined(LIBM_SCCS) && !defined(lint) 14static const char rcsid[] = 15 "$NetBSD: s_scalbn.c,v 1.8 1995/05/10 20:48:08 jtc Exp $"; 16#endif 17 18/* 19 * scalbn (double x, int n) 20 * scalbn(x,n) returns x* 2**n computed by exponent 21 * manipulation rather than by actually performing an 22 * exponentiation or a multiplication. 23 */ 24 25#include "math_libm.h" 26#include "math_private.h" 27 28libm_hidden_proto(copysign) 29#ifdef __STDC__ 30 static const double 31#else 32 static double 33#endif 34 two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */ 35 twom54 = 5.55111512312578270212e-17, /* 0x3C900000, 0x00000000 */ 36 huge_val = 1.0e+300, tiny = 1.0e-300; 37 38libm_hidden_proto(scalbn) 39#ifdef __STDC__ 40 double scalbn(double x, int n) 41#else 42 double scalbn(x, n) 43 double x; 44 int n; 45#endif 46{ 47 int32_t k, hx, lx; 48 EXTRACT_WORDS(hx, lx, x); 49 k = (hx & 0x7ff00000) >> 20; /* extract exponent */ 50 if (k == 0) { /* 0 or subnormal x */ 51 if ((lx | (hx & 0x7fffffff)) == 0) 52 return x; /* +-0 */ 53 x *= two54; 54 GET_HIGH_WORD(hx, x); 55 k = ((hx & 0x7ff00000) >> 20) - 54; 56 if (n < -50000) 57 return tiny * x; /* underflow */ 58 } 59 if (k == 0x7ff) 60 return x + x; /* NaN or Inf */ 61 k = k + n; 62 if (k > 0x7fe) 63 return huge_val * copysign(huge_val, x); /* overflow */ 64 if (k > 0) { /* normal result */ 65 SET_HIGH_WORD(x, (hx & 0x800fffff) | (k << 20)); 66 return x; 67 } 68 if (k <= -54) { 69 if (n > 50000) /* in case integer overflow in n+k */ 70 return huge_val * copysign(huge_val, x); /* overflow */ 71 else 72 return tiny * copysign(tiny, x); /* underflow */ 73 } 74 k += 54; /* subnormal result */ 75 SET_HIGH_WORD(x, (hx & 0x800fffff) | (k << 20)); 76 return x * twom54; 77} 78 79libm_hidden_def(scalbn)