[Git][cmucl/cmucl][rtoy-issue-78-unneded-code-code-in-complex-acos] Revert changes to atanh

[Raymond Toy]

Raymond Toy pushed to branch rtoy-issue-78-unneded-code-code-in-complex-acos at cmucl / cmucl


Commits:
6ed52aa0 by Raymond Toy at 2020-05-31T21:01:13-07:00
Revert changes to atanh

Still need the special case for real args to get the right branch
cuts.  Haven't yet figured out why the current code doesn't work when
the arg is a pure real.  Soemthing to do with computing y, but not
sure how to get that to work out without doing just explicitly
checking for real args.  In that case, it's just much easier to change
the arg to have the correct imaginary component to make it continuous
with the desired quadrant.

- - - - -


1 changed file:

- src/code/irrat.lisp


Changes:

=====================================
src/code/irrat.lisp
=====================================
@@ -1510,7 +1510,7 @@ Z may be any number, but the result is always a complex."
   ;; NOTE: this differs from what the CLHS says for the continuity.
   ;; Instead of the text in the CLHS, we choose to use the definition
   ;; to derive the correct values.
-  (if (and nil (realp z))
+  (if (realp z)
       (complex-atanh (complex (float z) (- (* 0 (float z)))))
       (let* ( ;; Constants
 	     (theta (/ (sqrt most-positive-double-float) 4.0d0))
@@ -1518,9 +1518,10 @@ Z may be any number, but the result is always a complex."
 	     (half-pi (/ pi 2.0d0))
 	     (rp (float (realpart z) 1.0d0))
 	     (beta (float-sign rp 1.0d0))
-	     (z* (conjugate z))
-	     (x (* beta (realpart z*)))
-	     (y (* beta (imagpart z*)))
+	     ;; x+iy = beta*conjugate(z), but being careful to produce
+	     ;; a signed-zero if z is  rational.
+	     (x (* beta (realpart z)))
+	     (y (* beta (- (float (imagpart z) 1d0))))
 	     (eta 0.0d0)
 	     (nu 0.0d0))
 	;; Shouldn't need this declare.
@@ -1562,8 +1563,7 @@ Z may be any number, but the result is always a complex."
 			       (atan (* 2.0d0 y)
 				     (- (* (- 1.0d0 x)
 					   (+ 1.0d0 x))
-					(square t1)))))
-		   (format t "eta = ~A nu ~A~%" eta nu))))
+					(square t1))))))))
 	  (coerce-to-complex-type (* beta eta)
 				  (- (* beta nu))
 				  z)))))
@@ -1715,16 +1715,13 @@ Z may be any number, but the result is always a complex."
   #+double-double
   (when (typep z '(or double-double-float (complex double-double-float)))
     (return-from complex-acos (dd-complex-acos z)))
-  (if (and nil (realp z) (> z 1))
-      ;; acos is continuous in quadrant IV in this case.
-      (complex-acos (complex z -0f0))
-      (let ((sqrt-1+z (complex-sqrt (1+z z)))
-	    (sqrt-1-z (complex-sqrt (1-z z))))
-	(with-float-traps-masked (:divide-by-zero)
-	  (complex (* 2 (atan (/ (realpart sqrt-1-z)
-				 (realpart sqrt-1+z))))
-		   (asinh (imagpart (* (conjugate sqrt-1+z)
-				       sqrt-1-z))))))))
+  (let ((sqrt-1+z (complex-sqrt (1+z z)))
+	(sqrt-1-z (complex-sqrt (1-z z))))
+    (with-float-traps-masked (:divide-by-zero)
+      (complex (* 2 (atan (/ (realpart sqrt-1-z)
+			     (realpart sqrt-1+z))))
+	       (asinh (imagpart (* (conjugate sqrt-1+z)
+				   sqrt-1-z)))))))
 
 ;; acosh(z) = 2*log(sqrt((x+1)/2) + sqrt((x-1)/2))
 ;;
@@ -1823,16 +1820,13 @@ Z may be any number, but the result is always a complex."
   #+double-double
   (when (typep z '(or double-double-float (complex double-double-float)))
     (return-from complex-asin (dd-complex-asin z)))
-  (if (and nil (realp z) (> z 1))
-      ;; asin is continuous in quadrant IV in this case.
-      (complex-asin (complex z -0f0))
-      (let ((sqrt-1-z (complex-sqrt (1-z z)))
-	    (sqrt-1+z (complex-sqrt (1+z z))))
-	(with-float-traps-masked (:divide-by-zero)
-	  (complex (atan (/ (realpart z)
-			    (realpart (* sqrt-1-z sqrt-1+z))))
-		   (asinh (imagpart (* (conjugate sqrt-1-z)
-				       sqrt-1+z))))))))
+  (let ((sqrt-1-z (complex-sqrt (1-z z)))
+	(sqrt-1+z (complex-sqrt (1+z z))))
+    (with-float-traps-masked (:divide-by-zero)
+      (complex (atan (/ (realpart z)
+			(realpart (* sqrt-1-z sqrt-1+z))))
+	       (asinh (imagpart (* (conjugate sqrt-1-z)
+				   sqrt-1+z)))))))
 
 (defun complex-asinh (z)
   "Compute asinh z = log(z + sqrt(1 + z*z))



View it on GitLab: https://gitlab.common-lisp.net/cmucl/cmucl/-/commit/6ed52aa06cd1b22f33b9417bf57585cc3235b6b4

-- 
View it on GitLab: https://gitlab.common-lisp.net/cmucl/cmucl/-/commit/6ed52aa06cd1b22f33b9417bf57585cc3235b6b4
You're receiving this email because of your account on gitlab.common-lisp.net.


-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://mailman.common-lisp.net/pipermail/cmucl-cvs/attachments/20200601/97c23fd0/attachment-0001.htm>