2 * IBM Accurate Mathematical Library
3 * written by International Business Machines Corp.
4 * Copyright (C) 2001-2017 Free Software Foundation, Inc.
6 * This program is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU Lesser General Public License as published by
8 * the Free Software Foundation; either version 2.1 of the License, or
9 * (at your option) any later version.
11 * This program is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU Lesser General Public License for more details.
16 * You should have received a copy of the GNU Lesser General Public License
17 * along with this program; if not, see <http://www.gnu.org/licenses/>.
19 /*********************************************************************/
20 /* MODULE_NAME: utan.c */
25 /* FILES NEEDED:dla.h endian.h mpa.h mydefs.h utan.h */
26 /* branred.c sincos32.c mptan.c */
29 /* An ultimate tan routine. Given an IEEE double machine number x */
30 /* it computes the correctly rounded (to nearest) value of tan(x). */
31 /* Assumption: Machine arithmetic operations are performed in */
32 /* round to nearest mode of IEEE 754 standard. */
34 /*********************************************************************/
43 #include <math_private.h>
44 #include <libm-alias-double.h>
46 #include <stap-probe.h>
52 static double tanMp (double);
53 void __mptan (double, mp_no
*, int);
63 double a
, da
, a2
, b
, db
, c
, dc
, c1
, cc1
, c2
, cc2
, c3
, cc3
, fi
, ffi
, gi
, pz
,
64 s
, sy
, t
, t1
, t2
, t3
, t4
, t7
, t8
, t9
, t10
, w
, x2
, xn
, xx2
, y
, ya
,
65 yya
, z0
, z
, zz
, z2
, zz2
;
71 mp_no mpa
, mpt1
, mpt2
;
75 int __branred (double, double *, double *);
76 int __mpranred (double, mp_no
*, int);
78 SET_RESTORE_ROUND_53BIT (FE_TONEAREST
);
82 ux
= num
.i
[HIGH_HALF
];
83 if ((ux
& 0x7ff00000) == 0x7ff00000)
85 if ((ux
& 0x7fffffff) == 0x7ff00000)
91 w
= (x
< 0.0) ? -x
: x
;
93 /* (I) The case abs(x) <= 1.259e-8 */
96 math_check_force_underflow_nonneg (w
);
101 /* (II) The case 1.259e-8 < abs(x) <= 0.0608 */
107 t2
= d9
.d
+ x2
* d11
.d
;
113 if ((y
= x
+ (t2
- u1
.d
* t2
)) == x
+ (t2
+ u1
.d
* t2
))
120 c1
= a25
.d
+ x2
* a27
.d
;
121 c1
= a23
.d
+ x2
* c1
;
122 c1
= a21
.d
+ x2
* c1
;
123 c1
= a19
.d
+ x2
* c1
;
124 c1
= a17
.d
+ x2
* c1
;
125 c1
= a15
.d
+ x2
* c1
;
128 EMULV (x
, x
, x2
, xx2
, t1
, t2
, t3
, t4
, t5
);
129 ADD2 (a13
.d
, aa13
.d
, c1
, 0.0, c2
, cc2
, t1
, t2
);
130 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
131 ADD2 (a11
.d
, aa11
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
132 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
133 ADD2 (a9
.d
, aa9
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
134 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
135 ADD2 (a7
.d
, aa7
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
136 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
137 ADD2 (a5
.d
, aa5
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
138 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
139 ADD2 (a3
.d
, aa3
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
140 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
141 MUL2 (x
, 0.0, c1
, cc1
, c2
, cc2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
142 ADD2 (x
, 0.0, c2
, cc2
, c1
, cc1
, t1
, t2
);
143 if ((y
= c1
+ (cc1
- u2
.d
* c1
)) == c1
+ (cc1
+ u2
.d
* c1
))
152 /* (III) The case 0.0608 < abs(x) <= 0.787 */
156 i
= ((int) (mfftnhf
.d
+ TWO8
* w
));
159 s
= (x
< 0.0) ? -1 : 1;
160 pz
= z
+ z
* z2
* (e0
.d
+ z2
* e1
.d
);
163 t2
= pz
* (gi
+ fi
) / (gi
- pz
);
164 if ((y
= fi
+ (t2
- fi
* u3
.d
)) == fi
+ (t2
+ fi
* u3
.d
))
169 t3
= (t2
< 0.0) ? -t2
: t2
;
170 t4
= fi
* ua3
.d
+ t3
* ub3
.d
;
171 if ((y
= fi
+ (t2
- t4
)) == fi
+ (t2
+ t4
))
179 c1
= z2
* (a7
.d
+ z2
* (a9
.d
+ z2
* a11
.d
));
180 EMULV (z
, z
, z2
, zz2
, t1
, t2
, t3
, t4
, t5
);
181 ADD2 (a5
.d
, aa5
.d
, c1
, 0.0, c2
, cc2
, t1
, t2
);
182 MUL2 (z2
, zz2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
183 ADD2 (a3
.d
, aa3
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
184 MUL2 (z2
, zz2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
185 MUL2 (z
, 0.0, c1
, cc1
, c2
, cc2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
186 ADD2 (z
, 0.0, c2
, cc2
, c1
, cc1
, t1
, t2
);
188 ADD2 (fi
, ffi
, c1
, cc1
, c2
, cc2
, t1
, t2
);
189 MUL2 (fi
, ffi
, c1
, cc1
, c3
, cc3
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
190 SUB2 (1.0, 0.0, c3
, cc3
, c1
, cc1
, t1
, t2
);
191 DIV2 (c2
, cc2
, c1
, cc1
, c3
, cc3
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
, t9
,
194 if ((y
= c3
+ (cc3
- u4
.d
* c3
)) == c3
+ (cc3
+ u4
.d
* c3
))
203 /* (---) The case 0.787 < abs(x) <= 25 */
206 /* Range reduction by algorithm i */
207 t
= (x
* hpinv
.d
+ toint
.d
);
210 t1
= (x
- xn
* mp1
.d
) - xn
* mp2
.d
;
211 n
= v
.i
[LOW_HALF
] & 0x00000001;
228 /* (IV),(V) The case 0.787 < abs(x) <= 25, abs(y) <= 1e-7 */
235 /* (VI) The case 0.787 < abs(x) <= 25, 1e-7 < abs(y) <= 0.0608 */
239 t2
= d9
.d
+ a2
* d11
.d
;
243 t2
= da
+ a
* a2
* t2
;
247 /* First stage -cot */
249 DIV2 (1.0, 0.0, b
, db
, c
, dc
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
,
251 if ((y
= c
+ (dc
- u6
.d
* c
)) == c
+ (dc
+ u6
.d
* c
))
259 /* First stage tan */
260 if ((y
= a
+ (t2
- u5
.d
* a
)) == a
+ (t2
+ u5
.d
* a
))
267 /* Range reduction by algorithm ii */
268 t
= (x
* hpinv
.d
+ toint
.d
);
271 t1
= (x
- xn
* mp1
.d
) - xn
* mp2
.d
;
272 n
= v
.i
[LOW_HALF
] & 0x00000001;
278 da
= ((t
- a
) - t1
) + da
;
281 EADD (a
, da
, t1
, t2
);
284 MUL2 (a
, da
, a
, da
, x2
, xx2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
286 c1
= a25
.d
+ x2
* a27
.d
;
287 c1
= a23
.d
+ x2
* c1
;
288 c1
= a21
.d
+ x2
* c1
;
289 c1
= a19
.d
+ x2
* c1
;
290 c1
= a17
.d
+ x2
* c1
;
291 c1
= a15
.d
+ x2
* c1
;
294 ADD2 (a13
.d
, aa13
.d
, c1
, 0.0, c2
, cc2
, t1
, t2
);
295 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
296 ADD2 (a11
.d
, aa11
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
297 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
298 ADD2 (a9
.d
, aa9
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
299 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
300 ADD2 (a7
.d
, aa7
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
301 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
302 ADD2 (a5
.d
, aa5
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
303 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
304 ADD2 (a3
.d
, aa3
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
305 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
306 MUL2 (a
, da
, c1
, cc1
, c2
, cc2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
307 ADD2 (a
, da
, c2
, cc2
, c1
, cc1
, t1
, t2
);
311 /* Second stage -cot */
312 DIV2 (1.0, 0.0, c1
, cc1
, c2
, cc2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
,
314 if ((y
= c2
+ (cc2
- u8
.d
* c2
)) == c2
+ (cc2
+ u8
.d
* c2
))
322 /* Second stage tan */
323 if ((y
= c1
+ (cc1
- u7
.d
* c1
)) == c1
+ (cc1
+ u7
.d
* c1
))
333 /* (VII) The case 0.787 < abs(x) <= 25, 0.0608 < abs(y) <= 0.787 */
336 i
= ((int) (mfftnhf
.d
+ TWO8
* ya
));
337 z
= (z0
= (ya
- xfg
[i
][0].d
)) + yya
;
339 pz
= z
+ z
* z2
* (e0
.d
+ z2
* e1
.d
);
346 t2
= pz
* (fi
+ gi
) / (fi
+ pz
);
347 if ((y
= gi
- (t2
- gi
* u10
.d
)) == gi
- (t2
+ gi
* u10
.d
))
352 t3
= (t2
< 0.0) ? -t2
: t2
;
353 t4
= gi
* ua10
.d
+ t3
* ub10
.d
;
354 if ((y
= gi
- (t2
- t4
)) == gi
- (t2
+ t4
))
363 t2
= pz
* (gi
+ fi
) / (gi
- pz
);
364 if ((y
= fi
+ (t2
- fi
* u9
.d
)) == fi
+ (t2
+ fi
* u9
.d
))
369 t3
= (t2
< 0.0) ? -t2
: t2
;
370 t4
= fi
* ua9
.d
+ t3
* ub9
.d
;
371 if ((y
= fi
+ (t2
- t4
)) == fi
+ (t2
+ t4
))
380 EADD (z0
, yya
, z
, zz
)
381 MUL2 (z
, zz
, z
, zz
, z2
, zz2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
382 c1
= z2
* (a7
.d
+ z2
* (a9
.d
+ z2
* a11
.d
));
383 ADD2 (a5
.d
, aa5
.d
, c1
, 0.0, c2
, cc2
, t1
, t2
);
384 MUL2 (z2
, zz2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
385 ADD2 (a3
.d
, aa3
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
386 MUL2 (z2
, zz2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
387 MUL2 (z
, zz
, c1
, cc1
, c2
, cc2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
388 ADD2 (z
, zz
, c2
, cc2
, c1
, cc1
, t1
, t2
);
390 ADD2 (fi
, ffi
, c1
, cc1
, c2
, cc2
, t1
, t2
);
391 MUL2 (fi
, ffi
, c1
, cc1
, c3
, cc3
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
392 SUB2 (1.0, 0.0, c3
, cc3
, c1
, cc1
, t1
, t2
);
397 DIV2 (c1
, cc1
, c2
, cc2
, c3
, cc3
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
, t9
,
399 if ((y
= c3
+ (cc3
- u12
.d
* c3
)) == c3
+ (cc3
+ u12
.d
* c3
))
408 DIV2 (c2
, cc2
, c1
, cc1
, c3
, cc3
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
, t9
,
410 if ((y
= c3
+ (cc3
- u11
.d
* c3
)) == c3
+ (cc3
+ u11
.d
* c3
))
421 /* (---) The case 25 < abs(x) <= 1e8 */
424 /* Range reduction by algorithm ii */
425 t
= (x
* hpinv
.d
+ toint
.d
);
428 t1
= (x
- xn
* mp1
.d
) - xn
* mp2
.d
;
429 n
= v
.i
[LOW_HALF
] & 0x00000001;
435 da
= ((t
- a
) - t1
) + da
;
436 EADD (a
, da
, t1
, t2
);
452 /* (+++) The case 25 < abs(x) <= 1e8, abs(y) <= 1e-7 */
459 /* (VIII) The case 25 < abs(x) <= 1e8, 1e-7 < abs(y) <= 0.0608 */
463 t2
= d9
.d
+ a2
* d11
.d
;
467 t2
= da
+ a
* a2
* t2
;
471 /* First stage -cot */
473 DIV2 (1.0, 0.0, b
, db
, c
, dc
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
,
475 if ((y
= c
+ (dc
- u14
.d
* c
)) == c
+ (dc
+ u14
.d
* c
))
483 /* First stage tan */
484 if ((y
= a
+ (t2
- u13
.d
* a
)) == a
+ (t2
+ u13
.d
* a
))
492 MUL2 (a
, da
, a
, da
, x2
, xx2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
493 c1
= a25
.d
+ x2
* a27
.d
;
494 c1
= a23
.d
+ x2
* c1
;
495 c1
= a21
.d
+ x2
* c1
;
496 c1
= a19
.d
+ x2
* c1
;
497 c1
= a17
.d
+ x2
* c1
;
498 c1
= a15
.d
+ x2
* c1
;
501 ADD2 (a13
.d
, aa13
.d
, c1
, 0.0, c2
, cc2
, t1
, t2
);
502 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
503 ADD2 (a11
.d
, aa11
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
504 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
505 ADD2 (a9
.d
, aa9
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
506 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
507 ADD2 (a7
.d
, aa7
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
508 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
509 ADD2 (a5
.d
, aa5
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
510 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
511 ADD2 (a3
.d
, aa3
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
512 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
513 MUL2 (a
, da
, c1
, cc1
, c2
, cc2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
514 ADD2 (a
, da
, c2
, cc2
, c1
, cc1
, t1
, t2
);
518 /* Second stage -cot */
519 DIV2 (1.0, 0.0, c1
, cc1
, c2
, cc2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
,
521 if ((y
= c2
+ (cc2
- u16
.d
* c2
)) == c2
+ (cc2
+ u16
.d
* c2
))
529 /* Second stage tan */
530 if ((y
= c1
+ (cc1
- u15
.d
* c1
)) == c1
+ (cc1
+ u15
.d
* c1
))
540 /* (IX) The case 25 < abs(x) <= 1e8, 0.0608 < abs(y) <= 0.787 */
542 i
= ((int) (mfftnhf
.d
+ TWO8
* ya
));
543 z
= (z0
= (ya
- xfg
[i
][0].d
)) + yya
;
545 pz
= z
+ z
* z2
* (e0
.d
+ z2
* e1
.d
);
552 t2
= pz
* (fi
+ gi
) / (fi
+ pz
);
553 if ((y
= gi
- (t2
- gi
* u18
.d
)) == gi
- (t2
+ gi
* u18
.d
))
558 t3
= (t2
< 0.0) ? -t2
: t2
;
559 t4
= gi
* ua18
.d
+ t3
* ub18
.d
;
560 if ((y
= gi
- (t2
- t4
)) == gi
- (t2
+ t4
))
569 t2
= pz
* (gi
+ fi
) / (gi
- pz
);
570 if ((y
= fi
+ (t2
- fi
* u17
.d
)) == fi
+ (t2
+ fi
* u17
.d
))
575 t3
= (t2
< 0.0) ? -t2
: t2
;
576 t4
= fi
* ua17
.d
+ t3
* ub17
.d
;
577 if ((y
= fi
+ (t2
- t4
)) == fi
+ (t2
+ t4
))
586 EADD (z0
, yya
, z
, zz
);
587 MUL2 (z
, zz
, z
, zz
, z2
, zz2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
588 c1
= z2
* (a7
.d
+ z2
* (a9
.d
+ z2
* a11
.d
));
589 ADD2 (a5
.d
, aa5
.d
, c1
, 0.0, c2
, cc2
, t1
, t2
);
590 MUL2 (z2
, zz2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
591 ADD2 (a3
.d
, aa3
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
592 MUL2 (z2
, zz2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
593 MUL2 (z
, zz
, c1
, cc1
, c2
, cc2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
594 ADD2 (z
, zz
, c2
, cc2
, c1
, cc1
, t1
, t2
);
596 ADD2 (fi
, ffi
, c1
, cc1
, c2
, cc2
, t1
, t2
);
597 MUL2 (fi
, ffi
, c1
, cc1
, c3
, cc3
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
598 SUB2 (1.0, 0.0, c3
, cc3
, c1
, cc1
, t1
, t2
);
603 DIV2 (c1
, cc1
, c2
, cc2
, c3
, cc3
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
, t9
,
605 if ((y
= c3
+ (cc3
- u20
.d
* c3
)) == c3
+ (cc3
+ u20
.d
* c3
))
614 DIV2 (c2
, cc2
, c1
, cc1
, c3
, cc3
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
, t9
,
616 if ((y
= c3
+ (cc3
- u19
.d
* c3
)) == c3
+ (cc3
+ u19
.d
* c3
))
626 /* (---) The case 1e8 < abs(x) < 2**1024 */
627 /* Range reduction by algorithm iii */
628 n
= (__branred (x
, &a
, &da
)) & 0x00000001;
629 EADD (a
, da
, t1
, t2
);
645 /* (+++) The case 1e8 < abs(x) < 2**1024, abs(y) <= 1e-7 */
652 /* (X) The case 1e8 < abs(x) < 2**1024, 1e-7 < abs(y) <= 0.0608 */
656 t2
= d9
.d
+ a2
* d11
.d
;
660 t2
= da
+ a
* a2
* t2
;
663 /* First stage -cot */
665 DIV2 (1.0, 0.0, b
, db
, c
, dc
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
, t9
,
667 if ((y
= c
+ (dc
- u22
.d
* c
)) == c
+ (dc
+ u22
.d
* c
))
675 /* First stage tan */
676 if ((y
= a
+ (t2
- u21
.d
* a
)) == a
+ (t2
+ u21
.d
* a
))
684 /* Reduction by algorithm iv */
686 n
= (__mpranred (x
, &mpa
, p
)) & 0x00000001;
687 __mp_dbl (&mpa
, &a
, p
);
688 __dbl_mp (a
, &mpt1
, p
);
689 __sub (&mpa
, &mpt1
, &mpt2
, p
);
690 __mp_dbl (&mpt2
, &da
, p
);
692 MUL2 (a
, da
, a
, da
, x2
, xx2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
694 c1
= a25
.d
+ x2
* a27
.d
;
695 c1
= a23
.d
+ x2
* c1
;
696 c1
= a21
.d
+ x2
* c1
;
697 c1
= a19
.d
+ x2
* c1
;
698 c1
= a17
.d
+ x2
* c1
;
699 c1
= a15
.d
+ x2
* c1
;
702 ADD2 (a13
.d
, aa13
.d
, c1
, 0.0, c2
, cc2
, t1
, t2
);
703 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
704 ADD2 (a11
.d
, aa11
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
705 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
706 ADD2 (a9
.d
, aa9
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
707 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
708 ADD2 (a7
.d
, aa7
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
709 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
710 ADD2 (a5
.d
, aa5
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
711 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
712 ADD2 (a3
.d
, aa3
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
713 MUL2 (x2
, xx2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
714 MUL2 (a
, da
, c1
, cc1
, c2
, cc2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
715 ADD2 (a
, da
, c2
, cc2
, c1
, cc1
, t1
, t2
);
719 /* Second stage -cot */
720 DIV2 (1.0, 0.0, c1
, cc1
, c2
, cc2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
,
722 if ((y
= c2
+ (cc2
- u24
.d
* c2
)) == c2
+ (cc2
+ u24
.d
* c2
))
730 /* Second stage tan */
731 if ((y
= c1
+ (cc1
- u23
.d
* c1
)) == c1
+ (cc1
+ u23
.d
* c1
))
741 /* (XI) The case 1e8 < abs(x) < 2**1024, 0.0608 < abs(y) <= 0.787 */
743 i
= ((int) (mfftnhf
.d
+ TWO8
* ya
));
744 z
= (z0
= (ya
- xfg
[i
][0].d
)) + yya
;
746 pz
= z
+ z
* z2
* (e0
.d
+ z2
* e1
.d
);
753 t2
= pz
* (fi
+ gi
) / (fi
+ pz
);
754 if ((y
= gi
- (t2
- gi
* u26
.d
)) == gi
- (t2
+ gi
* u26
.d
))
759 t3
= (t2
< 0.0) ? -t2
: t2
;
760 t4
= gi
* ua26
.d
+ t3
* ub26
.d
;
761 if ((y
= gi
- (t2
- t4
)) == gi
- (t2
+ t4
))
770 t2
= pz
* (gi
+ fi
) / (gi
- pz
);
771 if ((y
= fi
+ (t2
- fi
* u25
.d
)) == fi
+ (t2
+ fi
* u25
.d
))
776 t3
= (t2
< 0.0) ? -t2
: t2
;
777 t4
= fi
* ua25
.d
+ t3
* ub25
.d
;
778 if ((y
= fi
+ (t2
- t4
)) == fi
+ (t2
+ t4
))
787 EADD (z0
, yya
, z
, zz
);
788 MUL2 (z
, zz
, z
, zz
, z2
, zz2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
789 c1
= z2
* (a7
.d
+ z2
* (a9
.d
+ z2
* a11
.d
));
790 ADD2 (a5
.d
, aa5
.d
, c1
, 0.0, c2
, cc2
, t1
, t2
);
791 MUL2 (z2
, zz2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
792 ADD2 (a3
.d
, aa3
.d
, c1
, cc1
, c2
, cc2
, t1
, t2
);
793 MUL2 (z2
, zz2
, c2
, cc2
, c1
, cc1
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
794 MUL2 (z
, zz
, c1
, cc1
, c2
, cc2
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
795 ADD2 (z
, zz
, c2
, cc2
, c1
, cc1
, t1
, t2
);
797 ADD2 (fi
, ffi
, c1
, cc1
, c2
, cc2
, t1
, t2
);
798 MUL2 (fi
, ffi
, c1
, cc1
, c3
, cc3
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
);
799 SUB2 (1.0, 0.0, c3
, cc3
, c1
, cc1
, t1
, t2
);
804 DIV2 (c1
, cc1
, c2
, cc2
, c3
, cc3
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
, t9
,
806 if ((y
= c3
+ (cc3
- u28
.d
* c3
)) == c3
+ (cc3
+ u28
.d
* c3
))
815 DIV2 (c2
, cc2
, c1
, cc1
, c3
, cc3
, t1
, t2
, t3
, t4
, t5
, t6
, t7
, t8
, t9
,
817 if ((y
= c3
+ (cc3
- u27
.d
* c3
)) == c3
+ (cc3
+ u27
.d
* c3
))
830 /* multiple precision stage */
831 /* Convert x to multi precision number,compute tan(x) by mptan() routine */
832 /* and converts result back to double */
841 __mptan (x
, &mpy
, p
);
842 __mp_dbl (&mpy
, &y
, p
);
843 LIBC_PROBE (slowtan
, 2, &x
, &y
);
848 libm_alias_double (__tan
, tan
)