xref: /netbsd-src/external/gpl3/gcc.old/dist/libquadmath/math/ctanq.c (revision 627f7eb200a4419d89b531d55fccd2ee3ffdcde0)
1*627f7eb2Smrg /* Complex tangent function for a complex float type.
2*627f7eb2Smrg    Copyright (C) 1997-2018 Free Software Foundation, Inc.
3*627f7eb2Smrg    This file is part of the GNU C Library.
4*627f7eb2Smrg    Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
5*627f7eb2Smrg 
6*627f7eb2Smrg    The GNU C Library is free software; you can redistribute it and/or
7*627f7eb2Smrg    modify it under the terms of the GNU Lesser General Public
8*627f7eb2Smrg    License as published by the Free Software Foundation; either
9*627f7eb2Smrg    version 2.1 of the License, or (at your option) any later version.
10*627f7eb2Smrg 
11*627f7eb2Smrg    The GNU C Library is distributed in the hope that it will be useful,
12*627f7eb2Smrg    but WITHOUT ANY WARRANTY; without even the implied warranty of
13*627f7eb2Smrg    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
14*627f7eb2Smrg    Lesser General Public License for more details.
15*627f7eb2Smrg 
16*627f7eb2Smrg    You should have received a copy of the GNU Lesser General Public
17*627f7eb2Smrg    License along with the GNU C Library; if not, see
18*627f7eb2Smrg    <http://www.gnu.org/licenses/>.  */
19*627f7eb2Smrg 
20*627f7eb2Smrg #include "quadmath-imp.h"
21*627f7eb2Smrg 
22*627f7eb2Smrg __complex128
ctanq(__complex128 x)23*627f7eb2Smrg ctanq (__complex128 x)
24*627f7eb2Smrg {
25*627f7eb2Smrg   __complex128 res;
26*627f7eb2Smrg 
27*627f7eb2Smrg   if (__glibc_unlikely (!finiteq (__real__ x) || !finiteq (__imag__ x)))
28*627f7eb2Smrg     {
29*627f7eb2Smrg       if (isinfq (__imag__ x))
30*627f7eb2Smrg 	{
31*627f7eb2Smrg 	  if (finiteq (__real__ x) && fabsq (__real__ x) > 1)
32*627f7eb2Smrg 	    {
33*627f7eb2Smrg 	      __float128 sinrx, cosrx;
34*627f7eb2Smrg 	      sincosq (__real__ x, &sinrx, &cosrx);
35*627f7eb2Smrg 	      __real__ res = copysignq (0, sinrx * cosrx);
36*627f7eb2Smrg 	    }
37*627f7eb2Smrg 	  else
38*627f7eb2Smrg 	    __real__ res = copysignq (0, __real__ x);
39*627f7eb2Smrg 	  __imag__ res = copysignq (1, __imag__ x);
40*627f7eb2Smrg 	}
41*627f7eb2Smrg       else if (__real__ x == 0)
42*627f7eb2Smrg 	{
43*627f7eb2Smrg 	  res = x;
44*627f7eb2Smrg 	}
45*627f7eb2Smrg       else
46*627f7eb2Smrg 	{
47*627f7eb2Smrg 	  __real__ res = nanq ("");
48*627f7eb2Smrg 	  if (__imag__ x == 0)
49*627f7eb2Smrg 	    __imag__ res = __imag__ x;
50*627f7eb2Smrg 	  else
51*627f7eb2Smrg 	    __imag__ res = nanq ("");
52*627f7eb2Smrg 
53*627f7eb2Smrg 	  if (isinfq (__real__ x))
54*627f7eb2Smrg 	    feraiseexcept (FE_INVALID);
55*627f7eb2Smrg 	}
56*627f7eb2Smrg     }
57*627f7eb2Smrg   else
58*627f7eb2Smrg     {
59*627f7eb2Smrg       __float128 sinrx, cosrx;
60*627f7eb2Smrg       __float128 den;
61*627f7eb2Smrg       const int t = (int) ((FLT128_MAX_EXP - 1) * M_LN2q / 2);
62*627f7eb2Smrg 
63*627f7eb2Smrg       /* tan(x+iy) = (sin(2x) + i*sinh(2y))/(cos(2x) + cosh(2y))
64*627f7eb2Smrg 	 = (sin(x)*cos(x) + i*sinh(y)*cosh(y)/(cos(x)^2 + sinh(y)^2). */
65*627f7eb2Smrg 
66*627f7eb2Smrg       if (__glibc_likely (fabsq (__real__ x) > FLT128_MIN))
67*627f7eb2Smrg 	{
68*627f7eb2Smrg 	  sincosq (__real__ x, &sinrx, &cosrx);
69*627f7eb2Smrg 	}
70*627f7eb2Smrg       else
71*627f7eb2Smrg 	{
72*627f7eb2Smrg 	  sinrx = __real__ x;
73*627f7eb2Smrg 	  cosrx = 1;
74*627f7eb2Smrg 	}
75*627f7eb2Smrg 
76*627f7eb2Smrg       if (fabsq (__imag__ x) > t)
77*627f7eb2Smrg 	{
78*627f7eb2Smrg 	  /* Avoid intermediate overflow when the real part of the
79*627f7eb2Smrg 	     result may be subnormal.  Ignoring negligible terms, the
80*627f7eb2Smrg 	     imaginary part is +/- 1, the real part is
81*627f7eb2Smrg 	     sin(x)*cos(x)/sinh(y)^2 = 4*sin(x)*cos(x)/exp(2y).  */
82*627f7eb2Smrg 	  __float128 exp_2t = expq (2 * t);
83*627f7eb2Smrg 
84*627f7eb2Smrg 	  __imag__ res = copysignq (1, __imag__ x);
85*627f7eb2Smrg 	  __real__ res = 4 * sinrx * cosrx;
86*627f7eb2Smrg 	  __imag__ x = fabsq (__imag__ x);
87*627f7eb2Smrg 	  __imag__ x -= t;
88*627f7eb2Smrg 	  __real__ res /= exp_2t;
89*627f7eb2Smrg 	  if (__imag__ x > t)
90*627f7eb2Smrg 	    {
91*627f7eb2Smrg 	      /* Underflow (original imaginary part of x has absolute
92*627f7eb2Smrg 		 value > 2t).  */
93*627f7eb2Smrg 	      __real__ res /= exp_2t;
94*627f7eb2Smrg 	    }
95*627f7eb2Smrg 	  else
96*627f7eb2Smrg 	    __real__ res /= expq (2 * __imag__ x);
97*627f7eb2Smrg 	}
98*627f7eb2Smrg       else
99*627f7eb2Smrg 	{
100*627f7eb2Smrg 	  __float128 sinhix, coshix;
101*627f7eb2Smrg 	  if (fabsq (__imag__ x) > FLT128_MIN)
102*627f7eb2Smrg 	    {
103*627f7eb2Smrg 	      sinhix = sinhq (__imag__ x);
104*627f7eb2Smrg 	      coshix = coshq (__imag__ x);
105*627f7eb2Smrg 	    }
106*627f7eb2Smrg 	  else
107*627f7eb2Smrg 	    {
108*627f7eb2Smrg 	      sinhix = __imag__ x;
109*627f7eb2Smrg 	      coshix = 1;
110*627f7eb2Smrg 	    }
111*627f7eb2Smrg 
112*627f7eb2Smrg 	  if (fabsq (sinhix) > fabsq (cosrx) * FLT128_EPSILON)
113*627f7eb2Smrg 	    den = cosrx * cosrx + sinhix * sinhix;
114*627f7eb2Smrg 	  else
115*627f7eb2Smrg 	    den = cosrx * cosrx;
116*627f7eb2Smrg 	  __real__ res = sinrx * cosrx / den;
117*627f7eb2Smrg 	  __imag__ res = sinhix * coshix / den;
118*627f7eb2Smrg 	}
119*627f7eb2Smrg       math_check_force_underflow_complex (res);
120*627f7eb2Smrg     }
121*627f7eb2Smrg 
122*627f7eb2Smrg   return res;
123*627f7eb2Smrg }
124