Subversion Repositories Kolibri OS

Rev

Go to most recent revision | Blame | Last modification | View Log | Download | RSS feed

  1. /* 2009 for Newlib:  Sun's s_ilogb.c converted to be s_logb.c.  */
  2. /* @(#)s_ilogb.c 5.1 93/09/24 */
  3. /*
  4.  * ====================================================
  5.  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
  6.  *
  7.  * Developed at SunPro, a Sun Microsystems, Inc. business.
  8.  * Permission to use, copy, modify, and distribute this
  9.  * software is freely granted, provided that this notice
  10.  * is preserved.
  11.  * ====================================================
  12.  */
  13. /*
  14. FUNCTION
  15.        <<logb>>, <<logbf>>--get exponent of floating-point number
  16. INDEX
  17.         logb
  18. INDEX
  19.         logbf
  20.  
  21. ANSI_SYNOPSIS
  22.         #include <math.h>
  23.         double logb(double <[x]>);
  24.         float logbf(float <[x]>);
  25.  
  26. DESCRIPTION
  27. The <<logb>> functions extract the exponent of <[x]>, as a signed integer value
  28. in floating-point format.  If <[x]> is subnormal it is treated as though it were
  29. normalized; thus, for positive finite <[x]>,
  30. @ifnottex
  31. 1 <= (<[x]> * FLT_RADIX to the power (-logb(<[x]>))) < FLT_RADIX.
  32. @end ifnottex
  33. @tex
  34. $1 \leq ( x \cdot FLT\_RADIX ^ {-logb(x)} ) < FLT\_RADIX$.
  35. @end tex
  36. A domain error may occur if the argument is zero.
  37. In this floating-point implementation, FLT_RADIX is 2.  Which also means
  38. that for finite <[x]>, <<logb>>(<[x]>) = <<floor>>(<<log2>>(<<fabs>>(<[x]>))).
  39.  
  40. All nonzero, normal numbers can be described as
  41. @ifnottex
  42. <[m]> * 2**<[p]>, where 1.0 <= <[m]> < 2.0.
  43. @end ifnottex
  44. @tex
  45. $m \cdot 2^p$, where $1.0 \leq m < 2.0$.
  46. @end tex
  47. The <<logb>> functions examine the argument <[x]>, and return <[p]>.
  48. The <<frexp>> functions are similar to the <<logb>> functions, but
  49. returning <[m]> adjusted to the interval [.5, 1) or 0, and <[p]>+1.
  50.  
  51. RETURNS
  52. @comment Formatting note:  "$@" forces a new line
  53. When <[x]> is:@*
  54. +inf or -inf, +inf is returned;@*
  55. NaN, NaN is returned;@*
  56. 0, -inf is returned, and the divide-by-zero exception is raised;@*
  57. otherwise, the <<logb>> functions return the signed exponent of <[x]>.
  58.  
  59. PORTABILITY
  60. ANSI C, POSIX
  61.  
  62. SEEALSO
  63. frexp, ilogb
  64. */
  65.  
  66. /* double logb(double x)
  67.  * return the binary exponent of non-zero x
  68.  * logb(0) = -inf, raise divide-by-zero floating point exception
  69.  * logb(+inf|-inf) = +inf (no signal is raised)
  70.  * logb(NaN) = NaN (no signal is raised)
  71.  * Per C99 recommendation, a NaN argument is returned unchanged.
  72.  */
  73.  
  74. #include "fdlibm.h"
  75.  
  76. #ifndef _DOUBLE_IS_32BITS
  77.  
  78. double
  79. #ifdef __STDC__
  80. logb(double x)
  81. #else
  82. logb(x)
  83. double x;
  84. #endif
  85. {
  86.         __int32_t hx,lx,ix;
  87.  
  88.         EXTRACT_WORDS(hx,lx,x);
  89.         hx &= 0x7fffffff;               /* high |x| */
  90.         if(hx<0x00100000) {             /* 0 or subnormal */
  91.             if((hx|lx)==0)  {
  92.                 double  xx;
  93.                 /* arg==0:  return -inf and raise divide-by-zero exception */
  94.                 INSERT_WORDS(xx,hx,lx); /* +0.0 */
  95.                 return -1./xx;  /* logb(0) = -inf */
  96.                 }
  97.             else                        /* subnormal x */
  98.                 if(hx==0) {
  99.                     for (ix = -1043; lx>0; lx<<=1) ix -=1;
  100.                 } else {
  101.                     for (ix = -1022,hx<<=11; hx>0; hx<<=1) ix -=1;
  102.                 }
  103.             return (double) ix;
  104.         }
  105.         else if (hx<0x7ff00000) return (hx>>20)-1023;   /* normal # */
  106.         else if (hx>0x7ff00000 || lx)  return x;        /* x==NaN */
  107.         else  return HUGE_VAL;  /* x==inf (+ or -) */
  108. }
  109.  
  110. #endif /* _DOUBLE_IS_32BITS */
  111.