From: Taylor R Campbell <campbell@mumble.net>
Date: Fri, 30 Nov 2018 03:01:32 +0000 (+0000)
Subject: For +/-2i x, compute sqrt by sqrt(x) +/- i sqrt(x).
X-Git-Tag: mit-scheme-pucked-10.1.7~3^2~112
X-Git-Url: https://birchwood-abbey.net/git?a=commitdiff_plain;h=159cf0892afdc1b8be16437cb222f7b5b1ad0368;p=mit-scheme.git

For +/-2i x, compute sqrt by sqrt(x) +/- i sqrt(x).

Strike off some expected failures.
---

diff --git a/src/runtime/arith.scm b/src/runtime/arith.scm
index 4ebed69e1..c620b9068 100644
--- a/src/runtime/arith.scm
+++ b/src/runtime/arith.scm
@@ -1316,6 +1316,11 @@ USA.
 		 (rat:negative? y))
 	     (- (rat:abs x))
 	     (rat:abs x)))))
+
+(define (real:safe-negative? x)
+  (if (flonum? x)
+      (flo:safe-negative? x)
+      ((copy rat:negative?) x)))
 
 (define-syntax define-transcendental-unary
   (sc-macro-transformer
@@ -1963,8 +1968,16 @@ USA.
 
 (define (complex:sqrt z)
   (cond ((recnum? z)
-	 (complex:%make-polar (real:sqrt (complex:magnitude z))
-			      (real:/ (complex:angle z) 2)))
+	 (let ((x (rec:real-part z))
+	       (y (rec:imag-part z)))
+	   (if (real:zero? x)
+	       ;; sqrt(+/- 2i x) = sqrt(x) +/- sqrt(x)i
+	       (let ((sqrt-abs-y/2 (real:sqrt (real:/ (real:abs y) 2))))
+		 (complex:%make-rectangular
+		  sqrt-abs-y/2
+		  (real:copysign sqrt-abs-y/2 y)))
+	       (complex:%make-polar (real:sqrt (complex:magnitude z))
+                                    (real:/ (complex:angle z) 2)))))
 	((real:negative? z)
 	 (complex:%make-rectangular 0 (real:sqrt (real:negate z))))
 	(else
diff --git a/tests/runtime/test-arith.scm b/tests/runtime/test-arith.scm
index c93638458..4c6c984d0 100644
--- a/tests/runtime/test-arith.scm
+++ b/tests/runtime/test-arith.scm
@@ -667,22 +667,20 @@ USA.
      ,expect-failure)
     (,(make-rectangular 0. (* 2 flo:smallest-positive-subnormal))
      ,(make-rectangular (expt 2. (/ flo:subnormal-exponent-min 2))
-                        (expt 2. (/ flo:subnormal-exponent-min 2)))
-     ,expect-failure)
-    (+.125i .25+.25i ,expect-failure)
-    (+1/8i 1/4+1/4i ,expect-failure)
-    (+2i 1+1i ,expect-failure)
-    (+8i 2+2i ,expect-failure)
-    (+18i 3+3i ,expect-failure)
-    (+32i 4+4i ,expect-failure)
-    (+2.i 1.+1.i ,expect-failure)
-    (+8.i 2.+2.i ,expect-failure)
-    (+18.i 3.+3.i ,expect-failure)
-    (+32.i 4.+4.i ,expect-failure)
+                        (expt 2. (/ flo:subnormal-exponent-min 2))))
+    (+.125i .25+.25i)
+    (+1/8i 1/4+1/4i)
+    (+2i 1+1i)
+    (+8i 2+2i)
+    (+18i 3+3i)
+    (+32i 4+4i)
+    (+2.i 1.+1.i)
+    (+8.i 2.+2.i)
+    (+18.i 3.+3.i)
+    (+32.i 4.+4.i)
     (,(make-rectangular 0. (expt 2. flo:normal-exponent-max))
      ,(make-rectangular (expt 2. (/ (- flo:normal-exponent-max 1) 2))
-                        (expt 2. (/ (- flo:normal-exponent-max 1) 2)))
-     ,expect-failure)
+                        (expt 2. (/ (- flo:normal-exponent-max 1) 2))))
     (,(make-rectangular 0 (* 2 (expt 2 4000)))
      ,(make-rectangular (expt 2 2000) (expt 2 2000))
      ,expect-error)
@@ -693,30 +691,28 @@ USA.
     (,(make-rectangular 0. (- (* 2 flo:smallest-positive-subnormal)))
      ,(make-rectangular
        (expt 2. (/ flo:subnormal-exponent-min 2))
-       (- (expt 2. (/ flo:subnormal-exponent-min 2))))
-     ,expect-failure)
-    (-.125i .25-.25i ,expect-failure)
-    (-1/8i 1/4-1/4i ,expect-failure)
-    (-2i 1-1i ,expect-failure)
-    (-8i 2-2i ,expect-failure)
-    (-18i 3-3i ,expect-failure)
-    (-32i 4-4i ,expect-failure)
-    (-2.i 1.-1.i ,expect-failure)
-    (-8.i 2.-2.i ,expect-failure)
-    (-18.i 3.-3.i ,expect-failure)
-    (-32.i 4.-4.i ,expect-failure)
+       (- (expt 2. (/ flo:subnormal-exponent-min 2)))))
+    (-.125i .25-.25i)
+    (-1/8i 1/4-1/4i)
+    (-2i 1-1i)
+    (-8i 2-2i)
+    (-18i 3-3i)
+    (-32i 4-4i)
+    (-2.i 1.-1.i)
+    (-8.i 2.-2.i)
+    (-18.i 3.-3.i)
+    (-32.i 4.-4.i)
     (,(make-rectangular 0. (- (expt 2. flo:normal-exponent-max)))
      ,(make-rectangular (expt 2. (/ (- flo:normal-exponent-max 1) 2))
-                        (- (expt 2. (/ (- flo:normal-exponent-max 1) 2))))
-     ,expect-failure)
+                        (- (expt 2. (/ (- flo:normal-exponent-max 1) 2)))))
     (,(make-rectangular 0 (- (* 2 (expt 2 4000))))
      ,(make-rectangular (expt 2 2000) (- (expt 2 2000)))
      ,expect-error)
     ;; Handle signed zero carefully.
     (+0.i 0.+0.i)
-    (-0.i 0.-0.i ,expect-failure)
+    (-0.i 0.-0.i)
     (-0.+0.i +0.+0.i)
-    (-0.-0.i +0.-0.i ,expect-failure)
+    (-0.-0.i +0.-0.i)
     ;; Treat infinities carefully around branch cuts.
     (-inf.0 +inf.0i)
     (+inf.0 +inf.0)