[PATCH] Fix cexp when both real and imaginary parts of x are non-finite

[Mark H Weaver <mhw at netris dot org>]

The default implementations of cexp, cexpf, and cexpl are broken when
the real part of x is infinite and the imaginary part is non-finite
(infinite or NaN).

In these cases, the complex phase of the result is indeterminate, and
thus the only two sensible results are 0.0 + i*0.0 (when the real part
of x is negative infinity) or NaN + i*NaN in all other cases.

Currently, cexp(plusinf + i*nonfinite) => plusinf + i*NaN, which is
incorrect because it implies that the sign of the real part of the
result is positive, and cexp(minusinf + i*minusinf) => 0.0 - i*0.0,
which is incorrect because it implies that the sign of the imaginary
part of the result is negative.

The patch below fixes these problems.

   Thanks,
     Mark


>From 4ee69779ae8e525eaf1640b1b186d7029c64d2a9 Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Tue, 31 Jul 2012 16:02:06 -0400
Subject: [PATCH] Fix cexp when both real and imaginary parts of x are
 non-finite

* math/s_cexp.c (__cexp), math/s_cexpf.c (__cexpf), math/s_cexpl.c
  (__cexpl): For non-finite b, e^(+inf + b*i) == NaN + i*NaN,
  and e^(-inf + b*i) == 0.0 + i*0.0.
---
 math/s_cexp.c  |    4 ++--
 math/s_cexpf.c |    4 ++--
 math/s_cexpl.c |    4 ++--
 3 files changed, 6 insertions(+), 6 deletions(-)

diff --git a/math/s_cexp.c b/math/s_cexp.c
index 1d7a5a2..1b32da2 100644
--- a/math/s_cexp.c
+++ b/math/s_cexp.c
@@ -103,7 +103,7 @@ __cexp (__complex__ double x)
 	}
       else if (signbit (__real__ x) == 0)
 	{
-	  __real__ retval = HUGE_VAL;
+	  __real__ retval = __nan ("");
 	  __imag__ retval = __nan ("");
 
 	  if (icls == FP_INFINITE)
@@ -112,7 +112,7 @@ __cexp (__complex__ double x)
       else
 	{
 	  __real__ retval = 0.0;
-	  __imag__ retval = __copysign (0.0, __imag__ x);
+	  __imag__ retval = 0.0;
 	}
     }
   else
diff --git a/math/s_cexpf.c b/math/s_cexpf.c
index 4aa9765..888d8a0 100644
--- a/math/s_cexpf.c
+++ b/math/s_cexpf.c
@@ -103,7 +103,7 @@ __cexpf (__complex__ float x)
 	}
       else if (signbit (__real__ x) == 0)
 	{
-	  __real__ retval = HUGE_VALF;
+	  __real__ retval = __nanf ("");
 	  __imag__ retval = __nanf ("");
 
 	  if (icls == FP_INFINITE)
@@ -112,7 +112,7 @@ __cexpf (__complex__ float x)
       else
 	{
 	  __real__ retval = 0.0;
-	  __imag__ retval = __copysignf (0.0, __imag__ x);
+	  __imag__ retval = 0.0;
 	}
     }
   else
diff --git a/math/s_cexpl.c b/math/s_cexpl.c
index 2568249..269b5b4 100644
--- a/math/s_cexpl.c
+++ b/math/s_cexpl.c
@@ -103,7 +103,7 @@ __cexpl (__complex__ long double x)
 	}
       else if (signbit (__real__ x) == 0)
 	{
-	  __real__ retval = HUGE_VALL;
+	  __real__ retval = __nanl ("");
 	  __imag__ retval = __nanl ("");
 
 	  if (icls == FP_INFINITE)
@@ -112,7 +112,7 @@ __cexpl (__complex__ long double x)
       else
 	{
 	  __real__ retval = 0.0;
-	  __imag__ retval = __copysignl (0.0, __imag__ x);
+	  __imag__ retval = 0.0;
 	}
     }
   else
-- 
1.7.9.5