lib/math/asin.c

   1 /*
   2  * (c) copyright 1988 by the Vrije Universiteit, Amsterdam, The Netherlands.
   3  * See the copyright notice in the ACK home directory, in the file "Copyright".
   4  *
   5  * Author: Ceriel J.H. Jacobs
   6  */
   7 /* $Header$ */
   8
   9 #include        <math.h>
  10 #include        <errno.h>
  11 #include        "localmath.h"
  12
  13 static double
  14 asin_acos(double x, int cosfl)
  15 {
  16         int negative = x < 0;
  17         int     i;
  18         double  g;
  19         static double p[] = {
  20                 -0.27368494524164255994e+2,
  21                  0.57208227877891731407e+2,
  22                 -0.39688862997540877339e+2,
  23                  0.10152522233806463645e+2,
  24                 -0.69674573447350646411e+0
  25         };
  26         static double q[] = {
  27                 -0.16421096714498560795e+3,
  28                  0.41714430248260412556e+3,
  29                 -0.38186303361750149284e+3,
  30                  0.15095270841030604719e+3,
  31                 -0.23823859153670238830e+2,
  32                  1.0
  33         };
  34
  35         if (__IsNan(x)) {
  36                 errno = EDOM;
  37                 return x;
  38         }
  39
  40         if (negative) {
  41                 x = -x;
  42         }
  43         if (x > 0.5) {
  44                 i = 1;
  45                 if (x > 1) {
  46                         errno = EDOM;
  47                         return 0;
  48                 }
  49                 g = 0.5 - 0.5 * x;
  50                 x = - sqrt(g);
  51                 x += x;
  52         }
  53         else {
  54                 /* ??? avoid underflow ??? */
  55                 i = 0;
  56                 g = x * x;
  57         }
  58         x += x * g * POLYNOM4(g, p) / POLYNOM5(g, q);
  59         if (cosfl) {
  60                 if (! negative) x = -x;
  61         }
  62         if ((cosfl == 0) == (i == 1)) {
  63                 x = (x + M_PI_4) + M_PI_4;
  64         }
  65         else if (cosfl && negative && i == 1) {
  66                 x = (x + M_PI_2) + M_PI_2;
  67         }
  68         if (! cosfl && negative) x = -x;
  69         return x;
  70 }
  71
  72 double
  73 asin(double x)
  74 {
  75         return asin_acos(x, 0);
  76 }
  77
  78 double
  79 acos(double x)
  80 {
  81         return asin_acos(x, 1);
  82 }