xref: /plan9/sys/src/ape/lib/ap/math/tan.c (revision 3e12c5d1bb89fc02707907988834ef147769ddaf)
1 /*
2 	floating point tangent
3 
4 	A series is used after range reduction.
5 	Coefficients are #4285 from Hart & Cheney. (19.74D)
6  */
7 
8 #include <math.h>
9 #include <errno.h>
10 
11 static double invpi	  = 1.27323954473516268;
12 static double p0	 = -0.1306820264754825668269611177e+5;
13 static double p1	  = 0.1055970901714953193602353981e+4;
14 static double p2	 = -0.1550685653483266376941705728e+2;
15 static double p3	  = 0.3422554387241003435328470489e-1;
16 static double p4	  = 0.3386638642677172096076369e-4;
17 static double q0	 = -0.1663895238947119001851464661e+5;
18 static double q1	  = 0.4765751362916483698926655581e+4;
19 static double q2	 = -0.1555033164031709966900124574e+3;
20 
21 double
tan(double arg)22 tan(double arg)
23 {
24 	double sign, temp, e, x, xsq;
25 	int flag, i;
26 
27 	flag = 0;
28 	sign = 1;
29 	if(arg < 0){
30 		arg = -arg;
31 		sign = -1;
32 	}
33 	arg = arg*invpi;   /* overflow? */
34 	x = modf(arg, &e);
35 	i = e;
36 	switch(i%4) {
37 	case 1:
38 		x = 1 - x;
39 		flag = 1;
40 		break;
41 
42 	case 2:
43 		sign = - sign;
44 		flag = 1;
45 		break;
46 
47 	case 3:
48 		x = 1 - x;
49 		sign = - sign;
50 		break;
51 
52 	case 0:
53 		break;
54 	}
55 
56 	xsq = x*x;
57 	temp = ((((p4*xsq+p3)*xsq+p2)*xsq+p1)*xsq+p0)*x;
58 	temp = temp/(((xsq+q2)*xsq+q1)*xsq+q0);
59 
60 	if(flag == 1) {
61 		if(temp == 0) {
62 			errno = EDOM;
63 			if (sign > 0)
64 				return HUGE_VAL;
65 			return -HUGE_VAL;
66 		}
67 		temp = 1/temp;
68 	}
69 	return sign*temp;
70 }
71