unofficial mirror of guile-devel@gnu.org 
 help / color / mirror / Atom feed
* [PATCH] Complex numbers with inexact zero imaginary part, etc
@ 2011-02-02 11:25 Mark H Weaver
  2011-02-02 20:29 ` Andy Wingo
                   ` (2 more replies)
  0 siblings, 3 replies; 9+ messages in thread
From: Mark H Weaver @ 2011-02-02 11:25 UTC (permalink / raw)
  To: guile-devel

[-- Attachment #1: Type: text/plain, Size: 615 bytes --]

Here's another batch of numerics fixes and changes.  Most notably, the
final patch allows Guile to represent non-real complex numbers with
inexact zero imaginary part.  The first two patches fix bugs and improve
handling of signed zeroes.

Note that the patches are meant to be applied after my earlier patch
that makes trig functions return exact results in some cases, although
it _might_ work fine to apply these first.

There are two more patches coming soon: one to improve expt, and one to
add many test cases to numbers.test.  Hopefully these can make it into
the imminent prerelease.

    Best,
     Mark



[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #2: Fix `min' and `max' handling of NaNs, infinities, and signed zeroes --]
[-- Type: text/x-diff, Size: 19281 bytes --]

From 5e174213416b97e9446dccac70fd9689106a6fd6 Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Wed, 2 Feb 2011 01:02:49 -0500
Subject: [PATCH] Fix `min' and `max' handling of NaNs, infinities, and signed zeroes

* libguile/numbers.c (scm_min, scm_max): Properly order the real
  infinities and NaNs, per R6RS, and also take care to handle signed
  zeroes properly.  Note that this ordering is different than that of
  `<', `>', `<=', and `>=', which return #f if any argument is a real
  NaN, and consider the real zeroes to be equal.  The relevant real
  infinity (-inf.0 for min, +inf.0 for max) beats everything, including
  NaNs, and NaNs beat everything else.  Previously these were handled
  improperly in some cases, e.g.:
  (min 1/2 +nan.0) now returns +nan.0 (previously returned 0.5),
  (max 1/2 +nan.0) now returns +nan.0 (previously returned 0.5),
  (min -inf.0 +nan.0) now returns -inf.0 (previously returned +nan.0),
  (max +inf.0 +nan.0) now returns +inf.0 (previously returned +nan.0),
  (min -0.0  0.0) now returns -0.0 (previously returned  0.0),
  (max  0.0 -0.0) now returns  0.0 (previously returned -0.0),
  (max  0   -0.0) now returns  0.0 (previously returned -0.0),
  (max -0.0  0  ) now returns  0.0 (previously returned -0.0).

* test-suite/tests/numbers.test (min, max): Add many more test cases
  relating to NaNs, infinities, and signed zeroes.  Change most existing
  test cases to use `eqv?' instead of `=', in order to check exactness.
---
 libguile/numbers.c            |   93 ++++++++++++++----
 test-suite/tests/numbers.test |  226 ++++++++++++++++++++++++++++++-----------
 2 files changed, 241 insertions(+), 78 deletions(-)

diff --git a/libguile/numbers.c b/libguile/numbers.c
index df95c32..090fb75 100644
--- a/libguile/numbers.c
+++ b/libguile/numbers.c
@@ -498,6 +498,14 @@ scm_i_fraction2double (SCM z)
 					 SCM_FRACTION_DENOMINATOR (z)));
 }
 
+static int
+double_is_non_negative_zero (double x)
+{
+  static double zero = 0.0;
+
+  return !memcmp (&x, &zero, sizeof(double));
+}
+
 SCM_PRIMITIVE_GENERIC (scm_exact_p, "exact?", 1, 0, 0, 
 		       (SCM x),
 	    "Return @code{#t} if @var{x} is an exact number, @code{#f}\n"
@@ -5148,9 +5156,19 @@ scm_max (SCM x, SCM y)
 	}
       else if (SCM_REALP (y))
 	{
-	  double z = xx;
-	  /* if y==NaN then ">" is false and we return NaN */
-	  return (z > SCM_REAL_VALUE (y)) ? scm_from_double (z) : y;
+	  double xxd = xx;
+	  double yyd = SCM_REAL_VALUE (y);
+
+	  if (xxd > yyd)
+	    return scm_from_double (xxd);
+	  /* If y is a NaN, then "==" is false and we return the NaN */
+	  else if (SCM_LIKELY (!(xxd == yyd)))
+	    return y;
+	  /* Handle signed zeroes properly */
+	  else if (xx == 0)
+	    return flo0;
+	  else
+	    return y;
 	}
       else if (SCM_FRACTIONP (y))
 	{
@@ -5194,9 +5212,20 @@ scm_max (SCM x, SCM y)
     {
       if (SCM_I_INUMP (y))
 	{
-	  double z = SCM_I_INUM (y);
-	  /* if x==NaN then "<" is false and we return NaN */
-	  return (SCM_REAL_VALUE (x) < z) ? scm_from_double (z) : x;
+	  scm_t_inum yy = SCM_I_INUM (y);
+	  double xxd = SCM_REAL_VALUE (x);
+	  double yyd = yy;
+
+	  if (yyd > xxd)
+	    return scm_from_double (yyd);
+	  /* If x is a NaN, then "==" is false and we return the NaN */
+	  else if (SCM_LIKELY (!(xxd == yyd)))
+	    return x;
+	  /* Handle signed zeroes properly */
+	  else if (yy == 0)
+	    return flo0;
+	  else
+	    return x;
 	}
       else if (SCM_BIGP (y))
 	{
@@ -5205,12 +5234,25 @@ scm_max (SCM x, SCM y)
 	}
       else if (SCM_REALP (y))
 	{
-	  /* if x==NaN then our explicit check means we return NaN
-	     if y==NaN then ">" is false and we return NaN
-	     calling isnan is unavoidable, since it's the only way to know
-	     which of x or y causes any compares to be false */
 	  double xx = SCM_REAL_VALUE (x);
-	  return (isnan (xx) || xx > SCM_REAL_VALUE (y)) ? x : y;
+	  double yy = SCM_REAL_VALUE (y);
+
+	  /* For purposes of max: +inf.0 > nan > everything else, per R6RS */
+	  if (xx > yy)
+	    return x;
+	  else if (SCM_LIKELY (xx < yy))
+	    return y;
+	  /* If neither (xx > yy) nor (xx < yy), then
+	     either they're equal or one is a NaN */
+	  else if (SCM_UNLIKELY (isnan (xx)))
+	    return (isinf (yy) == 1) ? y : x;
+	  else if (SCM_UNLIKELY (isnan (yy)))
+	    return (isinf (xx) == 1) ? x : y;
+	  /* xx == yy, but handle signed zeroes properly */
+	  else if (double_is_non_negative_zero (yy))
+	    return y;
+	  else
+	    return x;
 	}
       else if (SCM_FRACTIONP (y))
 	{
@@ -5234,7 +5276,8 @@ scm_max (SCM x, SCM y)
       else if (SCM_REALP (y))
 	{
 	  double xx = scm_i_fraction2double (x);
-	  return (xx < SCM_REAL_VALUE (y)) ? y : scm_from_double (xx);
+	  /* if y==NaN then ">" is false, so we return the NaN y */
+	  return (xx > SCM_REAL_VALUE (y)) ? scm_from_double (xx) : y;
 	}
       else if (SCM_FRACTIONP (y))
 	{
@@ -5351,12 +5394,25 @@ scm_min (SCM x, SCM y)
 	}
       else if (SCM_REALP (y))
 	{
-	  /* if x==NaN then our explicit check means we return NaN
-	     if y==NaN then "<" is false and we return NaN
-	     calling isnan is unavoidable, since it's the only way to know
-	     which of x or y causes any compares to be false */
 	  double xx = SCM_REAL_VALUE (x);
-	  return (isnan (xx) || xx < SCM_REAL_VALUE (y)) ? x : y;
+	  double yy = SCM_REAL_VALUE (y);
+
+	  /* For purposes of min: -inf.0 < nan < everything else, per R6RS */
+	  if (xx < yy)
+	    return x;
+	  else if (SCM_LIKELY (xx > yy))
+	    return y;
+	  /* If neither (xx < yy) nor (xx > yy), then
+	     either they're equal or one is a NaN */
+	  else if (SCM_UNLIKELY (isnan (xx)))
+	    return (isinf (yy) == -1) ? y : x;
+	  else if (SCM_UNLIKELY (isnan (yy)))
+	    return (isinf (xx) == -1) ? x : y;
+	  /* xx == yy, but handle signed zeroes properly */
+	  else if (double_is_non_negative_zero (xx))
+	    return y;
+	  else
+	    return x;
 	}
       else if (SCM_FRACTIONP (y))
 	{
@@ -5380,7 +5436,8 @@ scm_min (SCM x, SCM y)
       else if (SCM_REALP (y))
 	{
 	  double xx = scm_i_fraction2double (x);
-	  return (SCM_REAL_VALUE (y) < xx) ? y : scm_from_double (xx);
+	  /* if y==NaN then "<" is false, so we return the NaN y */
+	  return (xx < SCM_REAL_VALUE (y)) ? scm_from_double (xx) : y;
 	}
       else if (SCM_FRACTIONP (y))
 	{
diff --git a/test-suite/tests/numbers.test b/test-suite/tests/numbers.test
index 9c01fa1..28db652 100644
--- a/test-suite/tests/numbers.test
+++ b/test-suite/tests/numbers.test
@@ -2471,26 +2471,79 @@
         (big*5 (* fixnum-max 5)))
 
     (with-test-prefix "inum / frac"
-      (pass-if (= 3 (max 3 5/2)))
-      (pass-if (= 5/2 (max 2 5/2))))
+      (pass-if (eqv? 3 (max 3 5/2)))
+      (pass-if (eqv? 5/2 (max 2 5/2))))
 
     (with-test-prefix "frac / inum"
-      (pass-if (= 3 (max 5/2 3)))
-      (pass-if (= 5/2 (max 5/2 2))))
-
-    (with-test-prefix "inum / real"
-      (pass-if (real-nan? (max 123 +nan.0))))
-
-    (with-test-prefix "real / inum"
-      (pass-if (real-nan? (max +nan.0 123))))
+      (pass-if (eqv? 3 (max 5/2 3)))
+      (pass-if (eqv? 5/2 (max 5/2 2))))
+
+    (with-test-prefix "infinities and NaNs"
+      ;; +inf.0 beats everything else, including NaNs
+      (pass-if (eqv?  +inf.0   (max   +inf.0    123   )))
+      (pass-if (eqv?  +inf.0   (max    123     +inf.0 )))
+      (pass-if (eqv?  +inf.0   (max   +inf.0   -123.3 )))
+      (pass-if (eqv?  +inf.0   (max   -123.3   +inf.0 )))
+      (pass-if (eqv?  +inf.0   (max   +inf.0   -7/2   )))
+      (pass-if (eqv?  +inf.0   (max   -7/2     +inf.0 )))
+      (pass-if (eqv?  +inf.0   (max   +inf.0  -1e20   )))
+      (pass-if (eqv?  +inf.0   (max  -1e20     +inf.0 )))
+      (pass-if (eqv?  +inf.0   (max   +inf.0 (- big*2))))
+      (pass-if (eqv?  +inf.0   (max (- big*2)  +inf.0 )))
+      (pass-if (eqv?  +inf.0   (max   +inf.0   +inf.0 )))
+      (pass-if (eqv?  +inf.0   (max   +inf.0   +inf.0 )))
+      (pass-if (eqv?  +inf.0   (max   +inf.0   +nan.0 )))
+      (pass-if (eqv?  +inf.0   (max   +nan.0   +inf.0 )))
+      (pass-if (eqv?  +inf.0   (max   +inf.0   +inf.0 )))
+
+      ;; NaNs beat everything except +inf.0
+      (pass-if (real-nan?      (max   +nan.0    123   )))
+      (pass-if (real-nan?      (max    123     +nan.0 )))
+      (pass-if (real-nan?      (max   +nan.0    123.3 )))
+      (pass-if (real-nan?      (max    123.3   +nan.0 )))
+      (pass-if (real-nan?      (max   +nan.0   -7/2   )))
+      (pass-if (real-nan?      (max   -7/2     +nan.0 )))
+      (pass-if (real-nan?      (max   +nan.0  -1e20   )))
+      (pass-if (real-nan?      (max  -1e20     +nan.0 )))
+      (pass-if (real-nan?      (max   +nan.0 (- big*2))))
+      (pass-if (real-nan?      (max (- big*2)  +nan.0 )))
+      (pass-if (real-nan?      (max   +nan.0   -inf.0 )))
+      (pass-if (real-nan?      (max   -inf.0   +nan.0 )))
+      (pass-if (real-nan?      (max   +nan.0   +nan.0 )))
+
+      ;; -inf.0 always loses, except against itself
+      (pass-if (eqv?   -inf.0  (max   -inf.0   -inf.0 )))
+      (pass-if (eqv?   -123.0  (max   -inf.0   -123   )))
+      (pass-if (eqv?   -123.0  (max   -123     -inf.0 )))
+      (pass-if (eqv?   -123.3  (max   -inf.0   -123.3 )))
+      (pass-if (eqv?   -123.3  (max   -123.3   -inf.0 )))
+      (pass-if (eqv?     -3.5  (max   -inf.0   -7/2   )))
+      (pass-if (eqv?     -3.5  (max   -7/2     -inf.0 )))
+      (pass-if (eqv?  -1.0e20  (max   -inf.0  -1e20   )))
+      (pass-if (eqv?  -1.0e20  (max  -1e20     -inf.0 )))
+      (pass-if (eqv? (exact->inexact (- big*2))
+                     (max   -inf.0 (- big*2))))
+      (pass-if (eqv? (exact->inexact (- big*2))
+                     (max (- big*2)  -inf.0 ))))
+
+    (with-test-prefix "signed zeroes"
+      (pass-if (eqv?  0.0  (max  0.0  0.0)))
+      (pass-if (eqv?  0.0  (max  0.0 -0.0)))
+      (pass-if (eqv?  0.0  (max -0.0  0.0)))
+      (pass-if (eqv? -0.0  (max -0.0 -0.0)))
+      (pass-if (eqv?  0.0  (max -0.0  0  )))
+      (pass-if (eqv?  0.0  (max  0.0  0  )))
+      (pass-if (eqv?  0.0  (max  0   -0.0)))
+      (pass-if (eqv?  0.0  (max  0    0.0)))
+      (pass-if (eqv?  0    (min  0    0  ))))
 
     (with-test-prefix "big / frac"
-      (pass-if (= big*2 (max big*2 5/2)))
-      (pass-if (= 5/2 (max (- big*2) 5/2))))
+      (pass-if (eqv? big*2 (max big*2 5/2)))
+      (pass-if (eqv? 5/2 (max (- big*2) 5/2))))
 
     (with-test-prefix "frac / big"
-      (pass-if (= big*2 (max 5/2 big*2)))
-      (pass-if (= 5/2 (max 5/2 (- big*2)))))
+      (pass-if (eqv? big*2 (max 5/2 big*2)))
+      (pass-if (eqv? 5/2 (max 5/2 (- big*2)))))
 
     (with-test-prefix "big / real"
       (pass-if (real-nan? (max big*5 +nan.0)))
@@ -2507,29 +2560,29 @@
       (pass-if (eqv? 1.0                     (max 1.0 (- big*5)))))
 
     (with-test-prefix "frac / frac"
-      (pass-if (= 2/3 (max 1/2 2/3)))
-      (pass-if (= 2/3 (max 2/3 1/2)))
-      (pass-if (= -1/2 (max -1/2 -2/3)))
-      (pass-if (= -1/2 (max -2/3 -1/2))))
+      (pass-if (eqv? 2/3 (max 1/2 2/3)))
+      (pass-if (eqv? 2/3 (max 2/3 1/2)))
+      (pass-if (eqv? -1/2 (max -1/2 -2/3)))
+      (pass-if (eqv? -1/2 (max -2/3 -1/2))))
 
     (with-test-prefix "real / real"
       (pass-if (real-nan? (max 123.0 +nan.0)))
       (pass-if (real-nan? (max +nan.0 123.0)))
       (pass-if (real-nan? (max +nan.0 +nan.0)))
-      (pass-if (= 456.0 (max 123.0 456.0)))
-      (pass-if (= 456.0 (max 456.0 123.0)))))
+      (pass-if (eqv? 456.0 (max 123.0 456.0)))
+      (pass-if (eqv? 456.0 (max 456.0 123.0)))))
 
   ;; in gmp prior to 4.2, mpz_cmp_d ended up treating Inf as 2^1024, make
   ;; sure we've avoided that
   (for-each (lambda (b)
 	      (pass-if (list b +inf.0)
-		(= +inf.0 (max b +inf.0)))
+		(eqv? +inf.0 (max b +inf.0)))
 	      (pass-if (list +inf.0 b)
-		(= +inf.0 (max b +inf.0)))
+		(eqv? +inf.0 (max b +inf.0)))
 	      (pass-if (list b -inf.0)
-		(= (exact->inexact b) (max b -inf.0)))
+		(eqv? (exact->inexact b) (max b -inf.0)))
 	      (pass-if (list -inf.0 b)
-		(= (exact->inexact b) (max b -inf.0))))
+		(eqv? (exact->inexact b) (max b -inf.0))))
 	    (list (1- (ash 1 1024))
 		  (ash 1 1024)
 		  (1+ (ash 1 1024))
@@ -2579,43 +2632,96 @@
         (big*5 (* fixnum-max 5)))
 
     (pass-if (documented? min))
