newlib/libm/mathfp/s_tan.c

   1
   2 /* @(#)z_tan.c 1.0 98/08/13 */
   3 /******************************************************************
   4  * The following routines are coded directly from the algorithms
   5  * and coefficients given in "Software Manual for the Elementary
   6  * Functions" by William J. Cody, Jr. and William Waite, Prentice
   7  * Hall, 1980.
   8  ******************************************************************/
   9
  10 /*
  11 FUNCTION
  12         <<tan>>, <<tanf>>---tangent
  13
  14 INDEX
  15 tan
  16 INDEX
  17 tanf
  18
  19 SYNOPSIS
  20         #include <math.h>
  21         double tan(double <[x]>);
  22         float tanf(float <[x]>);
  23
  24 DESCRIPTION
  25 <<tan>> computes the tangent of the argument <[x]>.
  26 Angles are specified in radians.
  27
  28 <<tanf>> is identical, save that it takes and returns <<float>> values.
  29
  30 RETURNS
  31 The tangent of <[x]> is returned.
  32
  33 PORTABILITY
  34 <<tan>> is ANSI. <<tanf>> is an extension.
  35 */
  36
  37 /******************************************************************
  38  * Tangent
  39  *
  40  * Input:
  41  *   x - floating point value
  42  *
  43  * Output:
  44  *   tangent of x
  45  *
  46  * Description:
  47  *   This routine calculates the tangent of x.
  48  *
  49  *****************************************************************/
  50
  51 #include "fdlibm.h"
  52 #include "zmath.h"
  53
  54 #ifndef _DOUBLE_IS_32BITS
  55
  56 static const double TWO_OVER_PI = 0.63661977236758134308;
  57 static const double p[] = { -0.13338350006421960681,
  58                              0.34248878235890589960e-2,
  59                             -0.17861707342254426711e-4 };
  60 static const double q[] = { -0.46671683339755294240,
  61                              0.25663832289440112864e-1,
  62                             -0.31181531907010027307e-3,
  63                              0.49819433993786512270e-6 };
  64
  65 double
  66 tan (double x)
  67 {
  68   double y, f, g, XN, xnum, xden, res;
  69   int N;
  70
  71   /* Check for special values. */
  72   switch (numtest (x))
  73     {
  74       case NAN:
  75         errno = EDOM;
  76         return (x);
  77       case INF:
  78         errno = EDOM;
  79         return (z_notanum.d);
  80     }
  81
  82   y = fabs (x);
  83
  84   /* Check for values that are out of our range. */
  85   if (y > 105414357.0)
  86     {
  87       errno = ERANGE;
  88       return (y);
  89     }
  90
  91   if (x < 0.0)
  92     N = (int) (x * TWO_OVER_PI - 0.5);
  93   else
  94     N = (int) (x * TWO_OVER_PI + 0.5);
  95
  96   XN = (double) N;
  97
  98   f = x - N * __PI_OVER_TWO;
  99
 100   /* Check for values that are too small. */
 101   if (-z_rooteps < f && f < z_rooteps)
 102     {
 103       xnum = f;
 104       xden = 1.0;
 105     }
 106
 107   /* Calculate the polynomial. */
 108   else
 109     {
 110       g = f * f;
 111
 112       xnum = f * ((p[2] * g + p[1]) * g + p[0]) * g + f;
 113       xden = (((q[3] * g + q[2]) * g + q[1]) * g + q[0]) * g + 1.0;
 114     }
 115
 116   if (N & 1)
 117     {
 118       xnum = -xnum;
 119       res = xden / xnum;
 120     }
 121   else
 122     {
 123       res = xnum / xden;
 124     }
 125
 126   return (res);
 127 }
 128
 129 #endif /* _DOUBLE_IS_32BITS */