1 | module math |
2 | |
3 | const ( |
4 | tan_p = [ |
5 | -1.30936939181383777646e+4, |
6 | 1.15351664838587416140e+6, |
7 | -1.79565251976484877988e+7, |
8 | ] |
9 | tan_q = [ |
10 | 1.00000000000000000000e+0, |
11 | 1.36812963470692954678e+4, |
12 | -1.32089234440210967447e+6, |
13 | 2.50083801823357915839e+7, |
14 | -5.38695755929454629881e+7, |
15 | ] |
16 | tan_dp1 = 7.853981554508209228515625e-1 |
17 | tan_dp2 = 7.94662735614792836714e-9 |
18 | tan_dp3 = 3.06161699786838294307e-17 |
19 | tan_lossth = 1.073741824e+9 |
20 | ) |
21 | |
22 | // tan calculates tangent of a number |
23 | pub fn tan(a f64) f64 { |
24 | mut x := a |
25 | if x == 0.0 || is_nan(x) { |
26 | return x |
27 | } |
28 | if is_inf(x, 0) { |
29 | return nan() |
30 | } |
31 | mut sign := 1 // make argument positive but save the sign |
32 | if x < 0 { |
33 | x = -x |
34 | sign = -1 |
35 | } |
36 | if x > math.tan_lossth { |
37 | return 0.0 |
38 | } |
39 | // compute x mod pi_4 |
40 | mut y := floor(x * 4.0 / pi) // strip high bits of integer part |
41 | mut z := ldexp(y, -3) |
42 | z = floor(z) // integer part of y/8 |
43 | z = y - ldexp(z, 3) // y - 16 * (y/16) // integer and fractional part modulo one octant |
44 | mut octant := int(z) // map zeros and singularities to origin |
45 | if (octant & 1) == 1 { |
46 | octant++ |
47 | y += 1.0 |
48 | } |
49 | z = ((x - y * math.tan_dp1) - y * math.tan_dp2) - y * math.tan_dp3 |
50 | zz := z * z |
51 | if zz > 1.0e-14 { |
52 | y = z + z * (zz * (((math.tan_p[0] * zz) + math.tan_p[1]) * zz + math.tan_p[2]) / ((((zz + |
53 | math.tan_q[1]) * zz + math.tan_q[2]) * zz + math.tan_q[3]) * zz + math.tan_q[4])) |
54 | } else { |
55 | y = z |
56 | } |
57 | if (octant & 2) == 2 { |
58 | y = -1.0 / y |
59 | } |
60 | if sign < 0 { |
61 | y = -y |
62 | } |
63 | return y |
64 | } |
65 | |
66 | // tanf calculates tangent. (float32) |
67 | [inline] |
68 | pub fn tanf(a f32) f32 { |
69 | return f32(tan(a)) |
70 | } |
71 | |
72 | // tan calculates cotangent of a number |
73 | pub fn cot(a f64) f64 { |
74 | mut x := a |
75 | if x == 0.0 { |
76 | return inf(1) |
77 | } |
78 | mut sign := 1 // make argument positive but save the sign |
79 | if x < 0 { |
80 | x = -x |
81 | sign = -1 |
82 | } |
83 | if x > math.tan_lossth { |
84 | return 0.0 |
85 | } |
86 | // compute x mod pi_4 |
87 | mut y := floor(x * 4.0 / pi) // strip high bits of integer part |
88 | mut z := ldexp(y, -3) |
89 | z = floor(z) // integer part of y/8 |
90 | z = y - ldexp(z, 3) // y - 16 * (y/16) // integer and fractional part modulo one octant |
91 | mut octant := int(z) // map zeros and singularities to origin |
92 | if (octant & 1) == 1 { |
93 | octant++ |
94 | y += 1.0 |
95 | } |
96 | z = ((x - y * math.tan_dp1) - y * math.tan_dp2) - y * math.tan_dp3 |
97 | zz := z * z |
98 | if zz > 1.0e-14 { |
99 | y = z + z * (zz * (((math.tan_p[0] * zz) + math.tan_p[1]) * zz + math.tan_p[2]) / ((((zz + |
100 | math.tan_q[1]) * zz + math.tan_q[2]) * zz + math.tan_q[3]) * zz + math.tan_q[4])) |
101 | } else { |
102 | y = z |
103 | } |
104 | if (octant & 2) == 2 { |
105 | y = -y |
106 | } else { |
107 | y = 1.0 / y |
108 | } |
109 | if sign < 0 { |
110 | y = -y |
111 | } |
112 | return y |
113 | } |