-    (pass-if (= 1 (min 7 3 1 5)))
-    (pass-if (= 1 (min 1 7 3 5)))
-    (pass-if (= 1 (min 7 3 5 1)))
-    (pass-if (= -7 (min 2 3 4 -2 5 -7 1 -1 4 2)))
-    (pass-if (= -7 (min -7 2 3 4 -2 5 1 -1 4 2)))
-    (pass-if (= -7 (min 2 3 4 -2 5 1 -1 4 2 -7)))
-    (pass-if (= big*2 (min big*3 big*5 big*2 big*4)))
-    (pass-if (= big*2 (min big*2 big*3 big*5 big*4)))
-    (pass-if (= big*2 (min big*3 big*5 big*4 big*2)))
+    (pass-if (eqv? 1 (min 7 3 1 5)))
+    (pass-if (eqv? 1 (min 1 7 3 5)))
+    (pass-if (eqv? 1 (min 7 3 5 1)))
+    (pass-if (eqv? -7 (min 2 3 4 -2 5 -7 1 -1 4 2)))
+    (pass-if (eqv? -7 (min -7 2 3 4 -2 5 1 -1 4 2)))
+    (pass-if (eqv? -7 (min 2 3 4 -2 5 1 -1 4 2 -7)))
+    (pass-if (eqv? big*2 (min big*3 big*5 big*2 big*4)))
+    (pass-if (eqv? big*2 (min big*2 big*3 big*5 big*4)))
+    (pass-if (eqv? big*2 (min big*3 big*5 big*4 big*2)))
     (pass-if
-        (= (- fixnum-min 1) (min 2 4 (- fixnum-min 1) 3 (* 2 fixnum-max))))
+        (eqv? (- fixnum-min 1) (min 2 4 (- fixnum-min 1) 3 (* 2 fixnum-max))))
     (pass-if
-        (= (- fixnum-min 1) (min (- fixnum-min 1) 2 4 3 (* 2 fixnum-max))))
+        (eqv? (- fixnum-min 1) (min (- fixnum-min 1) 2 4 3 (* 2 fixnum-max))))
     (pass-if
-        (= (- fixnum-min 1) (min 2 4 3 (* 2 fixnum-max) (- fixnum-min 1))))
+        (eqv? (- fixnum-min 1) (min 2 4 3 (* 2 fixnum-max) (- fixnum-min 1))))
 
     (with-test-prefix "inum / frac"
-      (pass-if (= 5/2 (min 3 5/2)))
-      (pass-if (= 2   (min 2 5/2))))
+      (pass-if (eqv? 5/2 (min 3 5/2)))
+      (pass-if (eqv? 2   (min 2 5/2))))
 
     (with-test-prefix "frac / inum"
-      (pass-if (= 5/2 (min 5/2 3)))
-      (pass-if (= 2   (min 5/2 2))))
-
-    (with-test-prefix "inum / real"
-      (pass-if (real-nan? (min 123 +nan.0))))
-
-    (with-test-prefix "real / inum"
-      (pass-if (real-nan? (min +nan.0 123))))
+      (pass-if (eqv? 5/2 (min 5/2 3)))
+      (pass-if (eqv? 2   (min 5/2 2))))
+
+    (with-test-prefix "infinities and NaNs"
+      ;; -inf.0 beats everything else, including NaNs
+      (pass-if (eqv?  -inf.0   (min   -inf.0    123   )))
+      (pass-if (eqv?  -inf.0   (min    123     -inf.0 )))
+      (pass-if (eqv?  -inf.0   (min   -inf.0   -123.3 )))
+      (pass-if (eqv?  -inf.0   (min   -123.3   -inf.0 )))
+      (pass-if (eqv?  -inf.0   (min   -inf.0   -7/2   )))
+      (pass-if (eqv?  -inf.0   (min   -7/2     -inf.0 )))
+      (pass-if (eqv?  -inf.0   (min   -inf.0  -1e20   )))
+      (pass-if (eqv?  -inf.0   (min  -1e20     -inf.0 )))
+      (pass-if (eqv?  -inf.0   (min   -inf.0 (- big*2))))
+      (pass-if (eqv?  -inf.0   (min (- big*2)  -inf.0 )))
+      (pass-if (eqv?  -inf.0   (min   -inf.0   +inf.0 )))
+      (pass-if (eqv?  -inf.0   (min   +inf.0   -inf.0 )))
+      (pass-if (eqv?  -inf.0   (min   -inf.0   +nan.0 )))
+      (pass-if (eqv?  -inf.0   (min   +nan.0   -inf.0 )))
+      (pass-if (eqv?  -inf.0   (min   -inf.0   -inf.0 )))
+
+      ;; NaNs beat everything except -inf.0
+      (pass-if (real-nan?      (min   +nan.0    123   )))
+      (pass-if (real-nan?      (min    123     +nan.0 )))
+      (pass-if (real-nan?      (min   +nan.0    123.3 )))
+      (pass-if (real-nan?      (min    123.3   +nan.0 )))
+      (pass-if (real-nan?      (min   +nan.0   -7/2   )))
+      (pass-if (real-nan?      (min   -7/2     +nan.0 )))
+      (pass-if (real-nan?      (min   +nan.0  -1e20   )))
+      (pass-if (real-nan?      (min  -1e20     +nan.0 )))
+      (pass-if (real-nan?      (min   +nan.0 (- big*2))))
+      (pass-if (real-nan?      (min (- big*2)  +nan.0 )))
+      (pass-if (real-nan?      (min   +nan.0   +inf.0 )))
+      (pass-if (real-nan?      (min   +inf.0   +nan.0 )))
+      (pass-if (real-nan?      (min   +nan.0   +nan.0 )))
+
+      ;; +inf.0 always loses, except against itself
+      (pass-if (eqv?   +inf.0  (min   +inf.0   +inf.0 )))
+      (pass-if (eqv?   -123.0  (min   +inf.0   -123   )))
+      (pass-if (eqv?   -123.0  (min   -123     +inf.0 )))
+      (pass-if (eqv?   -123.3  (min   +inf.0   -123.3 )))
+      (pass-if (eqv?   -123.3  (min   -123.3   +inf.0 )))
+      (pass-if (eqv?     -3.5  (min   +inf.0   -7/2   )))
+      (pass-if (eqv?     -3.5  (min   -7/2     +inf.0 )))
+      (pass-if (eqv?  -1.0e20  (min   +inf.0  -1e20   )))
+      (pass-if (eqv?  -1.0e20  (min  -1e20     +inf.0 )))
+      (pass-if (eqv? (exact->inexact (- big*2))
+                     (min   +inf.0 (- big*2))))
+      (pass-if (eqv? (exact->inexact (- big*2))
+                     (min (- big*2)  +inf.0 ))))
+
+    (with-test-prefix "signed zeroes"
+      (pass-if (eqv?  0.0  (min  0.0  0.0)))
+      (pass-if (eqv? -0.0  (min  0.0 -0.0)))
+      (pass-if (eqv? -0.0  (min -0.0  0.0)))
+      (pass-if (eqv? -0.0  (min -0.0 -0.0)))
+      (pass-if (eqv? -0.0  (min -0.0  0  )))
+      (pass-if (eqv?  0.0  (min  0.0  0  )))
+      (pass-if (eqv? -0.0  (min  0   -0.0)))
+      (pass-if (eqv?  0.0  (min  0    0.0)))
+      (pass-if (eqv?  0    (min  0    0  ))))
 
     (with-test-prefix "big / frac"
-      (pass-if (= 5/2       (min big*2 5/2)))
-      (pass-if (= (- big*2) (min (- big*2) 5/2))))
+      (pass-if (eqv? 5/2       (min big*2 5/2)))
+      (pass-if (eqv? (- big*2) (min (- big*2) 5/2))))
 
     (with-test-prefix "frac / big"
-      (pass-if (= 5/2       (min 5/2 big*2)))
-      (pass-if (= (- big*2) (min 5/2 (- big*2)))))
+      (pass-if (eqv? 5/2       (min 5/2 big*2)))
+      (pass-if (eqv? (- big*2) (min 5/2 (- big*2)))))
 
     (with-test-prefix "big / real"
       (pass-if (real-nan? (min big*5 +nan.0)))
@@ -2632,30 +2738,30 @@
       (pass-if (eqv? (exact->inexact (- big*5))  (min 1.0 (- big*5)))))
 
     (with-test-prefix "frac / frac"
-      (pass-if (= 1/2 (min 1/2 2/3)))
-      (pass-if (= 1/2 (min 2/3 1/2)))
-      (pass-if (= -2/3 (min -1/2 -2/3)))
-      (pass-if (= -2/3 (min -2/3 -1/2))))
+      (pass-if (eqv? 1/2 (min 1/2 2/3)))
+      (pass-if (eqv? 1/2 (min 2/3 1/2)))
+      (pass-if (eqv? -2/3 (min -1/2 -2/3)))
+      (pass-if (eqv? -2/3 (min -2/3 -1/2))))
 
     (with-test-prefix "real / real"
       (pass-if (real-nan? (min 123.0 +nan.0)))
       (pass-if (real-nan? (min +nan.0 123.0)))
       (pass-if (real-nan? (min +nan.0 +nan.0)))
-      (pass-if (= 123.0 (min 123.0 456.0)))
-      (pass-if (= 123.0 (min 456.0 123.0)))))
+      (pass-if (eqv? 123.0 (min 123.0 456.0)))
+      (pass-if (eqv? 123.0 (min 456.0 123.0)))))
 
 
   ;; in gmp prior to 4.2, mpz_cmp_d ended up treating Inf as 2^1024, make
   ;; sure we've avoided that
   (for-each (lambda (b)
 	      (pass-if (list b +inf.0)
-		(= (exact->inexact b) (min b +inf.0)))
+		(eqv? (exact->inexact b) (min b +inf.0)))
 	      (pass-if (list +inf.0 b)
-		(= (exact->inexact b) (min b +inf.0)))
+		(eqv? (exact->inexact b) (min b +inf.0)))
 	      (pass-if (list b -inf.0)
-		(= -inf.0 (min b -inf.0)))
+		(eqv? -inf.0 (min b -inf.0)))
 	      (pass-if (list -inf.0 b)
-		(= -inf.0 (min b -inf.0))))
+		(eqv? -inf.0 (min b -inf.0))))
 	    (list (1- (ash 1 1024))
 		  (ash 1 1024)
 		  (1+ (ash 1 1024))
-- 
1.5.6.5


[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #3: Improve handling of signed zeroes --]
[-- Type: text/x-diff, Size: 6658 bytes --]

From 0515f466960fd80f58a5a90b8bfdb530706b982a Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Wed, 2 Feb 2011 03:14:13 -0500
Subject: [PATCH] Improve handling of signed zeroes

* libguile/numbers.c (scm_abs): (abs -0.0) now returns 0.0.  Previously
  it returned -0.0.  Also move the REALP case above the BIGP case,
  and consider it SCM_LIKELY to be REALP if not INUMP.
  (scm_difference): (- 0 0.0) now returns -0.0.  Previously it returned
  0.0.  Also make sure that (- 0 0.0+0.0i) will return -0.0-0.0i.

* test-suite/tests/numbers.test (abs, -): Add test cases, and change
  some tests to use `eqv?' instead of `=', in order to test exactness
  and distinguish signed zeroes.
---
 libguile/numbers.c            |   49 +++++++++++++++++++------
 test-suite/tests/numbers.test |   77 +++++++++++++++++++++++++++++++++++------
 2 files changed, 103 insertions(+), 23 deletions(-)

diff --git a/libguile/numbers.c b/libguile/numbers.c
index 090fb75..9a7848a 100644
--- a/libguile/numbers.c
+++ b/libguile/numbers.c
@@ -745,6 +745,18 @@ SCM_PRIMITIVE_GENERIC (scm_abs, "abs", 1, 0, 0,
       else
 	return scm_i_inum2big (-xx);
     }
+  else if (SCM_LIKELY (SCM_REALP (x)))
+    {
+      double xx = SCM_REAL_VALUE (x);
+      /* If x is a NaN then xx<0 is false so we return x unchanged */
+      if (xx < 0.0)
+        return scm_from_double (-xx);
+      /* Handle signed zeroes properly */
+      else if (SCM_UNLIKELY (xx == 0.0))
+	return flo0;
+      else
+        return x;
+    }
   else if (SCM_BIGP (x))
     {
       const int sgn = mpz_sgn (SCM_I_BIG_MPZ (x));
@@ -753,15 +765,6 @@ SCM_PRIMITIVE_GENERIC (scm_abs, "abs", 1, 0, 0,
       else
 	return x;
     }
-  else if (SCM_REALP (x))
-    {
-      /* note that if x is a NaN then xx<0 is false so we return x unchanged */
-      double xx = SCM_REAL_VALUE (x);
-      if (xx < 0.0)
-        return scm_from_double (-xx);
-      else
-        return x;
-    }
   else if (SCM_FRACTIONP (x))
     {
       if (scm_is_false (scm_negative_p (SCM_FRACTION_NUMERATOR (x))))
@@ -5758,13 +5761,35 @@ scm_difference (SCM x, SCM y)
       else if (SCM_REALP (y))
 	{
 	  scm_t_inum xx = SCM_I_INUM (x);
-	  return scm_from_double (xx - SCM_REAL_VALUE (y));
+
+	  /*
+	   * We need to handle x == exact 0
+	   * specially because R6RS states that:
+	   *   (- 0.0)     ==> -0.0  and
+	   *   (- 0.0 0.0) ==>  0.0
+	   * and the scheme compiler changes
+	   *   (- 0.0) into (- 0 0.0)
+	   * So we need to treat (- 0 0.0) like (- 0.0).
+	   * At the C level, (-x) is different than (0.0 - x).
+	   * (0.0 - 0.0) ==> 0.0, but (- 0.0) ==> -0.0.
+	   */
+	  if (xx == 0)
+	    return scm_from_double (- SCM_REAL_VALUE (y));
+	  else
+	    return scm_from_double (xx - SCM_REAL_VALUE (y));
 	}
       else if (SCM_COMPLEXP (y))
 	{
 	  scm_t_inum xx = SCM_I_INUM (x);
-	  return scm_c_make_rectangular (xx - SCM_COMPLEX_REAL (y),
-				   - SCM_COMPLEX_IMAG (y));
+
+	  /* We need to handle x == exact 0 specially.
+	     See the comment above (for SCM_REALP (y)) */
+	  if (xx == 0)
+	    return scm_c_make_rectangular (- SCM_COMPLEX_REAL (y),
+					   - SCM_COMPLEX_IMAG (y));
+	  else
+	    return scm_c_make_rectangular (xx - SCM_COMPLEX_REAL (y),
+					      - SCM_COMPLEX_IMAG (y));
 	}
       else if (SCM_FRACTIONP (y))
 	/* a - b/c = (ac - b) / c */
diff --git a/test-suite/tests/numbers.test b/test-suite/tests/numbers.test
index 28db652..5a8b31b 100644
--- a/test-suite/tests/numbers.test
+++ b/test-suite/tests/numbers.test
@@ -423,17 +423,23 @@
 
 (with-test-prefix "abs"
   (pass-if (documented? abs))
-  (pass-if (zero? (abs 0)))
-  (pass-if (= 1 (abs 1)))
-  (pass-if (= 1 (abs -1)))
-  (pass-if (= (+ fixnum-max 1) (abs (+ fixnum-max 1))))
-  (pass-if (= (+ (- fixnum-min) 1) (abs (- fixnum-min 1))))  
-  (pass-if (= 0.0 (abs 0.0)))
-  (pass-if (= 1.0 (abs 1.0)))
-  (pass-if (= 1.0 (abs -1.0)))
-  (pass-if (real-nan? (abs +nan.0)))
-  (pass-if (= +inf.0 (abs +inf.0)))
-  (pass-if (= +inf.0 (abs -inf.0))))
+  (pass-if (eqv? 0 (abs 0)))
+  (pass-if (eqv? 1 (abs 1)))
+  (pass-if (eqv? 1 (abs -1)))
+
+  (with-test-prefix "double-negation of fixnum-min"
+    (pass-if (eqv?   fixnum-min (- (abs fixnum-min)))))
+
+  (pass-if (eqv? (+    fixnum-max  1) (abs (+ fixnum-max 1))))
+  (pass-if (eqv? (+ (- fixnum-min) 1) (abs (- fixnum-min 1))))
+
+  (pass-if (eqv? 0.0    (abs  0.0)))
+  (pass-if (eqv? 0.0    (abs -0.0)))
+  (pass-if (eqv? 1.0    (abs  1.0)))
+  (pass-if (eqv? 1.0    (abs -1.0)))
+  (pass-if (real-nan?   (abs +nan.0)))
+  (pass-if (eqv? +inf.0 (abs +inf.0)))
+  (pass-if (eqv? +inf.0 (abs -inf.0))))
 
 ;;;
 ;;; quotient
@@ -2814,6 +2820,55 @@
   (pass-if "binary double-negation of fixnum-min: equal?"
     (equal? fixnum-min (- 0 (- 0 fixnum-min))))
 
+  (pass-if "signed zeroes"
+    (and (eqv? +0.0 (- -0.0))
+         (eqv? -0.0 (- +0.0))
+         (eqv?  0.0 (-  0.0  0.0))
+         (eqv?  0.0 (-  0.0 -0.0))
+         (eqv?  0.0 (- -0.0 -0.0))
+         (eqv? -0.0 (- -0.0  0.0))))
+
+  (pass-if "exactness propagation"
+    (and (eqv?  3   (- 8 5))
+         (eqv?  3.0 (- 8 5.0))
+         (eqv?  3.0 (- 8.0 5))
+         (eqv?  3.0 (- 8.0 5.0))
+         (eqv? -1/6  (- 1/3 1/2))
+         (eqv? -4.5  (- 1/2 5.0))
+         (eqv?  2.75 (- 3.0 1/4))))
+
+  (pass-if "infinities"
+    (and (eqv? +inf.0 (- +inf.0 -inf.0))
+         (eqv? -inf.0 (- -inf.0 +inf.0))
+         (real-nan?   (- +inf.0 +inf.0))
+         (real-nan?   (- -inf.0 -inf.0))))
+
+  (pass-if "NaNs"
+    (and (real-nan?  (- +nan.0 +nan.0))
+         (real-nan?  (- 0 +nan.0))
+         (real-nan?  (- +nan.0 0))
+         (real-nan?  (- 1 +nan.0))
+         (real-nan?  (- +nan.0 1))
+         (real-nan?  (- -1 +nan.0))
+         (real-nan?  (- +nan.0 -1))
+         (real-nan?  (- -7/2 +nan.0))
+         (real-nan?  (- +nan.0 -7/2))
+         (real-nan?  (- 1e20 +nan.0))
+         (real-nan?  (- +nan.0 1e20))
+         (real-nan?  (- +inf.0 +nan.0))
+         (real-nan?  (- +nan.0 +inf.0))
+         (real-nan?  (- -inf.0 +nan.0))
+         (real-nan?  (- +nan.0 -inf.0))
+         (real-nan?  (- (* fixnum-max 2) +nan.0))
+         (real-nan?  (- +nan.0 (* fixnum-max 2)))))
+
+  (pass-if "(eqv? fixnum-min (- (- fixnum-min)))"
+    (eqv? fixnum-min (- (- fixnum-min))))
+  (pass-if "(eqv? fixnum-min (- 0 (- 0 fixnum-min)))"
+    (eqv? fixnum-min (- 0 (- 0 fixnum-min))))
+  (pass-if "(eqv? fixnum-num (apply - (list (apply - (list fixnum-min)))))"
+    (eqv? fixnum-min (apply - (list (apply - (list fixnum-min))))))
+
   (pass-if "-inum - +bignum"
     (= #x-100000000000000000000000000000001
        (- -1 #x100000000000000000000000000000000)))
-- 
1.5.6.5


[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #4: Support non-real complex numbers with inexact zero imaginary part --]
[-- Type: text/x-diff, Size: 24463 bytes --]

From ebe58755fab58390e508a1ac3e0ceaa66ab354b6 Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Wed, 2 Feb 2011 05:29:55 -0500
Subject: [PATCH] Support non-real complex numbers with inexact zero imaginary part

Add the ability to represent non-real complex numbers whose imaginary
part is an _inexact_ zero (0.0 or -0.0), per R6RS.  Previously, such
numbers were immediately changed into inexact reals.

* libguile/numbers.c: Remove from the list of `General assumptions' in
  numbers.c that objects satisfying SCM_COMPLEXP() have a non-zero
  complex component.  This is no longer true.  Also add a warning
  about another unrelated assumption that is not entirely correct
  (that floor(r) == r implies that mpz_set_d will DTRT; it won't
  if r is infinite).

  (icmplx2str): Always print the imaginary part, even if it is zero.
  Also handle a negative zero imaginary part more gracefully.  It
  now prints 0.0-0.0i, where previously it would print 0.0+-0.0i.

  (mem2ureal): Replace scm_from_double (0.0) with flo0.

  (scm_c_make_rectangular): Always create non-real complex numbers.
  Previously it would create inexact reals if the specified imaginary
  part was zero.

  (scm_make_rectangular): If the imaginary part is an _exact_ 0, return
  the real part unchanged (possibly exact), otherwise return a non-real
  complex number (possibly with an inexact zero imaginary part).
  Previously, it would return an inexact real number whenever the
  imaginary part was any kind of zero.

  (scm_make_polar): If the magnitude is an exact 0, return an exact 0.
  If the angle is an exact 0, return the magnitude unchanged (possibly
  exact).  Otherwise return a non-real complex number (possibly with an
  inexact zero imaginary part).  Previously, it would return a real
  number whenever the imaginary part was any kind of zero.

  (scm_imag_part): Return an exact 0 if applied to a real number.
  Previously it would return an inexact zero if applied to an inexact
  real number.

  (scm_inexact_to_exact): Accept complex numbers with inexact zero
  imaginary part.  In that case, simply use the real part and ignore the
  imaginary part.  Essentially we coerce the inexact zero imaginary part
  to an exact 0.

* test-suite/tests/numbers.test: Add many test cases, and modify
  existing tests as needed to reflect these changes.  Also add a new
  internal predicate: `almost-real-nan?' which tests for a non-real
  complex number with zero imaginary part whose real part is a NaN.

* doc/ref/api-data.texi (Complex Numbers): Update description of complex
  numbers to reflect these changes: non-real complex numbers in Guile
  need not have non-zero imaginary part.  Also, each part of a complex
  number may be any inexact real, not just rationals as was previously
  stated.  Explicitly mention that each part may be an infinity, a NaN,
  or a signed zero.

  (Complex Number Operations): Change the formal parameter names of
  `make-polar' from `x' and `y' to `mag' and `ang'.

* NEWS: Add news entries.
---
 NEWS                          |   42 +++++++++
 doc/ref/api-data.texi         |   20 +++--
 libguile/numbers.c            |   98 ++++++++++++---------
 test-suite/tests/numbers.test |  187 ++++++++++++++++++++++++++++++++++++-----
 4 files changed, 274 insertions(+), 73 deletions(-)

diff --git a/NEWS b/NEWS
index 64d2864..2d3afdc 100644
--- a/NEWS
+++ b/NEWS
@@ -115,6 +115,48 @@ Note that these operators are equivalent to the R6RS integer division
 operators `div', `mod', `div-and-mod', `div0', `mod0', and
 `div0-and-mod0'.
 
+*** Complex number changes
+
+Guile is now able to represent non-real complex numbers whose
+imaginary part is an _inexact_ zero (0.0 or -0.0), per R6RS.
+Previously, such numbers were immediately changed into inexact reals.
+
+(real? 0.0+0.0i) now returns #f, per R6RS, although (zero? 0.0+0.0i)
+still returns #t, per R6RS.  (= 0 0.0+0.0i) and (= 0.0 0.0+0.0i) are
+#t, but the same comparisons using `eqv?' or `equal?' are #f.
+
+Like other non-real numbers, these complex numbers with inexact zero
+imaginary part will raise exceptions is passed to procedures requiring
+reals, such as `<', `>', `<=', `>=', `min', `max', `positive?',
+`negative?', `inf?', `nan?', `finite?', etc.
+
+**** `make-rectangular' changes
+
+scm_make_rectangular `make-rectangular' now returns a real number only
+if the imaginary part is an _exact_ 0.  Previously, it would return a
+real number if the imaginary part was an inexact zero.
+
+scm_c_make_rectangular now always returns a non-real complex number,
+even if the imaginary part is zero.  Previously, it would return a
+real number if the imaginary part was zero.
+
+**** `make-polar' changes
+
+scm_make_polar `make-polar' now returns a real number only if the
+angle or magnitude is an _exact_ 0.  If the magnitude is an exact 0,
+it now returns an exact 0.  Previously, it would return a real
+number if the imaginary part was an inexact zero.
+
+scm_c_make_polar now always returns a non-real complex number, even if
+the imaginary part is 0.0.  Previously, it would return a real number
+if the imaginary part was 0.0.
+
+**** `imag-part' changes
+
+scm_imag_part `imag-part' now returns an exact 0 if applied to an
+inexact real number.  Previously it returned an inexact zero in this
+case.
+
 *** `eqv?' and `equal?' now compare numbers equivalently
 
 scm_equal_p `equal?' now behaves equivalently to scm_eqv_p `eqv?' for
diff --git a/doc/ref/api-data.texi b/doc/ref/api-data.texi
index 1ce9e1e..9b065a5 100755
--- a/doc/ref/api-data.texi
+++ b/doc/ref/api-data.texi
@@ -683,11 +683,15 @@ angle,
 -1@@1.57079 @result{} 0.0-1.0i  (approx)
 @end lisp
 
-Guile represents a complex number with a non-zero imaginary part as a
-pair of inexact rationals, so the real and imaginary parts of a
-complex number have the same properties of inexactness and limited
-precision as single inexact rational numbers.  Guile can not represent
-exact complex numbers with non-zero imaginary parts.
+Guile represents a complex number as a pair of inexact reals, so the
+real and imaginary parts of a complex number have the same properties of
+inexactness and limited precision as single inexact real numbers.
+
+Note that each part of a complex number may contain any inexact real
+value, including the special values @samp{+nan.0}, @samp{+inf.0} and
+@samp{-inf.0}, as well as either of the signed zeroes @samp{0.0} or
+@samp{-0.0}.
+
 
 @deffn {Scheme Procedure} complex? z
 @deffnx {C Function} scm_complex_p (z)
@@ -1077,10 +1081,10 @@ locale-dependent parsing).
 Return a complex number constructed of the given @var{real-part} and @var{imaginary-part} parts.
 @end deffn
 
-@deffn {Scheme Procedure} make-polar x y
-@deffnx {C Function} scm_make_polar (x, y)
+@deffn {Scheme Procedure} make-polar mag ang
+@deffnx {C Function} scm_make_polar (mag, ang)
 @cindex polar form
-Return the complex number @var{x} * e^(i * @var{y}).
+Return the complex number @var{mag} * e^(i * @var{ang}).
 @end deffn
 
 @c begin (texi-doc-string "guile" "real-part")
diff --git a/libguile/numbers.c b/libguile/numbers.c
index 9a7848a..18d5755 100644
--- a/libguile/numbers.c
+++ b/libguile/numbers.c
@@ -22,10 +22,10 @@
 
 \f
 /* General assumptions:
- * All objects satisfying SCM_COMPLEXP() have a non-zero complex component.
  * All objects satisfying SCM_BIGP() are too large to fit in a fixnum.
  * If an object satisfies integer?, it's either an inum, a bignum, or a real.
  * If floor (r) == r, r is an int, and mpz_set_d will DTRT.
+ *     XXX What about infinities?  They are equal to their own floor!  -mhw
  * All objects satisfying SCM_FRACTIONP are never an integer.
  */
 
@@ -3649,17 +3649,20 @@ static size_t
 icmplx2str (double real, double imag, char *str, int radix)
 {
   size_t i;
+  double sgn;
   
   i = idbl2str (real, str, radix);
-  if (imag != 0.0)
-    {
-      /* Don't output a '+' for negative numbers or for Inf and
-	 NaN.  They will provide their own sign. */
-      if (0 <= imag && !isinf (imag) && !isnan (imag))
-	str[i++] = '+';
-      i += idbl2str (imag, &str[i], radix);
-      str[i++] = 'i';
-    }
+#ifdef HAVE_COPYSIGN
+  sgn = copysign (1.0, imag);
+#else
+  sgn = imag;
+#endif
+  /* Don't output a '+' for negative numbers or for Inf and
+     NaN.  They will provide their own sign. */
+  if (sgn >= 0 && DOUBLE_IS_FINITE (imag))
+    str[i++] = '+';
+  i += idbl2str (imag, &str[i], radix);
+  str[i++] = 'i';
   return i;
 }
 
@@ -4206,7 +4209,7 @@ mem2ureal (SCM mem, unsigned int *p_idx,
      floating point value so that we can change its sign. 
   */
   if (scm_is_eq (result, SCM_INUM0) && *p_exactness == INEXACT)
-    result = scm_from_double (0.0);
+    result = flo0;
 
   return result;
 }
@@ -7107,19 +7110,14 @@ SCM_PRIMITIVE_GENERIC (scm_sys_atanh, "atanh", 1, 0, 0,
 SCM
 scm_c_make_rectangular (double re, double im)
 {
-  if (im == 0.0)
-    return scm_from_double (re);
-  else
-    {
-      SCM z;
+  SCM z;
 
-      z = PTR2SCM (scm_gc_malloc_pointerless (sizeof (scm_t_complex),
-					      "complex"));
-      SCM_SET_CELL_TYPE (z, scm_tc16_complex);
-      SCM_COMPLEX_REAL (z) = re;
-      SCM_COMPLEX_IMAG (z) = im;
-      return z;
-    }
+  z = PTR2SCM (scm_gc_malloc_pointerless (sizeof (scm_t_complex),
+					  "complex"));
+  SCM_SET_CELL_TYPE (z, scm_tc16_complex);
+  SCM_COMPLEX_REAL (z) = re;
+  SCM_COMPLEX_IMAG (z) = im;
+  return z;
 }
 
 SCM_DEFINE (scm_make_rectangular, "make-rectangular", 2, 0, 0,
@@ -7132,8 +7130,13 @@ SCM_DEFINE (scm_make_rectangular, "make-rectangular", 2, 0, 0,
                    SCM_ARG1, FUNC_NAME, "real");
   SCM_ASSERT_TYPE (scm_is_real (imaginary_part), imaginary_part,
                    SCM_ARG2, FUNC_NAME, "real");
-  return scm_c_make_rectangular (scm_to_double (real_part),
-                                 scm_to_double (imaginary_part));
+
+  /* Return a real if and only if the imaginary_part is an _exact_ 0 */
+  if (scm_is_eq (imaginary_part, SCM_INUM0))
+    return real_part;
+  else
+    return scm_c_make_rectangular (scm_to_double (real_part),
+				   scm_to_double (imaginary_part));
 }
 #undef FUNC_NAME
 
@@ -7156,13 +7159,21 @@ scm_c_make_polar (double mag, double ang)
 }
 
 SCM_DEFINE (scm_make_polar, "make-polar", 2, 0, 0,
-            (SCM x, SCM y),
-	    "Return the complex number @var{x} * e^(i * @var{y}).")
+            (SCM mag, SCM ang),
+	    "Return the complex number @var{mag} * e^(i * @var{ang}).")
 #define FUNC_NAME s_scm_make_polar
 {
-  SCM_ASSERT_TYPE (scm_is_real (x), x, SCM_ARG1, FUNC_NAME, "real");
-  SCM_ASSERT_TYPE (scm_is_real (y), y, SCM_ARG2, FUNC_NAME, "real");
-  return scm_c_make_polar (scm_to_double (x), scm_to_double (y));
+  SCM_ASSERT_TYPE (scm_is_real (mag), mag, SCM_ARG1, FUNC_NAME, "real");
+  SCM_ASSERT_TYPE (scm_is_real (ang), ang, SCM_ARG2, FUNC_NAME, "real");
+
+  /* If mag is exact0, return exact0 */
+  if (scm_is_eq (mag, SCM_INUM0))
+    return SCM_INUM0;
+  /* Return a real if ang is exact0 */
+  else if (scm_is_eq (ang, SCM_INUM0))
+    return mag;
+  else
+    return scm_c_make_polar (scm_to_double (mag), scm_to_double (ang));
 }
 #undef FUNC_NAME
 
@@ -7189,9 +7200,7 @@ SCM_PRIMITIVE_GENERIC (scm_imag_part, "imag-part", 1, 0, 0,
 {
   if (SCM_COMPLEXP (z))
     return scm_from_double (SCM_COMPLEX_IMAG (z));
-  else if (SCM_REALP (z))
-    return flo0;
-  else if (SCM_I_INUMP (z) || SCM_BIGP (z) || SCM_FRACTIONP (z))
+  else if (SCM_I_INUMP (z) || SCM_REALP (z) || SCM_BIGP (z) || SCM_FRACTIONP (z))
     return SCM_INUM0;
   else
     SCM_WTA_DISPATCH_1 (g_scm_imag_part, z, SCM_ARG1, s_scm_imag_part);
@@ -7344,11 +7353,20 @@ SCM_PRIMITIVE_GENERIC (scm_inexact_to_exact, "inexact->exact", 1, 0, 0,
 	"Return an exact number that is numerically closest to @var{z}.")
 #define FUNC_NAME s_scm_inexact_to_exact
 {
-  if (SCM_I_INUMP (z) || SCM_BIGP (z))
+  if (SCM_I_INUMP (z) || SCM_BIGP (z) || SCM_FRACTIONP (z))
     return z;
-  else if (SCM_REALP (z))
+  else
     {
-      if (!DOUBLE_IS_FINITE (SCM_REAL_VALUE (z)))
+      double val;
+
+      if (SCM_REALP (z))
+	val = SCM_REAL_VALUE (z);
+      else if (SCM_COMPLEXP (z) && SCM_COMPLEX_IMAG (z) == 0.0)
+	val = SCM_COMPLEX_REAL (z);
+      else
+	SCM_WTA_DISPATCH_1 (g_scm_inexact_to_exact, z, 1, s_scm_inexact_to_exact);
+
+      if (!SCM_LIKELY (DOUBLE_IS_FINITE (val)))
 	SCM_OUT_OF_RANGE (1, z);
       else
 	{
@@ -7356,9 +7374,9 @@ SCM_PRIMITIVE_GENERIC (scm_inexact_to_exact, "inexact->exact", 1, 0, 0,
 	  SCM q;
 	  
 	  mpq_init (frac);
-	  mpq_set_d (frac, SCM_REAL_VALUE (z));
+	  mpq_set_d (frac, val);
 	  q = scm_i_make_ratio (scm_i_mpz2num (mpq_numref (frac)),
-			      scm_i_mpz2num (mpq_denref (frac)));
+				scm_i_mpz2num (mpq_denref (frac)));
 
 	  /* When scm_i_make_ratio throws, we leak the memory allocated
 	     for frac...
@@ -7367,10 +7385,6 @@ SCM_PRIMITIVE_GENERIC (scm_inexact_to_exact, "inexact->exact", 1, 0, 0,
 	  return q;
 	}
     }
-  else if (SCM_FRACTIONP (z))
-    return z;
-  else
-    SCM_WTA_DISPATCH_1 (g_scm_inexact_to_exact, z, 1, s_scm_inexact_to_exact);
 }
 #undef FUNC_NAME
 
diff --git a/test-suite/tests/numbers.test b/test-suite/tests/numbers.test
index 5a8b31b..96f37df 100644
--- a/test-suite/tests/numbers.test
+++ b/test-suite/tests/numbers.test
@@ -125,6 +125,14 @@
   (and (real? obj)
        (nan? obj)))
 
+;; return true if OBJ is a non-real complex number
+;; whose real part is a nan, and whose imaginary
+;; part is an inexact zero.
+(define (almost-real-nan? obj)
+  (and (not (real? obj))
+       (nan?  (real-part obj))
+       (zero? (imag-part obj))))
+
 ;; return true if both the real and imaginary
 ;; parts of OBJ are NaNs
 (define (complex-nan? obj)
@@ -137,6 +145,12 @@
   (and (zero? (real-part obj))
        (nan?  (imag-part obj))))
 
+;; return true if OBJ is a non-real complex zero
+(define (complex-zero? obj)
+  (and (zero? obj)
+       (complex? obj)
+       (not (real? obj))))
+
 (define const-e    2.7182818284590452354)
 (define const-e^2  7.3890560989306502274)
 (define const-1/e  0.3678794411714423215)
@@ -1510,14 +1524,15 @@
                 ;; * <digit 10>+ #+ . #* <suffix>
                 ("3#." 30.0) ("3#.e0" 30.0) ("3#.#" 30.0) ("3#.#e0" 30.0)
                 ;; Complex:
-                ("1@0" 1.0) ("1@+0" 1.0) ("1@-0" 1.0)
+                ("1@0" 1) ("1@+0" 1) ("1@-0" 1)
+                ("1.0@0" 1.0+0i) ("1@+0.0" 1+0.0i) ("1.0@-0" 1.0-0i)
                 ("2+3i" ,(+ 2 (* 3 +i))) ("4-5i" ,(- 4 (* 5 +i)))
                 ("1+i" 1+1i) ("1-i" 1-1i) ("+1i" 0+1i) ("-1i" 0-1i)
                 ("+i" +1i) ("-i" -1i)
 		("1.0+.1i" 1.0+0.1i)
 		("1.0-.1i" 1.0-0.1i)
-		(".1+.0i" 0.1)
-		("1.+.0i" 1.0)
+		(".1+.0i" 0.1+0.0i)
+		("1.+.0i" 1.0+0.0i)
 		(".1+.1i" 0.1+0.1i)
 		("1e1+.1i" 10+0.1i)
 		))
@@ -1643,6 +1658,8 @@
   (pass-if (not (integer? +inf.0)))
   (pass-if (not (integer? -inf.0)))
   (pass-if (not (integer? +nan.0)))
+  (pass-if (not (integer? +inf.0-inf.0i)))
+  (pass-if (not (integer? +nan.0+nan.0i)))
   (pass-if (not (integer? 3+4i)))
   (pass-if (not (integer? #\a)))
   (pass-if (not (integer? "a")))
@@ -1708,12 +1725,19 @@
   (pass-if (equal? (- fixnum-min 1) (- fixnum-min 1)))
   (pass-if (equal?  0.0  0.0))
   (pass-if (equal? -0.0 -0.0))
+  (pass-if (equal?  0.0+0.0i  0.0+0.0i))
+  (pass-if (equal?  0.0-0.0i  0.0-0.0i))
+  (pass-if (equal? -0.0+0.0i -0.0+0.0i))
   (pass-if (not (equal? 0 1)))
   (pass-if (not (equal? 0 0.0)))
   (pass-if (not (equal? 1 1.0)))
   (pass-if (not (equal? 0.0 0)))
   (pass-if (not (equal? 1.0 1)))
   (pass-if (not (equal? -1.0 -1)))
+  (pass-if (not (equal? 1.0 1.0+0.0i)))
+  (pass-if (not (equal? 0.0 0.0+0.0i)))
+  (pass-if (not (equal? 0.0+0.0i  0.0-0.0i)))
+  (pass-if (not (equal? 0.0+0.0i -0.0+0.0i)))
   (pass-if (not (equal? fixnum-max (+ 1 fixnum-max))))
   (pass-if (not (equal? (+ 1 fixnum-max) fixnum-max)))
   (pass-if (not (equal? (+ 1 fixnum-max) (+ 2 fixnum-max))))
@@ -1778,12 +1802,21 @@
   (pass-if (eqv? (- fixnum-min 1) (- fixnum-min 1)))
   (pass-if (eqv?  0.0  0.0))
   (pass-if (eqv? -0.0 -0.0))
+  (pass-if (eqv?  0.0+0.0i  0.0+0.0i))
+  (pass-if (eqv?  0.0-0.0i  0.0-0.0i))
+  (pass-if (eqv? -0.0+0.0i -0.0+0.0i))
+  (pass-if (not (eqv? 0.0 -0.0)))
+  (pass-if (not (eqv? 0.0  0.0+0.0i)))
+  (pass-if (not (eqv? 0.0+0.0i  0.0-0.0i)))
+  (pass-if (not (eqv? 0.0+0.0i -0.0+0.0i)))
   (pass-if (not (eqv? 0 1)))
   (pass-if (not (eqv? 0 0.0)))
   (pass-if (not (eqv? 1 1.0)))
   (pass-if (not (eqv? 0.0 0)))
   (pass-if (not (eqv? 1.0 1)))
   (pass-if (not (eqv? -1.0 -1)))
+  (pass-if (not (eqv? 1.0 1.0+0.0i)))
+  (pass-if (not (eqv? 0.0 0.0+0.0i)))
   (pass-if (not (eqv? fixnum-max (+ 1 fixnum-max))))
   (pass-if (not (eqv? (+ 1 fixnum-max) fixnum-max)))
   (pass-if (not (eqv? (+ 1 fixnum-max) (+ 2 fixnum-max))))
@@ -1835,9 +1868,34 @@
 
 (with-test-prefix "="
   (pass-if (documented? =))
-  (pass-if (= 0 0))
   (pass-if (= 7 7))
   (pass-if (= -7 -7))
+  (pass-if (= 1.0 1))
+  (pass-if (= 1 1.0))
+  (pass-if (= -1 -1.0))
+  (pass-if (= 0.0  0.0))
+  (pass-if (= 0.0 -0.0))
+  (pass-if (= 1 1.0+0.0i))
+
+  (pass-if (= 0  0))
+  (pass-if (= 0  0.0))
+  (pass-if (= 0 -0.0))
+  (pass-if (= 0  0.0+0.0i))
+  (pass-if (= 0  0.0-0.0i))
+  (pass-if (= 0  0.0+0.0i))
+  (pass-if (= 0 -0.0-0.0i))
+
+  (pass-if (=  0        0))
+  (pass-if (=  0.0      0))
+  (pass-if (= -0.0      0))
+  (pass-if (=  0.0+0.0i 0))
+  (pass-if (=  0.0-0.0i 0))
+  (pass-if (=  0.0+0.0i 0))
+  (pass-if (= -0.0-0.0i 0))
+
+  (pass-if (= 0.0+0.0i  0.0-0.0i))
+  (pass-if (= 0.0+0.0i -0.0+0.0i))
+
   (pass-if (= (+ 1 fixnum-max) (+ 1 fixnum-max)))
   (pass-if (= (- fixnum-min 1) (- fixnum-min 1)))
   (pass-if (not (= 0 1)))
@@ -2406,13 +2464,28 @@
 
 (with-test-prefix "zero?"
   (pass-if (documented? zero?))
-  (pass-if (zero? 0))
-  (pass-if (not (zero? 7)))
+
+  (pass-if (zero?  0))
+  (pass-if (zero?  0.0))
+  (pass-if (zero? -0.0))
+
+  (pass-if (zero?  0.0+0.0i))
+  (pass-if (zero?  0.0-0.0i))
+  (pass-if (zero?  0.0+0.0i))
+  (pass-if (zero? -0.0-0.0i))
+
+  (pass-if (not (zero?  7)))
   (pass-if (not (zero? -7)))
+  (pass-if (not (zero? 1/7)))
+  (pass-if (not (zero? -inf.0)))
+  (pass-if (not (zero? +inf.0)))
+  (pass-if (not (zero? +nan.0)))
   (pass-if (not (zero? (+ 1 fixnum-max))))
   (pass-if (not (zero? (- 1 fixnum-min))))
   (pass-if (not (zero? 1.3)))
-  (pass-if (not (zero? 3.1+4.2i))))
+  (pass-if (not (zero? 3.1+4.2i)))
+  (pass-if (not (zero? 1.0+0.0i)))
+  (pass-if (not (zero? 0.0-1.0i))))
 
 ;;;
 ;;; positive?
@@ -2789,6 +2862,54 @@
   (pass-if "documented?"
     (documented? +))
 
+  (pass-if "simple"
+    (and (eqv?  7 (+ 3 4))
+         (eqv?  3 (+ 3))
+         (eqv?  0 (+))))
+
+  (pass-if "exactness propagation"
+    (and (eqv?  8   (+ 3 5))
+         (eqv?  8.0 (+ 3 5.0))
+         (eqv?  8.0 (+ 3.0 5))
+         (eqv?  8.0 (+ 3.0 5.0))
+
+         (eqv?  5/6  (+ 1/2 1/3))
+         (eqv?  5.5  (+ 1/2 5.0))
+         (eqv?  3.25 (+ 3.0 1/4))))
+
+  (pass-if "signed zeroes"
+    (and (eqv?  0.0 (+  0.0))
+         (eqv? -0.0 (+ -0.0))
+         (eqv?  0.0 (+  0.0  0.0))
+         (eqv?  0.0 (+  0.0 -0.0))
+         (eqv?  0.0 (+ -0.0  0.0))
+         (eqv? -0.0 (+ -0.0 -0.0))))
+
+  (pass-if "NaNs"
+    (and (real-nan?  (+ +nan.0 +nan.0))
+         (real-nan?  (+ 0 +nan.0))
+         (real-nan?  (+ +nan.0 0))
+         (real-nan?  (+ 1 +nan.0))
+         (real-nan?  (+ +nan.0 1))
+         (real-nan?  (+ -1 +nan.0))
+         (real-nan?  (+ +nan.0 -1))
+         (real-nan?  (+ -7/2 +nan.0))
+         (real-nan?  (+ +nan.0 -7/2))
+         (real-nan?  (+ 1e20 +nan.0))
+         (real-nan?  (+ +nan.0 1e20))
+         (real-nan?  (+ +inf.0 +nan.0))
+         (real-nan?  (+ +nan.0 +inf.0))
+         (real-nan?  (+ -inf.0 +nan.0))
+         (real-nan?  (+ +nan.0 -inf.0))
+         (real-nan?  (+ (* fixnum-max 2) +nan.0))
+         (real-nan?  (+ +nan.0 (* fixnum-max 2)))))
+
+  (pass-if "infinities"
+    (and (eqv? +inf.0 (+ +inf.0 +inf.0))
+         (eqv? -inf.0 (+ -inf.0 -inf.0))
+         (real-nan?   (+ +inf.0 -inf.0))
+         (real-nan?   (+ -inf.0 +inf.0))))
+
   ;; The maximum fixnum on a 32-bit architecture: 2^29 - 1.
   (pass-if "fixnum + fixnum = bignum (32-bit)"
     (eqv? 536870912 (+ 536870910 2)))
@@ -2906,6 +3027,16 @@
     (pass-if (eqv?   fixnum-min (* (* fixnum-min -1) -1)))
     (pass-if (equal? fixnum-min (* (* fixnum-min -1) -1))))
 
+  (with-test-prefix "signed zeroes"
+    (pass-if (eqv? +0.0 (* +0.0 +0.0)))
+    (pass-if (eqv? -0.0 (* -0.0 +0.0)))
+    (pass-if (eqv? +0.0 (* -0.0 -0.0)))
+    (pass-if (eqv? -0.0 (* +0.0 -0.0)))
+    (pass-if (eqv? +0.0+0.0i (* +i +0.0)))
+    (pass-if (eqv? +0.0-0.0i (* -i +0.0)))
+    (pass-if (eqv? -0.0-0.0i (* +i -0.0)))
+    (pass-if (eqv? -0.0+0.0i (* -i -0.0))))
+
   (with-test-prefix "exactness propagation"
     (pass-if (eqv? -0.0 (*  0 -1.0 )))
     (pass-if (eqv?  0.0 (*  0  1.0 )))
@@ -2941,10 +3072,10 @@
     (pass-if (real-nan?      (* (* fixnum-max 2) +nan.0)))
     (pass-if (real-nan?      (* +nan.0 (* fixnum-max 2))))
 
-    (pass-if (real-nan?      (*     0     +nan.0  )))
-    (pass-if (real-nan?      (*  +nan.0      0    )))
-    (pass-if (real-nan?      (*     0     +nan.0+i)))
-    (pass-if (real-nan?      (*  +nan.0+i    0    )))
+    (pass-if (real-nan?        (*     0     +nan.0  )))
+    (pass-if (real-nan?        (*  +nan.0      0    )))
+    (pass-if (almost-real-nan? (*     0     +nan.0+i)))
+    (pass-if (almost-real-nan? (*  +nan.0+i    0    )))
 
     (pass-if (imaginary-nan? (*     0     +nan.0i )))
     (pass-if (imaginary-nan? (*  +nan.0i     0    )))
@@ -2991,13 +3122,13 @@
 
     (pass-if (real-nan?      (*     0     +inf.0  )))
     (pass-if (real-nan?      (*  +inf.0      0    )))
-    (pass-if (real-nan?      (*     0     +inf.0+i)))
-    (pass-if (real-nan?      (*  +inf.0+i    0    )))
-
     (pass-if (real-nan?      (*     0     -inf.0  )))
     (pass-if (real-nan?      (*  -inf.0      0    )))
-    (pass-if (real-nan?      (*     0     -inf.0+i)))
-    (pass-if (real-nan?      (*  -inf.0+i    0    )))
+
+    (pass-if (almost-real-nan? (*     0     +inf.0+i)))
+    (pass-if (almost-real-nan? (*  +inf.0+i    0    )))
+    (pass-if (almost-real-nan? (*     0     -inf.0+i)))
+    (pass-if (almost-real-nan? (*  -inf.0+i    0    )))
 
     (pass-if (imaginary-nan? (*     0     +inf.0i )))
     (pass-if (imaginary-nan? (*  +inf.0i     0    )))
@@ -3466,7 +3597,7 @@
   (pass-if (eqv? -0.125 (expt -2 -3.0)))
   (pass-if (eqv? -0.125 (expt -2.0 -3.0)))
   (pass-if (eqv? 0.25 (expt 2.0 -2.0)))
-  (pass-if (eqv? (* -1.0 12398 12398) (expt +12398i 2.0)))
+  (pass-if (test-eqv? (* -1.0+0.0i 12398 12398) (expt +12398i 2.0)))
   (pass-if (eqv-loosely? +i (expt -1 0.5)))
   (pass-if (eqv-loosely? +i (expt -1 1/2)))
   (pass-if (eqv-loosely? 1.0+1.7320508075688i (expt -8 1/3)))
@@ -3596,6 +3727,11 @@
 ;;;
 ;;; make-rectangular
 ;;;
+  
+(with-test-prefix "make-rectangular"
+  (pass-if (real?      (make-rectangular 5.0  0  )))
+  (pass-if (not (real? (make-rectangular 5.0  0.0))))
+  (pass-if (not (real? (make-rectangular 5.0 -0.0)))))
 
 ;;;
 ;;; make-polar
@@ -3606,10 +3742,15 @@
   (define (almost= x y)
     (> 0.01 (magnitude (- x y))))
   
-  (pass-if (= 0 (make-polar 0 0)))
-  (pass-if (= 0 (make-polar 0 123.456)))
-  (pass-if (= 1 (make-polar 1 0)))
-  (pass-if (= -1 (make-polar -1 0)))
+  (pass-if (real?      (make-polar 0    1.0)))
+  (pass-if (real?      (make-polar 5.0  0  )))
+  (pass-if (not (real? (make-polar 5.0  0.0))))
+  (pass-if (not (real? (make-polar 5.0 -0.0))))
+
+  (pass-if (eqv?  0 (make-polar  0 0)))
+  (pass-if (eqv?  0 (make-polar  0 123.456)))
+  (pass-if (eqv?  1 (make-polar  1 0)))
+  (pass-if (eqv? -1 (make-polar -1 0)))
   
   (pass-if (almost= 0+i (make-polar 1 (* 0.5 pi))))
   (pass-if (almost= -1  (make-polar 1 (* 1.0 pi))))
@@ -3634,7 +3775,7 @@
 
 (with-test-prefix "imag-part"
   (pass-if (documented? imag-part))
-  (pass-if (eqv? 0.0 (imag-part  5.0)))
+  (pass-if (eqv? 0   (imag-part  5.0)))
   (pass-if (eqv? 5.0 (imag-part +5.0i)))
   (pass-if (eqv? 0   (imag-part  5)))
   (pass-if (eqv? 0   (imag-part  1/5)))
@@ -3726,7 +3867,7 @@
   (pass-if (eqv? -1/8 (integer-expt -2 -3)))
   (pass-if (eqv? -0.125 (integer-expt -2.0 -3)))
   (pass-if (eqv? 0.25 (integer-expt 2.0 -2)))
-  (pass-if (eqv? (* -1.0 12398 12398) (integer-expt +12398.0i 2))))
+  (pass-if (test-eqv? (* -1.0+0.0i 12398 12398) (integer-expt +12398.0i 2))))
 
 
 ;;;
-- 
1.5.6.5


^ permalink raw reply related	[flat|nested] 9+ messages in thread

* Re: [PATCH] Complex numbers with inexact zero imaginary part, etc
  2011-02-02 11:25 [PATCH] Complex numbers with inexact zero imaginary part, etc Mark H Weaver
@ 2011-02-02 20:29 ` Andy Wingo
  2011-02-02 20:30 ` Andy Wingo
  2011-02-02 20:35 ` [PATCH] Complex numbers with inexact zero imaginary part, etc Andy Wingo
  2 siblings, 0 replies; 9+ messages in thread
From: Andy Wingo @ 2011-02-02 20:29 UTC (permalink / raw)
  To: Mark H Weaver; +Cc: guile-devel

On Wed 02 Feb 2011 12:25, Mark H Weaver <mhw@netris.org> writes:

> * libguile/numbers.c (scm_abs): (abs -0.0) now returns 0.0.  Previously
>   it returned -0.0.

I applied this, but is it right?  I can convince myself both ways.

Andy
-- 
http://wingolog.org/



^ permalink raw reply	[flat|nested] 9+ messages in thread

* Re: [PATCH] Complex numbers with inexact zero imaginary part, etc
  2011-02-02 11:25 [PATCH] Complex numbers with inexact zero imaginary part, etc Mark H Weaver
  2011-02-02 20:29 ` Andy Wingo
@ 2011-02-02 20:30 ` Andy Wingo
  2011-02-02 21:36   ` Mark H Weaver
  2011-02-02 20:35 ` [PATCH] Complex numbers with inexact zero imaginary part, etc Andy Wingo
  2 siblings, 1 reply; 9+ messages in thread
From: Andy Wingo @ 2011-02-02 20:30 UTC (permalink / raw)
  To: Mark H Weaver; +Cc: guile-devel

Hi,

I just sent a mail I didn't mean to send, I said:

On Wed 02 Feb 2011 12:25, Mark H Weaver <mhw@netris.org> writes:

> * libguile/numbers.c (scm_abs): (abs -0.0) now returns 0.0.
> Previously it returned -0.0.

I questioned this, but I think it's pretty fine, obviously; I meant to
ask about:

>   (scm_difference): (- 0 0.0) now returns -0.0.  Previously it returned
>   0.0.  Also make sure that (- 0 0.0+0.0i) will return -0.0-0.0i.

Is this right?  I can convince myself both ways.

Regards,

Andy
-- 
http://wingolog.org/



^ permalink raw reply	[flat|nested] 9+ messages in thread

* Re: [PATCH] Complex numbers with inexact zero imaginary part, etc
  2011-02-02 11:25 [PATCH] Complex numbers with inexact zero imaginary part, etc Mark H Weaver
  2011-02-02 20:29 ` Andy Wingo
  2011-02-02 20:30 ` Andy Wingo
@ 2011-02-02 20:35 ` Andy Wingo
  2 siblings, 0 replies; 9+ messages in thread
From: Andy Wingo @ 2011-02-02 20:35 UTC (permalink / raw)
  To: Mark H Weaver; +Cc: guile-devel

On Wed 02 Feb 2011 12:25, Mark H Weaver <mhw@netris.org> writes:

> Here's another batch of numerics fixes and changes.

Applied, thanks!

> There are two more patches coming soon: one to improve expt, and one to
> add many test cases to numbers.test.  Hopefully these can make it into
> the imminent prerelease.

They'll make it in shortly after.

Cheers,

Andy
-- 
http://wingolog.org/



^ permalink raw reply	[flat|nested] 9+ messages in thread

* Re: [PATCH] Complex numbers with inexact zero imaginary part, etc
  2011-02-02 20:30 ` Andy Wingo
@ 2011-02-02 21:36   ` Mark H Weaver
  2011-02-02 22:42     ` Andy Wingo
  0 siblings, 1 reply; 9+ messages in thread
From: Mark H Weaver @ 2011-02-02 21:36 UTC (permalink / raw)
  To: Andy Wingo; +Cc: guile-devel

Andy Wingo <wingo@pobox.com> writes:
>>   (scm_difference): (- 0 0.0) now returns -0.0.  Previously it returned
>>   0.0.  Also make sure that (- 0 0.0+0.0i) will return -0.0-0.0i.
>
> Is this right?  I can convince myself both ways.

I'm not 100% confident, but I'm pretty sure it's the right thing.

As far as I can tell, the semantics of the signed zeroes are based on
one-sided limits, so (f 0.0) aka (f +0.0) should be the limit of (f x)
as x approaches zero from the right, similarly for (f -0.0) except
approaching from the left.  (f +inf.0) and (f -inf.0) should also be
defined in terms of limits.  I don't know how this is handled when more
than one of these special values is present in the operand list.

As for the rules of when to return 0.0 vs -0.0, I'm not sure.  IEEE 754
certainly produces -0.0 when a negative number underflows, but -0.0 is
returned in more cases than that.  For example, there is no underflow
for (+ -0.0 -0.0) or (- 0.0).  Based on the rules I know about, the
pattern seems to be that if the operation is defined in terms of a limit
(i.e. if +0.0, -0.0, +inf.0, or -inf.0 appears as an operand), and the
result of the operation is known to approach zero from the left, then
-0.0 is returned, otherwise 0.0 is returned.

In this case, (- 0 0.0) can be defined as the limit of (- 0 x) as x
approaches zero from the right, in which case the result approaches zero
from the left.

Another argument is that some Schemes apparently take the position that
an inexact zero is "weaker" than an exact 0, for example by evaluating
(/ 0 0.0) => 0, as pointed out by Aubrey Jaffer in the "Division by
zero" section of:

  http://people.csail.mit.edu/jaffer/III/RAWI

In any case, the R6RS explicitly gives these examples, which are widely
accepted rules for the signed zeroes:

  (- 0.0)      ==>  -0.0
  (- 0.0 0.0)  ==>   0.0

Therefore, if we decide that (- 0 0.0) should be 0.0 instead of -0.0,
then you'll have to change the compiler to no longer change (- x) to
(- x 0).

      Best,
       Mark



^ permalink raw reply	[flat|nested] 9+ messages in thread

* Re: [PATCH] Complex numbers with inexact zero imaginary part, etc
  2011-02-02 21:36   ` Mark H Weaver
@ 2011-02-02 22:42     ` Andy Wingo
  2011-02-03  0:38       ` [PATCH] Fix non-portable usage of `isinf' in `max' and `min' Mark H Weaver
  0 siblings, 1 reply; 9+ messages in thread
From: Andy Wingo @ 2011-02-02 22:42 UTC (permalink / raw)
  To: Mark H Weaver; +Cc: guile-devel

On Wed 02 Feb 2011 22:36, Mark H Weaver <mhw@netris.org> writes:

>   http://people.csail.mit.edu/jaffer/III/RAWI

Fascinating link, thanks.

I'm OK with the way things are.

BTW: did you see the failures on darwin?

  http://hydra.nixos.org/build/882506/nixlog/1

Seems there were errors in:

    FAIL: numbers.test: max: infinities and NaNs: (real-nan? (max +nan.0 -inf.0))
    FAIL: numbers.test: max: infinities and NaNs: (real-nan? (max -inf.0 +nan.0))
    FAIL: numbers.test: min: infinities and NaNs: (eqv? -inf.0 (min -inf.0 +nan.0))
    FAIL: numbers.test: min: infinities and NaNs: (eqv? -inf.0 (min +nan.0 -inf.0))

Cheers,

Andy
-- 
http://wingolog.org/



^ permalink raw reply	[flat|nested] 9+ messages in thread

* [PATCH] Fix non-portable usage of `isinf' in `max' and `min'
  2011-02-02 22:42     ` Andy Wingo
@ 2011-02-03  0:38       ` Mark H Weaver
  2011-02-03  9:43         ` Andy Wingo
  2011-02-03 13:07         ` Ludovic Courtès
  0 siblings, 2 replies; 9+ messages in thread
From: Mark H Weaver @ 2011-02-03  0:38 UTC (permalink / raw)
  To: Andy Wingo; +Cc: guile-devel

[-- Attachment #1: Type: text/plain, Size: 651 bytes --]

Andy Wingo <wingo@pobox.com> writes:
> BTW: did you see the failures on darwin?
>
>   http://hydra.nixos.org/build/882506/nixlog/1
>
> Seems there were errors in:
>
>     FAIL: numbers.test: max: infinities and NaNs: (real-nan? (max +nan.0 -inf.0))
>     FAIL: numbers.test: max: infinities and NaNs: (real-nan? (max -inf.0 +nan.0))
>     FAIL: numbers.test: min: infinities and NaNs: (eqv? -inf.0 (min -inf.0 +nan.0))
>     FAIL: numbers.test: min: infinities and NaNs: (eqv? -inf.0 (min +nan.0 -inf.0))

Ah, I was using a non-portable extension of isinf(x) to determine the
sign of the infinity.  This patch should fix it.

    Thanks,
      Mark



[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #2: Fix non-portable usage of `isinf' in `max' and `min' --]
[-- Type: text/x-diff, Size: 2252 bytes --]

From 6801f4c8503be81c03f503520c8e2d70944f371d Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Wed, 2 Feb 2011 19:32:16 -0500
Subject: [PATCH] Fix non-portable usage of `isinf' in `max' and `min'

* numbers.c: Add new macros DOUBLE_IS_POSITIVE_INFINITY and
  DOUBLE_IS_NEGATIVE_INFINITY.
  (scm_max, scm_min): Use the new macros to detect particular
  infinities.  Previously we checked the return value of `isinf'
  to determine the sign of the infinity, but that is not portable.
---
 libguile/numbers.c |   13 +++++++++----
 1 files changed, 9 insertions(+), 4 deletions(-)

diff --git a/libguile/numbers.c b/libguile/numbers.c
index 18d5755..3be4478 100644
--- a/libguile/numbers.c
+++ b/libguile/numbers.c
@@ -83,6 +83,11 @@ typedef scm_t_signed_bits scm_t_inum;
    TODO: if it's available, use C99's isfinite(x) instead */
 #define DOUBLE_IS_FINITE(x) (!isinf(x) && !isnan(x))
 
+/* On some platforms, isinf(x) returns 0, 1 or -1, indicating the sign
+   of the infinity, but other platforms return a boolean only. */
+#define DOUBLE_IS_POSITIVE_INFINITY(x) (isinf(x) && ((x) > 0))
+#define DOUBLE_IS_NEGATIVE_INFINITY(x) (isinf(x) && ((x) < 0))
+
 \f
 
 /*
@@ -5251,9 +5256,9 @@ scm_max (SCM x, SCM y)
 	  /* If neither (xx > yy) nor (xx < yy), then
 	     either they're equal or one is a NaN */
 	  else if (SCM_UNLIKELY (isnan (xx)))
-	    return (isinf (yy) == 1) ? y : x;
+	    return DOUBLE_IS_POSITIVE_INFINITY (yy) ? y : x;
 	  else if (SCM_UNLIKELY (isnan (yy)))
-	    return (isinf (xx) == 1) ? x : y;
+	    return DOUBLE_IS_POSITIVE_INFINITY (xx) ? x : y;
 	  /* xx == yy, but handle signed zeroes properly */
 	  else if (double_is_non_negative_zero (yy))
 	    return y;
@@ -5411,9 +5416,9 @@ scm_min (SCM x, SCM y)
 	  /* If neither (xx < yy) nor (xx > yy), then
 	     either they're equal or one is a NaN */
 	  else if (SCM_UNLIKELY (isnan (xx)))
-	    return (isinf (yy) == -1) ? y : x;
+	    return DOUBLE_IS_NEGATIVE_INFINITY (yy) ? y : x;
 	  else if (SCM_UNLIKELY (isnan (yy)))
-	    return (isinf (xx) == -1) ? x : y;
+	    return DOUBLE_IS_NEGATIVE_INFINITY (xx) ? x : y;
 	  /* xx == yy, but handle signed zeroes properly */
 	  else if (double_is_non_negative_zero (xx))
 	    return y;
-- 
1.5.6.5


^ permalink raw reply related	[flat|nested] 9+ messages in thread

* Re: [PATCH] Fix non-portable usage of `isinf' in `max' and `min'
  2011-02-03  0:38       ` [PATCH] Fix non-portable usage of `isinf' in `max' and `min' Mark H Weaver
@ 2011-02-03  9:43         ` Andy Wingo
  2011-02-03 13:07         ` Ludovic Courtès
  1 sibling, 0 replies; 9+ messages in thread
From: Andy Wingo @ 2011-02-03  9:43 UTC (permalink / raw)
  To: Mark H Weaver; +Cc: guile-devel

On Thu 03 Feb 2011 01:38, Mark H Weaver <mhw@netris.org> writes:

> Ah, I was using a non-portable extension of isinf(x) to determine the
> sign of the infinity.  This patch should fix it.

Applied, thanks.

Andy
-- 
http://wingolog.org/



^ permalink raw reply	[flat|nested] 9+ messages in thread

* Re: [PATCH] Fix non-portable usage of `isinf' in `max' and `min'
  2011-02-03  0:38       ` [PATCH] Fix non-portable usage of `isinf' in `max' and `min' Mark H Weaver
  2011-02-03  9:43         ` Andy Wingo
@ 2011-02-03 13:07         ` Ludovic Courtès
  1 sibling, 0 replies; 9+ messages in thread
From: Ludovic Courtès @ 2011-02-03 13:07 UTC (permalink / raw)
  To: guile-devel

Hello!

And thanks for the patch stream and quick fixes!  :-)

Mark H Weaver <mhw@netris.org> writes:

> Ah, I was using a non-portable extension of isinf(x) to determine the
> sign of the infinity.  This patch should fix it.

We use Gnulib’s ‘isinf’ module.  It could be that it’s buggy, or it
could be improved to handle this specific case, in which it would be
nice to share with bug-gnulib@gnu.org.

Ludo’.




^ permalink raw reply	[flat|nested] 9+ messages in thread

end of thread, other threads:[~2011-02-03 13:07 UTC | newest]

Thread overview: 9+ messages (download: mbox.gz follow: Atom feed
-- links below jump to the message on this page --
2011-02-02 11:25 [PATCH] Complex numbers with inexact zero imaginary part, etc Mark H Weaver
2011-02-02 20:29 ` Andy Wingo
2011-02-02 20:30 ` Andy Wingo
2011-02-02 21:36   ` Mark H Weaver
2011-02-02 22:42     ` Andy Wingo
2011-02-03  0:38       ` [PATCH] Fix non-portable usage of `isinf' in `max' and `min' Mark H Weaver
2011-02-03  9:43         ` Andy Wingo
2011-02-03 13:07         ` Ludovic Courtès
2011-02-02 20:35 ` [PATCH] Complex numbers with inexact zero imaginary part, etc Andy Wingo

This is a public inbox, see mirroring instructions
for how to clone and mirror all data and code used for this inbox;
as well as URLs for read-only IMAP folder(s) and NNTP newsgroup(s).