Home | History | Annotate | Line # | Download | only in math
      1 /* Return value of complex exponential function for a float type.
      2    Copyright (C) 1997-2018 Free Software Foundation, Inc.
      3    This file is part of the GNU C Library.
      4    Contributed by Ulrich Drepper <drepper (at) 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 __complex128
     23 cexpq (__complex128 x)
     24 {
     25   __complex128 retval;
     26   int rcls = fpclassifyq (__real__ x);
     27   int icls = fpclassifyq (__imag__ x);
     28 
     29   if (__glibc_likely (rcls >= QUADFP_ZERO))
     30     {
     31       /* Real part is finite.  */
     32       if (__glibc_likely (icls >= QUADFP_ZERO))
     33 	{
     34 	  /* Imaginary part is finite.  */
     35 	  const int t = (int) ((FLT128_MAX_EXP - 1) * M_LN2q);
     36 	  __float128 sinix, cosix;
     37 
     38 	  if (__glibc_likely (fabsq (__imag__ x) > FLT128_MIN))
     39 	    {
     40 	      sincosq (__imag__ x, &sinix, &cosix);
     41 	    }
     42 	  else
     43 	    {
     44 	      sinix = __imag__ x;
     45 	      cosix = 1;
     46 	    }
     47 
     48 	  if (__real__ x > t)
     49 	    {
     50 	      __float128 exp_t = expq (t);
     51 	      __real__ x -= t;
     52 	      sinix *= exp_t;
     53 	      cosix *= exp_t;
     54 	      if (__real__ x > t)
     55 		{
     56 		  __real__ x -= t;
     57 		  sinix *= exp_t;
     58 		  cosix *= exp_t;
     59 		}
     60 	    }
     61 	  if (__real__ x > t)
     62 	    {
     63 	      /* Overflow (original real part of x > 3t).  */
     64 	      __real__ retval = FLT128_MAX * cosix;
     65 	      __imag__ retval = FLT128_MAX * sinix;
     66 	    }
     67 	  else
     68 	    {
     69 	      __float128 exp_val = expq (__real__ x);
     70 	      __real__ retval = exp_val * cosix;
     71 	      __imag__ retval = exp_val * sinix;
     72 	    }
     73 	  math_check_force_underflow_complex (retval);
     74 	}
     75       else
     76 	{
     77 	  /* If the imaginary part is +-inf or NaN and the real part
     78 	     is not +-inf the result is NaN + iNaN.  */
     79 	  __real__ retval = nanq ("");
     80 	  __imag__ retval = nanq ("");
     81 
     82 	  feraiseexcept (FE_INVALID);
     83 	}
     84     }
     85   else if (__glibc_likely (rcls == QUADFP_INFINITE))
     86     {
     87       /* Real part is infinite.  */
     88       if (__glibc_likely (icls >= QUADFP_ZERO))
     89 	{
     90 	  /* Imaginary part is finite.  */
     91 	  __float128 value = signbitq (__real__ x) ? 0 : HUGE_VALQ;
     92 
     93 	  if (icls == QUADFP_ZERO)
     94 	    {
     95 	      /* Imaginary part is 0.0.  */
     96 	      __real__ retval = value;
     97 	      __imag__ retval = __imag__ x;
     98 	    }
     99 	  else
    100 	    {
    101 	      __float128 sinix, cosix;
    102 
    103 	      if (__glibc_likely (fabsq (__imag__ x) > FLT128_MIN))
    104 		{
    105 		  sincosq (__imag__ x, &sinix, &cosix);
    106 		}
    107 	      else
    108 		{
    109 		  sinix = __imag__ x;
    110 		  cosix = 1;
    111 		}
    112 
    113 	      __real__ retval = copysignq (value, cosix);
    114 	      __imag__ retval = copysignq (value, sinix);
    115 	    }
    116 	}
    117       else if (signbitq (__real__ x) == 0)
    118 	{
    119 	  __real__ retval = HUGE_VALQ;
    120 	  __imag__ retval = __imag__ x - __imag__ x;
    121 	}
    122       else
    123 	{
    124 	  __real__ retval = 0;
    125 	  __imag__ retval = copysignq (0, __imag__ x);
    126 	}
    127     }
    128   else
    129     {
    130       /* If the real part is NaN the result is NaN + iNaN unless the
    131 	 imaginary part is zero.  */
    132       __real__ retval = nanq ("");
    133       if (icls == QUADFP_ZERO)
    134 	__imag__ retval = __imag__ x;
    135       else
    136 	{
    137 	  __imag__ retval = nanq ("");
    138 
    139 	  if (rcls != QUADFP_NAN || icls != QUADFP_NAN)
    140 	    feraiseexcept (FE_INVALID);
    141 	}
    142     }
    143 
    144   return retval;
    145 }
    146