Module Name:    src
Committed By:   maya
Date:           Sat Dec 31 20:01:15 UTC 2016

Modified Files:
        src/lib/libm/complex: csqrt.c csqrtf.c

Log Message:
csqrt has a branch cut on the negative real axis, and this requires
delicacy in order to maintain continuity around it.

we have an initial case to deal with a fairly common case: getting
a real number. Avoid dealing with the branch cut in this case by
checking if the real part is negative.

later, -0.0 < 0 is not met, so instead, test for a negative number
using signbit, so negative zero is also treated as a negative number.

Fixes last part of PR lib/51427: libm issues triggered by py-numpy

ok riastradh


To generate a diff of this commit:
cvs rdiff -u -r1.2 -r1.3 src/lib/libm/complex/csqrt.c
cvs rdiff -u -r1.1 -r1.2 src/lib/libm/complex/csqrtf.c

Please note that diffs are not public domain; they are subject to the
copyright notices on the relevant files.

Modified files:

Index: src/lib/libm/complex/csqrt.c
diff -u src/lib/libm/complex/csqrt.c:1.2 src/lib/libm/complex/csqrt.c:1.3
--- src/lib/libm/complex/csqrt.c:1.2	Sat Dec 31 15:33:03 2016
+++ src/lib/libm/complex/csqrt.c	Sat Dec 31 20:01:15 2016
@@ -1,4 +1,4 @@
-/* $NetBSD: csqrt.c,v 1.2 2016/12/31 15:33:03 maya Exp $ */
+/* $NetBSD: csqrt.c,v 1.3 2016/12/31 20:01:15 maya Exp $ */
 
 /*-
  * Copyright (c) 2007 The NetBSD Foundation, Inc.
@@ -41,7 +41,8 @@ csqrt(double complex z)
 	x = creal (z);
 	y = cimag (z);
 
-	if (y == 0.0) {
+	/* Input is a real number that isn't on the branch cut */
+	if ((y == 0.0) && !signbit(y)) {
 		if (x == 0.0) {
 			w = 0.0 + y * I;
 		} else {
@@ -92,7 +93,7 @@ csqrt(double complex z)
 		t = scale * fabs((0.5 * y) / r);
 		r *= scale;
 	}
-	if (y < 0)
+	if (signbit(y))
 		w = t - r * I;
 	else
 		w = t + r * I;

Index: src/lib/libm/complex/csqrtf.c
diff -u src/lib/libm/complex/csqrtf.c:1.1 src/lib/libm/complex/csqrtf.c:1.2
--- src/lib/libm/complex/csqrtf.c:1.1	Mon Aug 20 16:01:37 2007
+++ src/lib/libm/complex/csqrtf.c	Sat Dec 31 20:01:15 2016
@@ -1,4 +1,4 @@
-/* $NetBSD: csqrtf.c,v 1.1 2007/08/20 16:01:37 drochner Exp $ */
+/* $NetBSD: csqrtf.c,v 1.2 2016/12/31 20:01:15 maya Exp $ */
 
 /*-
  * Copyright (c) 2007 The NetBSD Foundation, Inc.
@@ -41,7 +41,8 @@ csqrtf(float complex z)
 	x = crealf (z);
 	y = cimagf (z);
 
-	if (y == 0.0f) {
+	/* Input is a real number that isn't on the branch cut */
+	if ((y == 0.0f) && !signbit(y)) {
 		if (x < 0.0f) {
 			w = 0.0f + sqrtf(-x) * I;
 			return w;
@@ -91,7 +92,7 @@ csqrtf(float complex z)
 		r *= scale;
 	}
 
-	if (y < 0)
+	if (signbit(y))
 		w = t - r * I;
 	else
 		w = t + r * I;

Reply via email to