Print this page
11210 libm should be cstyle(1ONBLD) clean

Split Close
Expand all
Collapse all
          --- old/usr/src/lib/libm/common/R/atanf.c
          +++ new/usr/src/lib/libm/common/R/atanf.c
↓ open down ↓ 14 lines elided ↑ open up ↑
  15   15   * If applicable, add the following below this CDDL HEADER, with the
  16   16   * fields enclosed by brackets "[]" replaced with your own identifying
  17   17   * information: Portions Copyright [yyyy] [name of copyright owner]
  18   18   *
  19   19   * CDDL HEADER END
  20   20   */
  21   21  
  22   22  /*
  23   23   * Copyright 2011 Nexenta Systems, Inc.  All rights reserved.
  24   24   */
       25 +
  25   26  /*
  26   27   * Copyright 2006 Sun Microsystems, Inc.  All rights reserved.
  27   28   * Use is subject to license terms.
  28   29   */
  29   30  
  30   31  #pragma weak __atanf = atanf
  31   32  
  32      -/* INDENT OFF */
       33 +
  33   34  /*
  34   35   * float atanf(float x);
  35   36   * Table look-up algorithm
  36   37   * By K.C. Ng, March 9, 1989
  37   38   *
  38   39   * Algorithm.
  39   40   *
  40   41   * The algorithm is based on atan(x)=atan(y)+atan((x-y)/(1+x*y)).
  41   42   * We use poly1(x) to approximate atan(x) for x in [0,1/8] with
  42   43   * error (relative)
  43      - *      |(atan(x)-poly1(x))/x|<= 2^-115.94      long double
  44      - *      |(atan(x)-poly1(x))/x|<= 2^-58.85       double
  45      - *      |(atan(x)-poly1(x))/x|<= 2^-25.53       float
       44 + *      |(atan(x)-poly1(x))/x|<= 2^-115.94      long double
       45 + *      |(atan(x)-poly1(x))/x|<= 2^-58.85       double
       46 + *      |(atan(x)-poly1(x))/x|<= 2^-25.53       float
  46   47   * and use poly2(x) to approximate atan(x) for x in [0,1/65] with
  47   48   * error (absolute)
  48   49   *      |atan(x)-poly2(x)|<= 2^-122.15          long double
  49   50   *      |atan(x)-poly2(x)|<= 2^-64.79           double
  50   51   *      |atan(x)-poly2(x)|<= 2^-35.36           float
  51   52   * and use poly3(x) to approximate atan(x) for x in [1/8,7/16] with
  52   53   * error (relative, on for single precision)
  53      - *      |(atan(x)-poly1(x))/x|<= 2^-25.53       float
       54 + *      |(atan(x)-poly1(x))/x|<= 2^-25.53       float
  54   55   *
  55   56   * Here poly1-3 are odd polynomial with the following form:
  56   57   *              x + x^3*(a1+x^2*(a2+...))
  57   58   *
  58   59   * (0). Purge off Inf and NaN and 0
  59   60   * (1). Reduce x to positive by atan(x) = -atan(-x).
  60   61   * (2). For x <= 1/8, use
  61   62   *      (2.1) if x < 2^(-prec/2-2), atan(x) =  x  with inexact
  62   63   *      (2.2) Otherwise
  63   64   *              atan(x) = poly1(x)
↓ open down ↓ 15 lines elided ↑ open up ↑
  79   80   *      atan(x) = atan(y) + poly1(s)
  80   81   *              = _TBL_r_atan_hi[j] + (_TBL_r_atan_lo[j] + poly2(s) )
  81   82   *
  82   83   *      Note. |s| <= 1.5384615385e-02 = 1/65. Maxium occurs at x = 1.03125
  83   84   *
  84   85   */
  85   86  
  86   87  #include "libm.h"
  87   88  
  88   89  extern const float _TBL_r_atan_hi[], _TBL_r_atan_lo[];
  89      -static const float
  90      -        big     =   1.0e37F,
  91      -        one     =   1.0F,
  92      -        p1      =  -3.333185951111688247225368498733544672172e-0001F,
  93      -        p2      =   1.969352894213455405211341983203180636021e-0001F,
  94      -        q1      =  -3.332921964095646819563419704110132937456e-0001F,
  95      -        a1      =  -3.333323465223893614063523351509338934592e-0001F,
  96      -        a2      =   1.999425625935277805494082274808174062403e-0001F,
  97      -        a3      =  -1.417547090509737780085769846290301788559e-0001F,
  98      -        a4      =   1.016250813871991983097273733227432685084e-0001F,
  99      -        a5      =  -5.137023693688358515753093811791755221805e-0002F,
 100      -        pio2hi  =   1.570796371e+0000F,
 101      -        pio2lo  =  -4.371139000e-0008F;
 102      -/* INDENT ON */
       90 +static const float big = 1.0e37F,
       91 +        one = 1.0F,
       92 +        p1 = -3.333185951111688247225368498733544672172e-0001F,
       93 +        p2 = 1.969352894213455405211341983203180636021e-0001F,
       94 +        q1 = -3.332921964095646819563419704110132937456e-0001F,
       95 +        a1 = -3.333323465223893614063523351509338934592e-0001F,
       96 +        a2 = 1.999425625935277805494082274808174062403e-0001F,
       97 +        a3 = -1.417547090509737780085769846290301788559e-0001F,
       98 +        a4 = 1.016250813871991983097273733227432685084e-0001F,
       99 +        a5 = -5.137023693688358515753093811791755221805e-0002F,
      100 +        pio2hi = 1.570796371e+0000F,
      101 +        pio2lo = -4.371139000e-0008F;
      102 +
 103  103  
 104  104  float
 105      -atanf(float xx) {
      105 +atanf(float xx)
      106 +{
 106  107          float x, y, z, r, p, s;
 107  108          volatile double dummy __unused;
 108  109          int ix, iy, sign, j;
 109  110  
 110  111          x = xx;
 111      -        ix = *(int *) &x;
      112 +        ix = *(int *)&x;
 112  113          sign = ix & 0x80000000;
 113  114          ix ^= sign;
 114  115  
 115  116          /* for |x| < 1/8 */
 116  117          if (ix < 0x3e000000) {
 117      -                if (ix < 0x38800000) {  /* if |x| < 2**(-prec/2-2) */
      118 +                if (ix < 0x38800000) {          /* if |x| < 2**(-prec/2-2) */
 118  119                          dummy = big + x;        /* get inexact flag if x != 0 */
 119  120  #ifdef lint
 120  121                          dummy = dummy;
 121  122  #endif
 122  123                          return (x);
 123  124                  }
      125 +
 124  126                  z = x * x;
      127 +
 125  128                  if (ix < 0x3c000000) {  /* if |x| < 2**(-prec/4-1) */
 126  129                          x = x + (x * z) * p1;
 127  130                          return (x);
 128  131                  } else {
 129  132                          x = x + (x * z) * (p1 + z * p2);
 130  133                          return (x);
 131  134                  }
 132  135          }
 133  136  
 134  137          /* for |x| >= 8.0 */
 135  138          if (ix >= 0x41000000) {
 136      -                *(int *) &x = ix;
 137      -                if (ix < 0x42820000) {  /* x <  65 */
      139 +                *(int *)&x = ix;
      140 +
      141 +                if (ix < 0x42820000) {                          /* x <  65 */
 138  142                          r = one / x;
 139  143                          z = r * r;
 140  144                          y = r * (one + z * (p1 + z * p2));      /* poly1 */
 141  145                          y -= pio2lo;
 142      -                } else if (ix < 0x44800000) {   /* x <  2**(prec/3+2) */
      146 +                } else if (ix < 0x44800000) { /* x <  2**(prec/3+2) */
 143  147                          r = one / x;
 144  148                          z = r * r;
 145  149                          y = r * (one + z * q1); /* poly2 */
 146  150                          y -= pio2lo;
 147      -                } else if (ix < 0x4c800000) {   /* x <  2**(prec+2) */
      151 +                } else if (ix < 0x4c800000) { /* x <  2**(prec+2) */
 148  152                          y = one / x - pio2lo;
 149      -                } else if (ix < 0x7f800000) {   /* x <  inf */
      153 +                } else if (ix < 0x7f800000) { /* x <  inf */
 150  154                          y = -pio2lo;
 151      -                } else {                /* x is inf or NaN */
 152      -                        if (ix > 0x7f800000) {
      155 +                } else {        /* x is inf or NaN */
      156 +                        if (ix > 0x7f800000)
 153  157                                  return (x * x); /* - -> * for Cheetah */
 154      -                        }
      158 +
 155  159                          y = -pio2lo;
 156  160                  }
 157  161  
 158  162                  if (sign == 0)
 159  163                          x = pio2hi - y;
 160  164                  else
 161  165                          x = y - pio2hi;
      166 +
 162  167                  return (x);
 163  168          }
 164  169  
 165      -
 166  170          /* now x is between 1/8 and 8 */
 167      -        if (ix < 0x3f000000) {  /* between 1/8 and 1/2 */
      171 +        if (ix < 0x3f000000) {          /* between 1/8 and 1/2 */
 168  172                  z = x * x;
 169      -                x = x + (x * z) * (a1 + z * (a2 + z * (a3 + z * (a4 +
 170      -                        z * a5))));
      173 +                x = x + (x * z) * (a1 + z * (a2 + z * (a3 + z * (a4 + z *
      174 +                    a5))));
 171  175                  return (x);
 172  176          }
 173      -        *(int *) &x = ix;
      177 +
      178 +        *(int *)&x = ix;
 174  179          iy = (ix + 0x00040000) & 0x7ff80000;
 175      -        *(int *) &y = iy;
      180 +        *(int *)&y = iy;
 176  181          j = (iy - 0x3f000000) >> 19;
 177  182  
 178      -        if (ix == iy)
 179      -                p = x - y;      /* p=0.0 */
 180      -        else {
      183 +        if (ix == iy) {
      184 +                p = x - y;              /* p=0.0 */
      185 +        } else {
 181  186                  if (sign == 0)
 182  187                          s = (x - y) / (one + x * y);
 183  188                  else
 184  189                          s = (y - x) / (one + x * y);
      190 +
 185  191                  z = s * s;
 186  192                  p = s * (one + z * q1);
 187  193          }
      194 +
 188  195          if (sign == 0) {
 189  196                  r = p + _TBL_r_atan_lo[j];
 190  197                  x = r + _TBL_r_atan_hi[j];
 191  198          } else {
 192  199                  r = p - _TBL_r_atan_lo[j];
 193  200                  x = r - _TBL_r_atan_hi[j];
 194  201          }
      202 +
 195  203          return (x);
 196  204  }
    
XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX