1/* Return arc tangent of complex long double value.
2 Copyright (C) 1997-2016 Free Software Foundation, Inc.
3 This file is part of the GNU C Library.
4 Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
5
6 The GNU C Library is free software; you can redistribute it and/or
7 modify it under the terms of the GNU Lesser General Public
8 License as published by the Free Software Foundation; either
9 version 2.1 of the License, or (at your option) any later version.
10
11 The GNU C Library 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 GNU
14 Lesser General Public License for more details.
15
16 You should have received a copy of the GNU Lesser General Public
17 License along with the GNU C Library; if not, see
18 <http://www.gnu.org/licenses/>. */
19
20#include <complex.h>
21#include <math.h>
22#include <math_private.h>
23#include <float.h>
24
25/* To avoid spurious overflows, use this definition to treat IBM long
26 double as approximating an IEEE-style format. */
27#if LDBL_MANT_DIG == 106
28# undef LDBL_EPSILON
29# define LDBL_EPSILON 0x1p-106L
30#endif
31
32__complex__ long double
33__catanl (__complex__ long double x)
34{
35 __complex__ long double res;
36 int rcls = fpclassify (__real__ x);
37 int icls = fpclassify (__imag__ x);
38
39 if (__glibc_unlikely (rcls <= FP_INFINITE || icls <= FP_INFINITE))
40 {
41 if (rcls == FP_INFINITE)
42 {
43 __real__ res = __copysignl (M_PI_2l, __real__ x);
44 __imag__ res = __copysignl (0.0, __imag__ x);
45 }
46 else if (icls == FP_INFINITE)
47 {
48 if (rcls >= FP_ZERO)
49 __real__ res = __copysignl (M_PI_2l, __real__ x);
50 else
51 __real__ res = __nanl ("");
52 __imag__ res = __copysignl (0.0, __imag__ x);
53 }
54 else if (icls == FP_ZERO || icls == FP_INFINITE)
55 {
56 __real__ res = __nanl ("");
57 __imag__ res = __copysignl (0.0, __imag__ x);
58 }
59 else
60 {
61 __real__ res = __nanl ("");
62 __imag__ res = __nanl ("");
63 }
64 }
65 else if (__glibc_unlikely (rcls == FP_ZERO && icls == FP_ZERO))
66 {
67 res = x;
68 }
69 else
70 {
71 if (fabsl (__real__ x) >= 16.0L / LDBL_EPSILON
72 || fabsl (__imag__ x) >= 16.0L / LDBL_EPSILON)
73 {
74 __real__ res = __copysignl (M_PI_2l, __real__ x);
75 if (fabsl (__real__ x) <= 1.0L)
76 __imag__ res = 1.0L / __imag__ x;
77 else if (fabsl (__imag__ x) <= 1.0L)
78 __imag__ res = __imag__ x / __real__ x / __real__ x;
79 else
80 {
81 long double h = __ieee754_hypotl (__real__ x / 2.0L,
82 __imag__ x / 2.0L);
83 __imag__ res = __imag__ x / h / h / 4.0L;
84 }
85 }
86 else
87 {
88 long double den, absx, absy;
89
90 absx = fabsl (__real__ x);
91 absy = fabsl (__imag__ x);
92 if (absx < absy)
93 {
94 long double t = absx;
95 absx = absy;
96 absy = t;
97 }
98
99 if (absy < LDBL_EPSILON / 2.0L)
100 {
101 den = (1.0L - absx) * (1.0L + absx);
102 if (den == -0.0L)
103 den = 0.0L;
104 }
105 else if (absx >= 1.0L)
106 den = (1.0L - absx) * (1.0L + absx) - absy * absy;
107 else if (absx >= 0.75L || absy >= 0.5L)
108 den = -__x2y2m1l (absx, absy);
109 else
110 den = (1.0L - absx) * (1.0L + absx) - absy * absy;
111
112 __real__ res = 0.5L * __ieee754_atan2l (2.0L * __real__ x, den);
113
114 if (fabsl (__imag__ x) == 1.0L
115 && fabsl (__real__ x) < LDBL_EPSILON * LDBL_EPSILON)
116 __imag__ res = (__copysignl (0.5L, __imag__ x)
117 * (M_LN2l - __ieee754_logl (fabsl (__real__ x))));
118 else
119 {
120 long double r2 = 0.0L, num, f;
121
122 if (fabsl (__real__ x) >= LDBL_EPSILON * LDBL_EPSILON)
123 r2 = __real__ x * __real__ x;
124
125 num = __imag__ x + 1.0L;
126 num = r2 + num * num;
127
128 den = __imag__ x - 1.0L;
129 den = r2 + den * den;
130
131 f = num / den;
132 if (f < 0.5L)
133 __imag__ res = 0.25L * __ieee754_logl (f);
134 else
135 {
136 num = 4.0L * __imag__ x;
137 __imag__ res = 0.25L * __log1pl (num / den);
138 }
139 }
140 }
141
142 math_check_force_underflow_complex (res);
143 }
144
145 return res;
146}
147weak_alias (__catanl, catanl)
148