1 /* Complex hyperbolic tangent for float types.
2 Copyright (C) 1997-2022 Free Software Foundation, Inc.
3 This file is part of the GNU C Library.
4
5 The GNU C Library is free software; you can redistribute it and/or
6 modify it under the terms of the GNU Lesser General Public
7 License as published by the Free Software Foundation; either
8 version 2.1 of the License, or (at your option) any later version.
9
10 The GNU C Library is distributed in the hope that it will be useful,
11 but WITHOUT ANY WARRANTY; without even the implied warranty of
12 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13 Lesser General Public License for more details.
14
15 You should have received a copy of the GNU Lesser General Public
16 License along with the GNU C Library; if not, see
17 <https://www.gnu.org/licenses/>. */
18
19 #include <complex.h>
20 #include <fenv.h>
21 #include <math.h>
22 #include <math_private.h>
23 #include <math-underflow.h>
24 #include <float.h>
25
26 CFLOAT
M_DECL_FUNC(__ctanh)27 M_DECL_FUNC (__ctanh) (CFLOAT x)
28 {
29 CFLOAT res;
30
31 if (__glibc_unlikely (!isfinite (__real__ x) || !isfinite (__imag__ x)))
32 {
33 if (isinf (__real__ x))
34 {
35 __real__ res = M_COPYSIGN (1, __real__ x);
36 if (isfinite (__imag__ x) && M_FABS (__imag__ x) > 1)
37 {
38 FLOAT sinix, cosix;
39 M_SINCOS (__imag__ x, &sinix, &cosix);
40 __imag__ res = M_COPYSIGN (0, sinix * cosix);
41 }
42 else
43 __imag__ res = M_COPYSIGN (0, __imag__ x);
44 }
45 else if (__imag__ x == 0)
46 {
47 res = x;
48 }
49 else
50 {
51 if (__real__ x == 0)
52 __real__ res = __real__ x;
53 else
54 __real__ res = M_NAN;
55 __imag__ res = M_NAN;
56
57 if (isinf (__imag__ x))
58 feraiseexcept (FE_INVALID);
59 }
60 }
61 else
62 {
63 FLOAT sinix, cosix;
64 FLOAT den;
65 const int t = (int) ((M_MAX_EXP - 1) * M_MLIT (M_LN2) / 2);
66
67 /* tanh(x+iy) = (sinh(2x) + i*sin(2y))/(cosh(2x) + cos(2y))
68 = (sinh(x)*cosh(x) + i*sin(y)*cos(y))/(sinh(x)^2 + cos(y)^2). */
69
70 if (__glibc_likely (M_FABS (__imag__ x) > M_MIN))
71 {
72 M_SINCOS (__imag__ x, &sinix, &cosix);
73 }
74 else
75 {
76 sinix = __imag__ x;
77 cosix = 1;
78 }
79
80 if (M_FABS (__real__ x) > t)
81 {
82 /* Avoid intermediate overflow when the imaginary part of
83 the result may be subnormal. Ignoring negligible terms,
84 the real part is +/- 1, the imaginary part is
85 sin(y)*cos(y)/sinh(x)^2 = 4*sin(y)*cos(y)/exp(2x). */
86 FLOAT exp_2t = M_EXP (2 * t);
87
88 __real__ res = M_COPYSIGN (1, __real__ x);
89 __imag__ res = 4 * sinix * cosix;
90 __real__ x = M_FABS (__real__ x);
91 __real__ x -= t;
92 __imag__ res /= exp_2t;
93 if (__real__ x > t)
94 {
95 /* Underflow (original real part of x has absolute value
96 > 2t). */
97 __imag__ res /= exp_2t;
98 }
99 else
100 __imag__ res /= M_EXP (2 * __real__ x);
101 }
102 else
103 {
104 FLOAT sinhrx, coshrx;
105 if (M_FABS (__real__ x) > M_MIN)
106 {
107 sinhrx = M_SINH (__real__ x);
108 coshrx = M_COSH (__real__ x);
109 }
110 else
111 {
112 sinhrx = __real__ x;
113 coshrx = 1;
114 }
115
116 if (M_FABS (sinhrx) > M_FABS (cosix) * M_EPSILON)
117 den = sinhrx * sinhrx + cosix * cosix;
118 else
119 den = cosix * cosix;
120 __real__ res = sinhrx * coshrx / den;
121 __imag__ res = sinix * cosix / den;
122 }
123 math_check_force_underflow_complex (res);
124 }
125
126 return res;
127 }
128
129 declare_mgen_alias (__ctanh, ctanh)
130