1 /* Complex sine function for complex __float128.
2    Copyright (C) 1997-2012 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 "quadmath-imp.h"
21 
22 #ifdef HAVE_FENV_H
23 # include <fenv.h>
24 #endif
25 
26 
27 __complex128
csinq(__complex128 x)28 csinq (__complex128 x)
29 {
30   __complex128 retval;
31   int negate = signbitq (__real__ x);
32   int rcls = fpclassifyq (__real__ x);
33   int icls = fpclassifyq (__imag__ x);
34 
35   __real__ x = fabsq (__real__ x);
36 
37   if (__builtin_expect (icls >= QUADFP_ZERO, 1))
38     {
39       /* Imaginary part is finite.  */
40       if (__builtin_expect (rcls >= QUADFP_ZERO, 1))
41 	{
42 	  /* Real part is finite.  */
43 	  const int t = (int) ((FLT128_MAX_EXP - 1) * M_LN2q);
44 	  __float128 sinix, cosix;
45 
46 	  if (__builtin_expect (rcls != QUADFP_SUBNORMAL, 1))
47 	    {
48 	      sincosq (__real__ x, &sinix, &cosix);
49 	    }
50 	  else
51 	    {
52 	      sinix = __real__ x;
53 	      cosix = 1.0Q;
54 	    }
55 
56 	  if (fabsq (__imag__ x) > t)
57 	    {
58 	      __float128 exp_t = expq (t);
59 	      __float128 ix = fabsq (__imag__ x);
60 	      if (signbitq (__imag__ x))
61 		cosix = -cosix;
62 	      ix -= t;
63 	      sinix *= exp_t / 2.0Q;
64 	      cosix *= exp_t / 2.0Q;
65 	      if (ix > t)
66 		{
67 		  ix -= t;
68 		  sinix *= exp_t;
69 		  cosix *= exp_t;
70 		}
71 	      if (ix > t)
72 		{
73 		  /* Overflow (original imaginary part of x > 3t).  */
74 		  __real__ retval = FLT128_MAX * sinix;
75 		  __imag__ retval = FLT128_MAX * cosix;
76 		}
77 	      else
78 		{
79 		  __float128 exp_val = expq (ix);
80 		  __real__ retval = exp_val * sinix;
81 		  __imag__ retval = exp_val * cosix;
82 		}
83 	    }
84 	  else
85 	    {
86 	      __real__ retval = coshq (__imag__ x) * sinix;
87 	      __imag__ retval = sinhq (__imag__ x) * cosix;
88 	    }
89 
90 	  if (negate)
91 	    __real__ retval = -__real__ retval;
92 	}
93       else
94 	{
95 	  if (icls == QUADFP_ZERO)
96 	    {
97 	      /* Imaginary part is 0.0.  */
98 	      __real__ retval = nanq ("");
99 	      __imag__ retval = __imag__ x;
100 
101 #ifdef HAVE_FENV_H
102 	      if (rcls == QUADFP_INFINITE)
103 		feraiseexcept (FE_INVALID);
104 #endif
105 	    }
106 	  else
107 	    {
108 	      __real__ retval = nanq ("");
109 	      __imag__ retval = nanq ("");
110 
111 #ifdef HAVE_FENV_H
112 	      feraiseexcept (FE_INVALID);
113 #endif
114 	    }
115 	}
116     }
117   else if (icls == QUADFP_INFINITE)
118     {
119       /* Imaginary part is infinite.  */
120       if (rcls == QUADFP_ZERO)
121 	{
122 	  /* Real part is 0.0.  */
123 	  __real__ retval = copysignq (0.0Q, negate ? -1.0Q : 1.0Q);
124 	  __imag__ retval = __imag__ x;
125 	}
126       else if (rcls > QUADFP_ZERO)
127 	{
128 	  /* Real part is finite.  */
129 	  __float128 sinix, cosix;
130 
131 	  if (__builtin_expect (rcls != QUADFP_SUBNORMAL, 1))
132 	    {
133 	      sincosq (__real__ x, &sinix, &cosix);
134 	    }
135 	  else
136 	    {
137 	      sinix = __real__ x;
138 	      cosix = 1.0;
139 	    }
140 
141 	  __real__ retval = copysignq (HUGE_VALQ, sinix);
142 	  __imag__ retval = copysignq (HUGE_VALQ, cosix);
143 
144 	  if (negate)
145 	    __real__ retval = -__real__ retval;
146 	  if (signbitq (__imag__ x))
147 	    __imag__ retval = -__imag__ retval;
148 	}
149       else
150 	{
151 	  /* The addition raises the invalid exception.  */
152 	  __real__ retval = nanq ("");
153 	  __imag__ retval = HUGE_VALQ;
154 
155 #ifdef HAVE_FENV_H
156 	  if (rcls == QUADFP_INFINITE)
157 	    feraiseexcept (FE_INVALID);
158 #endif
159 	}
160     }
161   else
162     {
163       if (rcls == QUADFP_ZERO)
164 	__real__ retval = copysignq (0.0Q, negate ? -1.0Q : 1.0Q);
165       else
166 	__real__ retval = nanq ("");
167       __imag__ retval = nanq ("");
168     }
169 
170   return retval;
171 }
172