The following reply was made to PR bin/170206; it has been noted by GNATS.
From: Stephen Montgomery-Smith <[email protected]>
To: Bruce Evans <[email protected]>
Cc: [email protected], [email protected],
Stephen Montgomery-Smith <[email protected]>
Subject: Re: bin/170206: complex arcsinh, log, etc.
Date: Sat, 28 Jul 2012 21:29:42 -0500
This is a multi-part message in MIME format.
--------------050105020207030702040007
Content-Type: text/plain; charset=ISO-8859-1; format=flowed
Content-Transfer-Encoding: 7bit
On 07/28/2012 04:46 PM, Stephen Montgomery-Smith wrote:
> Here are some diffs to catrig.c so that it completely passes Peter
> Jeremy's program www.rulingia.com/~peter/ctest.c. That is, it seems to
> get all the infs and nans correct.
And I think I messed up these diffs as well. Can we try this instead?
--------------050105020207030702040007
Content-Type: text/x-diff;
name="ca.diff"
Content-Transfer-Encoding: 7bit
Content-Disposition: attachment;
filename="ca.diff"
--- catrig.c-old1 2012-07-28 15:00:46.000000000 -0500
+++ catrig.c 2012-07-28 21:10:20.000000000 -0500
@@ -89,7 +89,7 @@
S = hypot(x,y-1);
A = 0.5*(R + S);
- if (A < 10) {
+ if (A < 10 && isfinite(A)) {
fp = f(x,1+y,&fpuf);
fm = f(x,1-y,&fmuf);
if (fpuf == 1 && fmuf == 1) {
@@ -108,9 +108,24 @@
} else {
*rx = log1p(fp + fm + sqrt((fp+fm)*(A+1)));
}
- } else
+ } else if (isinf(y))
+ *rx = y;
+ else
*rx = log(A + sqrt(A*A-1));
+ if (!isfinite(y)) {
+ *B_good = 0;
+ if (isfinite(x)) *A2my2 = 0;
+ else if (isnan(x)) *A2my2 = x;
+ else *A2my2 = y;
+ return;
+ }
+ if (isnan(x) && y == 0) {
+ *B_good = 0;
+ *A2my2 = 1;
+ return;
+ }
+
*B = y/A; /* = 0.5*(R - S) */
*B_good = 1;
@@ -147,7 +162,7 @@
x = fabs(x);
y = fabs(y);
- if (cabs(z) > 1e20) {
+ if (cabs(z) > 1e20 && isfinite(x) && isfinite(y)) {
if (huge+x>one) { /* set inexact flag. */
if (sx == 0) return clog(2*z);
if (sx == 1) return -clog(-2*z);
@@ -164,8 +179,8 @@
else
ry = atan2(y,A2my2);
- if (sx == 1) rx = -rx;
- if (sy == 1) ry = -ry;
+ if (sx == 1) rx = copysign(rx, -1);
+ if (sy == 1) ry = copysign(ry, -1);
return cpack(rx,ry);
}
@@ -206,7 +221,7 @@
x = fabs(x);
y = fabs(y);
- if (cabs(z) > 1e20) {
+ if (cabs(z) > 1e20 && isfinite(x) && isfinite(y)) {
if (huge+x>one) { /* set inexact flag. */
w = clog(2*z);
if (signbit(cimag(w)) == 0)
@@ -232,7 +247,7 @@
rx = atan2(A2my2,-x);
}
- if (sy==0) ry = -ry;
+ if (sy==0) ry = copysign(ry, -1);
return cpack(rx,ry);
}
@@ -271,11 +286,23 @@
if (huge+x>one) /* set inexact flag. */
return z;
- if (cabs(z) > 1e20)
- if (huge+x>one) { /* set inexact flag. */
- if (signbit(x) == 0)
- return cpack(0,M_PI_2);
- return cpack(0,-M_PI_2);
+ if (isinf(x) && isnan(y))
+ return cpack(copysign(0, x), y);
+
+ if (isnan(x) && isinf(y))
+ return cpack(copysign(0, x), copysign(M_PI_2,y));
+
+ if (x == 0 && isnan(y))
+ return cpack(x, y);
+
+ if (isnan(x) || isnan(y))
+ return clog(z);
+
+ if (cabs(z) > 1e20) {
+ if (isinf(x) || isinf(y))
+ return cpack(copysign(0,x),copysign(M_PI_2,y));
+ if (huge+x>one) /* set inexact flag. */
+ return cpack(copysign(0,x),copysign(M_PI_2,y));
}
if (fabs(y) < 1e-100) {
@@ -290,13 +317,17 @@
if (hp < 0.5 || hm < 0.5)
rx = 0.25*(log(hp/hm));
+ else if (x == 0)
+ rx = x;
else if (x > 0)
rx = 0.25*log1p(4*x/hm);
else
rx = -0.25*log1p(-4*x/hp);
if (x==1 || x==-1) {
- if (signbit(y) == 0)
+ if (y==0)
+ ry = y;
+ else if (signbit(y) == 0)
ry = atan2(2, -y)/2;
else
ry = atan2(-2, y)/2;
--------------050105020207030702040007--
_______________________________________________
[email protected] mailing list
http://lists.freebsd.org/mailman/listinfo/freebsd-bugs
To unsubscribe, send any mail to "[email protected]"
