[PATCH] Handle products with exact 0 differently, etc

[PATCH] Handle products with exact 0 differently, etc

[Mark H Weaver]

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

Here's another batch of numerics patches.  The most important one
changes the way products involving exact 0 are handled:

* libguile/numbers.c (scm_product): Handle exact 0 differently.  A
  product containing an exact 0 now returns an exact 0 if and only if
  the other arguments are all exact.  An inexact zero is returned if and
  only if the other arguments are all finite but not all exact.  If an
  infinite or NaN value is present, a NaN value is returned.
  Previously, any product containing an exact 0 yielded an exact 0,
  regardless of the other arguments.

  A note on the rationale for (* 0 0.0) returning 0.0 and not exact 0:
  The exactness propagation rules allow us to return an exact result in
  the presence of inexact arguments only if the values of the inexact
  arguments do not affect the result.  In this case, the value of the
  inexact argument _does_ affect the result, because an infinite or NaN
  value causes the result to be a NaN.

  A note on the rationale for (* 0 +inf.0) being a NaN and not exact 0:
  The R6RS requires that (/ 0 0.0) return a NaN value, and that (/ 0.0)
  return +inf.0.  We would like (/ x y) to be the same as (* x (/ y)),
  and in particular, for (/ 0 0.0) to be the same as (* 0 (/ 0.0)),
  which reduces to (* 0 +inf.0).  Therefore (* 0 +inf.0) should return
  a NaN.

     Best,
      Mark



[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #2: Fix bugs in `rationalize' --]
[-- Type: text/x-diff, Size: 6212 bytes --]

From d4a0dfbaa775f6268a20fde2161911c5ce12e9a9 Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Tue, 1 Feb 2011 05:19:24 -0500
Subject: [PATCH] Fix bugs in `rationalize'

* libguile/numbers.c (scm_rationalize): Fix bugs.  Previously, it
  returned exact integers unmodified, although that was incorrect if
  the epsilon was at least 1 or inexact, e.g. (rationalize 4 1) should
  return 3 per R5RS and R6RS, but previously it returned 4.  Also
  handle cases involving infinities and NaNs properly, per R6RS.

* test-suite/tests/numbers.test: Add test cases for `rationalize'.

* NEWS: Add NEWS entry
---
 NEWS                          |    8 ++++++
 libguile/numbers.c            |   52 +++++++++++++++++++++++++++++++---------
 test-suite/tests/numbers.test |   37 +++++++++++++++++++++++++++++
 3 files changed, 85 insertions(+), 12 deletions(-)

diff --git a/NEWS b/NEWS
index 2ba79a6..3769b81 100644
--- a/NEWS
+++ b/NEWS
@@ -169,6 +169,14 @@ an error when a non-real number or non-number is passed to these
 procedures.  (Note that NaNs _are_ considered numbers by scheme, despite
 their name).
 
+*** `rationalize' bugfixes and changes
+
+Fixed bugs in scm_rationalize `rationalize'.  Previously, it returned
+exact integers unmodified, although that was incorrect if the epsilon
+was at least 1 or inexact, e.g. (rationalize 4 1) should return 3 per
+R5RS and R6RS, but previously it returned 4.  It also now handles
+cases involving infinities and NaNs properly, per R6RS.
+
 *** New procedure: `finite?'
 
 Add scm_finite_p `finite?' from R6RS to guile core, which returns #t
diff --git a/libguile/numbers.c b/libguile/numbers.c
index d08d15f..d4380dd 100644
--- a/libguile/numbers.c
+++ b/libguile/numbers.c
@@ -7267,11 +7267,46 @@ SCM_DEFINE (scm_rationalize, "rationalize", 2, 0, 0,
 	    "@end lisp")
 #define FUNC_NAME s_scm_rationalize
 {
-  if (SCM_I_INUMP (x))
-    return x;
-  else if (SCM_BIGP (x))
+  SCM_ASSERT_TYPE (scm_is_real (x), x, SCM_ARG1, FUNC_NAME, "real");
+  SCM_ASSERT_TYPE (scm_is_real (eps), eps, SCM_ARG2, FUNC_NAME, "real");
+  eps = scm_abs (eps);
+  if (scm_is_false (scm_positive_p (eps)))
+    {
+      /* eps is either zero or a NaN */
+      if (scm_is_true (scm_nan_p (eps)))
+	return scm_nan ();
+      else if (SCM_INEXACTP (eps))
+	return scm_exact_to_inexact (x);
+      else
+	return x;
+    }
+  else if (scm_is_false (scm_finite_p (eps)))
+    {
+      if (scm_is_true (scm_finite_p (x)))
+	return flo0;
+      else
+	return scm_nan ();
+    }
+  else if (scm_is_false (scm_finite_p (x))) /* checks for both inf and nan */
     return x;
-  else if ((SCM_REALP (x)) || SCM_FRACTIONP (x)) 
+  else if (scm_is_false (scm_less_p (scm_floor (scm_sum (x, eps)),
+				     scm_ceiling (scm_difference (x, eps)))))
+    {
+      /* There's an integer within range; we want the one closest to zero */
+      if (scm_is_false (scm_less_p (eps, scm_abs (x))))
+	{
+	  /* zero is within range */
+	  if (SCM_INEXACTP (x) || SCM_INEXACTP (eps))
+	    return flo0;
+	  else
+	    return SCM_INUM0;
+	}
+      else if (scm_is_true (scm_positive_p (x)))
+	return scm_ceiling (scm_difference (x, eps));
+      else
+	return scm_floor (scm_sum (x, eps));
+    }
+  else
     {
       /* Use continued fractions to find closest ratio.  All
 	 arithmetic is done with exact numbers.
@@ -7285,9 +7320,6 @@ SCM_DEFINE (scm_rationalize, "rationalize", 2, 0, 0,
       SCM rx;
       int i = 0;
 
-      if (scm_is_true (scm_num_eq_p (ex, int_part)))
-	return ex;
-      
       ex = scm_difference (ex, int_part);            /* x = x-int_part */
       rx = scm_divide (ex, SCM_UNDEFINED); 	       /* rx = 1/x */
 
@@ -7296,7 +7328,6 @@ SCM_DEFINE (scm_rationalize, "rationalize", 2, 0, 0,
 	 converges after less than a dozen iterations.
       */
 
-      eps = scm_abs (eps);
       while (++i < 1000000)
 	{
 	  a = scm_sum (scm_product (a1, tt), a2);    /* a = a1*tt + a2 */
@@ -7307,8 +7338,7 @@ SCM_DEFINE (scm_rationalize, "rationalize", 2, 0, 0,
 			 eps)))                      /* abs(x-a/b) <= eps */
 	    {
 	      SCM res = scm_sum (int_part, scm_divide (a, b));
-	      if (scm_is_false (scm_exact_p (x))
-		  || scm_is_false (scm_exact_p (eps)))
+	      if (SCM_INEXACTP (x) || SCM_INEXACTP (eps))
 		return scm_exact_to_inexact (res);
 	      else
 		return res;
@@ -7323,8 +7353,6 @@ SCM_DEFINE (scm_rationalize, "rationalize", 2, 0, 0,
 	}
       scm_num_overflow (s_scm_rationalize);
     }
-  else
-    SCM_WRONG_TYPE_ARG (1, x);
 }
 #undef FUNC_NAME
 
diff --git a/test-suite/tests/numbers.test b/test-suite/tests/numbers.test
index d85e44c..5619bf0 100644
--- a/test-suite/tests/numbers.test
+++ b/test-suite/tests/numbers.test
@@ -1328,6 +1328,43 @@
     (pass-if (= lcm-of-big-n-and-11 (lcm 11 big-n 11)))))
 
 ;;;
+;;; rationalize
+;;;
+(with-test-prefix "rationalize"
+  (pass-if (documented? rationalize))
+  (pass-if (eqv?  2     (rationalize  4   2  )))
+  (pass-if (eqv? -2     (rationalize -4   2  )))
+  (pass-if (eqv?  2.0   (rationalize  4   2.0)))
+  (pass-if (eqv? -2.0   (rationalize -4.0 2  )))
+
+  (pass-if (eqv?  0     (rationalize  4   8  )))
+  (pass-if (eqv?  0     (rationalize -4   8  )))
+  (pass-if (eqv?  0.0   (rationalize  4   8.0)))
+  (pass-if (eqv?  0.0   (rationalize -4.0 8  )))
+
+  (pass-if (eqv?  0.0   (rationalize  3   +inf.0)))
+  (pass-if (eqv?  0.0   (rationalize -3   +inf.0)))
+
+  (pass-if (nan?        (rationalize +inf.0 +inf.0)))
+  (pass-if (nan?        (rationalize +nan.0 +inf.0)))
+  (pass-if (nan?        (rationalize +nan.0 4)))
+  (pass-if (eqv? +inf.0 (rationalize +inf.0 3)))
+
+  (pass-if (eqv?  3/10  (rationalize  3/10 0)))
+  (pass-if (eqv? -3/10  (rationalize -3/10 0)))
+
+  (pass-if (eqv?  1/3   (rationalize  3/10 1/10)))
+  (pass-if (eqv? -1/3   (rationalize -3/10 1/10)))
+
+  (pass-if (eqv?  1/3   (rationalize  3/10 -1/10)))
+  (pass-if (eqv? -1/3   (rationalize -3/10 -1/10)))
+
+  (pass-if (test-eqv? (/  1.0 3) (rationalize  0.3  1/10)))
+  (pass-if (test-eqv? (/ -1.0 3) (rationalize -0.3  1/10)))
+  (pass-if (test-eqv? (/  1.0 3) (rationalize  0.3 -1/10)))
+  (pass-if (test-eqv? (/ -1.0 3) (rationalize -0.3 -1/10))))
+
+;;;
 ;;; number->string
 ;;;
 
-- 
1.5.6.5


[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #3: More discriminating NaN predicates for numbers.test --]
[-- Type: text/x-diff, Size: 7593 bytes --]

From ab106861f0bf59f4a71535a745ae5770d4830e3d Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Tue, 1 Feb 2011 05:22:40 -0500
Subject: [PATCH] More discriminating NaN predicates for numbers.test

* test-suite/tests/numbers.test: (real-nan?, complex-nan?,
  imaginary-nan?): Add more discriminating NaN testing predicates
  internal to numbers.test, and convert several uses of `nan?'
  to use these instead:
   * `real-nan?' checks that its argument is real and a NaN.
   * `complex-nan?' checks that both the real and imaginary
                    parts of its argument are NaNs.
   * `imaginary-nan?' checks that its argument's real part
                      is zero and the imaginary part is a NaN.
---
 test-suite/tests/numbers.test |   73 +++++++++++++++++++++++++----------------
 1 files changed, 45 insertions(+), 28 deletions(-)

diff --git a/test-suite/tests/numbers.test b/test-suite/tests/numbers.test
index 5619bf0..75d3790 100644
--- a/test-suite/tests/numbers.test
+++ b/test-suite/tests/numbers.test
@@ -120,6 +120,23 @@
 	 (eqv? x y))
 	(else (and (inexact? y) (> test-epsilon (abs (- x y)))))))
 
+;; return true if OBJ is a real NaN
+(define (real-nan? obj)
+  (and (real? obj)
+       (nan? obj)))
+
+;; return true if both the real and imaginary
+;; parts of OBJ are NaNs
+(define (complex-nan? obj)
+  (and (nan? (real-part obj))
+       (nan? (imag-part obj))))
+
+;; return true if the real part of OBJ is zero
+;; and the imaginary part is a NaN.
+(define (imaginary-nan? obj)
+  (and (zero? (real-part obj))
+       (nan?  (imag-part obj))))
+
 (define const-e    2.7182818284590452354)
 (define const-e^2  7.3890560989306502274)
 (define const-1/e  0.3678794411714423215)
@@ -414,7 +431,7 @@
   (pass-if (= 0.0 (abs 0.0)))
   (pass-if (= 1.0 (abs 1.0)))
   (pass-if (= 1.0 (abs -1.0)))
-  (pass-if (nan? (abs +nan.0)))
+  (pass-if (real-nan? (abs +nan.0)))
   (pass-if (= +inf.0 (abs +inf.0)))
   (pass-if (= +inf.0 (abs -inf.0))))
 
@@ -1345,9 +1362,9 @@
   (pass-if (eqv?  0.0   (rationalize  3   +inf.0)))
   (pass-if (eqv?  0.0   (rationalize -3   +inf.0)))
 
-  (pass-if (nan?        (rationalize +inf.0 +inf.0)))
-  (pass-if (nan?        (rationalize +nan.0 +inf.0)))
-  (pass-if (nan?        (rationalize +nan.0 4)))
+  (pass-if (real-nan?   (rationalize +inf.0 +inf.0)))
+  (pass-if (real-nan?   (rationalize +nan.0 +inf.0)))
+  (pass-if (real-nan?   (rationalize +nan.0 4)))
   (pass-if (eqv? +inf.0 (rationalize +inf.0 3)))
 
   (pass-if (eqv?  3/10  (rationalize  3/10 0)))
@@ -2462,10 +2479,10 @@
       (pass-if (= 5/2 (max 5/2 2))))
 
     (with-test-prefix "inum / real"
-      (pass-if (nan? (max 123 +nan.0))))
+      (pass-if (real-nan? (max 123 +nan.0))))
 
     (with-test-prefix "real / inum"
-      (pass-if (nan? (max +nan.0 123))))
+      (pass-if (real-nan? (max +nan.0 123))))
 
     (with-test-prefix "big / frac"
       (pass-if (= big*2 (max big*2 5/2)))
@@ -2476,14 +2493,14 @@
       (pass-if (= 5/2 (max 5/2 (- big*2)))))
 
     (with-test-prefix "big / real"
-      (pass-if (nan? (max big*5 +nan.0)))
+      (pass-if (real-nan? (max big*5 +nan.0)))
       (pass-if (eqv? (exact->inexact big*5)  (max big*5 -inf.0)))
       (pass-if (eqv? (exact->inexact big*5)  (max big*5 1.0)))
       (pass-if (eqv? +inf.0                  (max big*5 +inf.0)))
       (pass-if (eqv? 1.0                     (max (- big*5) 1.0))))
 
     (with-test-prefix "real / big"
-      (pass-if (nan? (max +nan.0 big*5)))
+      (pass-if (real-nan? (max +nan.0 big*5)))
       (pass-if (eqv? (exact->inexact big*5)  (max -inf.0 big*5)))
       (pass-if (eqv? (exact->inexact big*5)  (max 1.0 big*5)))
       (pass-if (eqv? +inf.0                  (max +inf.0 big*5)))
@@ -2496,9 +2513,9 @@
       (pass-if (= -1/2 (max -2/3 -1/2))))
 
     (with-test-prefix "real / real"
-      (pass-if (nan? (max 123.0 +nan.0)))
-      (pass-if (nan? (max +nan.0 123.0)))
-      (pass-if (nan? (max +nan.0 +nan.0)))
+      (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)))))
 
@@ -2522,8 +2539,8 @@
 
   ;; in gmp prior to 4.2, mpz_cmp_d ended up treating NaN as 3*2^1023, make
   ;; sure we've avoided that
-  (pass-if (nan? (max (ash 1 2048) +nan.0)))
-  (pass-if (nan? (max +nan.0 (ash 1 2048)))))
+  (pass-if (real-nan? (max (ash 1 2048) +nan.0)))
+  (pass-if (real-nan? (max +nan.0 (ash 1 2048)))))
 
 ;;;
 ;;; min
@@ -2587,10 +2604,10 @@
       (pass-if (= 2   (min 5/2 2))))
 
     (with-test-prefix "inum / real"
-      (pass-if (nan? (min 123 +nan.0))))
+      (pass-if (real-nan? (min 123 +nan.0))))
 
     (with-test-prefix "real / inum"
-      (pass-if (nan? (min +nan.0 123))))
+      (pass-if (real-nan? (min +nan.0 123))))
 
     (with-test-prefix "big / frac"
       (pass-if (= 5/2       (min big*2 5/2)))
@@ -2601,14 +2618,14 @@
       (pass-if (= (- big*2) (min 5/2 (- big*2)))))
 
     (with-test-prefix "big / real"
-      (pass-if (nan? (min big*5 +nan.0)))
+      (pass-if (real-nan? (min big*5 +nan.0)))
       (pass-if (eqv? (exact->inexact big*5)      (min big*5  +inf.0)))
       (pass-if (eqv? -inf.0                      (min big*5  -inf.0)))
       (pass-if (eqv? 1.0                         (min big*5 1.0)))
       (pass-if (eqv? (exact->inexact (- big*5))  (min (- big*5) 1.0))))
 
     (with-test-prefix "real / big"
-      (pass-if (nan? (min +nan.0 big*5)))
+      (pass-if (real-nan? (min +nan.0 big*5)))
       (pass-if (eqv? (exact->inexact big*5)      (min +inf.0 big*5)))
       (pass-if (eqv? -inf.0                      (min -inf.0 big*5)))
       (pass-if (eqv? 1.0                         (min 1.0 big*5)))
@@ -2621,9 +2638,9 @@
       (pass-if (= -2/3 (min -2/3 -1/2))))
 
     (with-test-prefix "real / real"
-      (pass-if (nan? (min 123.0 +nan.0)))
-      (pass-if (nan? (min +nan.0 123.0)))
-      (pass-if (nan? (min +nan.0 +nan.0)))
+      (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)))))
 
@@ -2648,8 +2665,8 @@
 
   ;; in gmp prior to 4.2, mpz_cmp_d ended up treating NaN as 3*2^1023, make
   ;; sure we've avoided that
-  (pass-if (nan? (min (- (ash 1 2048)) (- +nan.0))))
-  (pass-if (nan? (min (- +nan.0) (- (ash 1 2048))))))
+  (pass-if (real-nan? (min (- (ash 1 2048)) (- +nan.0))))
+  (pass-if (real-nan? (min (- +nan.0) (- (ash 1 2048))))))
 
 ;;;
 ;;; +
@@ -3166,10 +3183,10 @@
   (pass-if (eqv? 1 (expt 0.0 0)))
   (pass-if (eqv? 1.0 (expt 0 0.0)))
   (pass-if (eqv? 1.0 (expt 0.0 0.0)))
-  (pass-if (nan? (expt 0 -1)))
-  (pass-if (nan? (expt 0 -1.0)))
-  (pass-if (nan? (expt 0.0 -1)))
-  (pass-if (nan? (expt 0.0 -1.0)))
+  (pass-if (real-nan? (expt 0 -1)))
+  (pass-if (real-nan? (expt 0 -1.0)))
+  (pass-if (real-nan? (expt 0.0 -1)))
+  (pass-if (real-nan? (expt 0.0 -1.0)))
   (pass-if (eqv? 0 (expt 0 3)))
   (pass-if (= 0 (expt 0 4.0)))
   (pass-if (eqv? 0.0 (expt 0.0 5)))
@@ -3336,8 +3353,8 @@
 
   (pass-if (eqv? 1 (integer-expt 0 0)))
   (pass-if (eqv? 1 (integer-expt 0.0 0)))
-  (pass-if (nan? (integer-expt 0 -1)))
-  (pass-if (nan? (integer-expt 0.0 -1)))
+  (pass-if (real-nan? (integer-expt 0 -1)))
+  (pass-if (real-nan? (integer-expt 0.0 -1)))
   (pass-if (eqv? 0 (integer-expt 0 3)))
   (pass-if (eqv? 0.0 (integer-expt 0.0 5)))
   (pass-if (eqv? -2742638075.5 (integer-expt -2742638075.5 1)))
-- 
1.5.6.5


[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #4: Handle products with exact 0 differently --]
[-- Type: text/x-diff, Size: 11715 bytes --]

From d4c1abe9bfd5397602f4d2c00ffd66a8ff133b01 Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Tue, 1 Feb 2011 06:30:29 -0500
Subject: [PATCH] Handle products with exact 0 differently

* libguile/numbers.c (scm_product): Handle exact 0 differently.  A
  product containing an exact 0 now returns an exact 0 if and only if
  the other arguments are all exact.  An inexact zero is returned if and
  only if the other arguments are all finite but not all exact.  If an
  infinite or NaN value is present, a NaN value is returned.
  Previously, any product containing an exact 0 yielded an exact 0,
  regardless of the other arguments.

  A note on the rationale for (* 0 0.0) returning 0.0 and not exact 0:
  The exactness propagation rules allow us to return an exact result in
  the presence of inexact arguments only if the values of the inexact
  arguments do not affect the result.  In this case, the value of the
  inexact argument _does_ affect the result, because an infinite or NaN
  value causes the result to be a NaN.

  A note on the rationale for (* 0 +inf.0) being a NaN and not exact 0:
  The R6RS requires that (/ 0 0.0) return a NaN value, and that (/ 0.0)
  return +inf.0.  We would like (/ x y) to be the same as (* x (/ y)),
  and in particular, for (/ 0 0.0) to be the same as (* 0 (/ 0.0)),
  which reduces to (* 0 +inf.0).  Therefore (* 0 +inf.0) should return
  a NaN.

* test-suite/tests/numbers.test: Add many multiplication tests.

* NEWS: Add NEWS entry.
---
 NEWS                          |   10 +++
 libguile/numbers.c            |   56 +++++++++++------
 test-suite/tests/numbers.test |  129 +++++++++++++++++++++++++++++++++++++----
 3 files changed, 163 insertions(+), 32 deletions(-)

diff --git a/NEWS b/NEWS
index 3769b81..63df7db 100644
--- a/NEWS
+++ b/NEWS
@@ -130,6 +130,16 @@ Previously, `(equal? +nan.0 +nan.0)' returned #f, although
 both returned #t.  R5RS requires that `equal?' behave like
 `eqv?' when comparing numbers.
 
+*** Change in handling products `*' involving exact 0
+
+scm_product `*' now handles exact 0 differently.  A product containing
+an exact 0 now returns an exact 0 if and only if the other arguments
+are all exact.  An inexact zero is returned if and only if the other
+arguments are all finite but not all exact.  If an infinite or NaN
+value is present, a NaN value is returned.  Previously, any product
+containing an exact 0 yielded an exact 0, regardless of the other
+arguments.
+
 *** `expt' and `integer-expt' changes when the base is 0
 
 While `(expt 0 0)' is still 1, and `(expt 0 N)' for N > 0 is still
diff --git a/libguile/numbers.c b/libguile/numbers.c
index d4380dd..9ba340f 100644
--- a/libguile/numbers.c
+++ b/libguile/numbers.c
@@ -5900,22 +5900,43 @@ scm_product (SCM x, SCM y)
     {
       scm_t_inum xx;
 
-    intbig:
+    xinum:
       xx = SCM_I_INUM (x);
 
       switch (xx)
 	{
-        case 0: return x; break;
-        case 1: return y; break;
+        case 1:
+	  /* exact1 is the universal multiplicative identity */
+	  return y;
+	  break;
+        case 0:
+	  /* exact0 times a fixnum is exact0: optimize this case */
+	  if (SCM_LIKELY (SCM_I_INUMP (y)))
+	    return SCM_INUM0;
+	  /* if the other argument is inexact, the result is inexact,
+	     and we must do the multiplication in order to handle
+	     infinities and NaNs properly. */
+	  else if (SCM_REALP (y))
+	    return scm_from_double (0.0 * SCM_REAL_VALUE (y));
+	  else if (SCM_COMPLEXP (y))
+	    return scm_c_make_rectangular (0.0 * SCM_COMPLEX_REAL (y),
+					   0.0 * SCM_COMPLEX_IMAG (y));
+	  /* we've already handled inexact numbers,
+	     so y must be exact, and we return exact0 */
+	  else if (SCM_NUMP (y))
+	    return SCM_INUM0;
+	  else
+	    SCM_WTA_DISPATCH_2 (g_product, x, y, SCM_ARGn, s_product);
+	  break;
+        case -1:
 	  /*
-	   * The following case (x = -1) is important for more than
-	   * just optimization.  It handles the case of negating
+	   * This case is important for more than just optimization.
+	   * It handles the case of negating
 	   * (+ 1 most-positive-fixnum) aka (- most-negative-fixnum),
 	   * which is a bignum that must be changed back into a fixnum.
 	   * Failure to do so will cause the following to return #f:
 	   * (= most-negative-fixnum (* -1 (- most-negative-fixnum)))
 	   */
-        case -1:
 	  return scm_difference(y, SCM_UNDEFINED);
 	  break;
 	}
@@ -5957,7 +5978,7 @@ scm_product (SCM x, SCM y)
       if (SCM_I_INUMP (y))
 	{
 	  SCM_SWAP (x, y);
-	  goto intbig;
+	  goto xinum;
 	}
       else if (SCM_BIGP (y))
 	{
@@ -5990,12 +6011,10 @@ scm_product (SCM x, SCM y)
   else if (SCM_REALP (x))
     {
       if (SCM_I_INUMP (y))
-        {
-          /* inexact*exact0 => exact 0, per R5RS "Exactness" section */
-          if (scm_is_eq (y, SCM_INUM0))
-            return y;
-          return scm_from_double (SCM_I_INUM (y) * SCM_REAL_VALUE (x));
-        }
+	{
+	  SCM_SWAP (x, y);
+	  goto xinum;
+	}
       else if (SCM_BIGP (y))
 	{
 	  double result = mpz_get_d (SCM_I_BIG_MPZ (y)) * SCM_REAL_VALUE (x);
@@ -6015,13 +6034,10 @@ scm_product (SCM x, SCM y)
   else if (SCM_COMPLEXP (x))
     {
       if (SCM_I_INUMP (y))
-        {
-          /* inexact*exact0 => exact 0, per R5RS "Exactness" section */
-          if (scm_is_eq (y, SCM_INUM0))
-            return y;
-          return scm_c_make_rectangular (SCM_I_INUM (y) * SCM_COMPLEX_REAL (x),
-                                         SCM_I_INUM (y) * SCM_COMPLEX_IMAG (x));
-        }
+	{
+	  SCM_SWAP (x, y);
+	  goto xinum;
+	}
       else if (SCM_BIGP (y))
 	{
 	  double z = mpz_get_d (SCM_I_BIG_MPZ (y));
diff --git a/test-suite/tests/numbers.test b/test-suite/tests/numbers.test
index 75d3790..96fb6d9 100644
--- a/test-suite/tests/numbers.test
+++ b/test-suite/tests/numbers.test
@@ -2745,6 +2745,115 @@
     (pass-if (eqv?   fixnum-min (* (* fixnum-min -1) -1)))
     (pass-if (equal? fixnum-min (* (* fixnum-min -1) -1))))
 
+  (with-test-prefix "exactness propagation"
+    (pass-if (eqv? -0.0 (*  0 -1.0 )))
+    (pass-if (eqv?  0.0 (*  0  1.0 )))
+    (pass-if (eqv? -0.0 (* -1.0  0 )))
+    (pass-if (eqv?  0.0 (*  1.0  0 )))
+    (pass-if (eqv?  0   (*  0  1/2 )))
+    (pass-if (eqv?  0   (*  1/2  0 )))
+    (pass-if (eqv?  0.0+0.0i (*  0  1+i )))
+    (pass-if (eqv?  0.0+0.0i (*  1+i  0 )))
+    (pass-if (eqv? -1.0 (*  1 -1.0 )))
+    (pass-if (eqv?  1.0 (*  1  1.0 )))
+    (pass-if (eqv? -1.0 (* -1.0  1 )))
+    (pass-if (eqv?  1.0 (*  1.0  1 )))
+    (pass-if (eqv?  1/2 (*  1  1/2 )))
+    (pass-if (eqv?  1/2 (*  1/2  1 )))
+    (pass-if (eqv?  1+i (*  1  1+i )))
+    (pass-if (eqv?  1+i (*  1+i  1 ))))
+
+  (with-test-prefix "propagation of NaNs"
+    (pass-if (real-nan?      (* +nan.0 +nan.0)))
+    (pass-if (real-nan?      (* +nan.0    1  )))
+    (pass-if (real-nan?      (* +nan.0   -1  )))
+    (pass-if (real-nan?      (* +nan.0 -7/2  )))
+    (pass-if (real-nan?      (* +nan.0 1e20  )))
+    (pass-if (real-nan?      (*  1     +nan.0)))
+    (pass-if (real-nan?      (* -1     +nan.0)))
+    (pass-if (real-nan?      (* -7/2   +nan.0)))
+    (pass-if (real-nan?      (* 1e20   +nan.0)))
+    (pass-if (real-nan?      (* +inf.0 +nan.0)))
+    (pass-if (real-nan?      (* +nan.0 +inf.0)))
+    (pass-if (real-nan?      (* -inf.0 +nan.0)))
+    (pass-if (real-nan?      (* +nan.0 -inf.0)))
+    (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 (imaginary-nan? (*     0     +nan.0i )))
+    (pass-if (imaginary-nan? (*  +nan.0i     0    )))
+    (pass-if (imaginary-nan? (*     0    1+nan.0i )))
+    (pass-if (imaginary-nan? (* 1+nan.0i     0    )))
+
+    (pass-if (complex-nan?   (* 0   +nan.0+nan.0i )))
+    (pass-if (complex-nan?   (* +nan.0+nan.0i   0 ))))
+
+  (with-test-prefix "infinities"
+    (pass-if (eqv?   +inf.0  (* +inf.0  5  )))
+    (pass-if (eqv?   -inf.0  (* +inf.0 -5  )))
+    (pass-if (eqv?   +inf.0  (* +inf.0 73.1)))
+    (pass-if (eqv?   -inf.0  (* +inf.0 -9.2)))
+    (pass-if (eqv?   +inf.0  (* +inf.0  5/2)))
+    (pass-if (eqv?   -inf.0  (* +inf.0 -5/2)))
+    (pass-if (eqv?   -inf.0  (* -5   +inf.0)))
+    (pass-if (eqv?   +inf.0  (* 73.1 +inf.0)))
+    (pass-if (eqv?   -inf.0  (* -9.2 +inf.0)))
+    (pass-if (eqv?   +inf.0  (*  5/2 +inf.0)))
+    (pass-if (eqv?   -inf.0  (* -5/2 +inf.0)))
+
+    (pass-if (eqv?   -inf.0  (* -inf.0  5  )))
+    (pass-if (eqv?   +inf.0  (* -inf.0 -5  )))
+    (pass-if (eqv?   -inf.0  (* -inf.0 73.1)))
+    (pass-if (eqv?   +inf.0  (* -inf.0 -9.2)))
+    (pass-if (eqv?   -inf.0  (* -inf.0  5/2)))
+    (pass-if (eqv?   +inf.0  (* -inf.0 -5/2)))
+    (pass-if (eqv?   +inf.0  (* -5   -inf.0)))
+    (pass-if (eqv?   -inf.0  (* 73.1 -inf.0)))
+    (pass-if (eqv?   +inf.0  (* -9.2 -inf.0)))
+    (pass-if (eqv?   -inf.0  (* 5/2  -inf.0)))
+    (pass-if (eqv?   +inf.0  (* -5/2 -inf.0)))
+
+    (pass-if (real-nan?      (*    0.0 +inf.0)))
+    (pass-if (real-nan?      (*   -0.0 +inf.0)))
+    (pass-if (real-nan?      (* +inf.0    0.0)))
+    (pass-if (real-nan?      (* +inf.0   -0.0)))
+
+    (pass-if (real-nan?      (*    0.0 -inf.0)))
+    (pass-if (real-nan?      (*   -0.0 -inf.0)))
+    (pass-if (real-nan?      (* -inf.0    0.0)))
+    (pass-if (real-nan?      (* -inf.0   -0.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 (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 (imaginary-nan? (*     0     +inf.0i )))
+    (pass-if (imaginary-nan? (*  +inf.0i     0    )))
+    (pass-if (imaginary-nan? (*     0    1+inf.0i )))
+    (pass-if (imaginary-nan? (* 1+inf.0i     0    )))
+
+    (pass-if (imaginary-nan? (*     0     -inf.0i )))
+    (pass-if (imaginary-nan? (*  -inf.0i     0    )))
+    (pass-if (imaginary-nan? (*     0    1-inf.0i )))
+    (pass-if (imaginary-nan? (* 1-inf.0i     0    )))
+
+    (pass-if (complex-nan?   (* 0   +inf.0+inf.0i )))
+    (pass-if (complex-nan?   (* +inf.0+inf.0i   0 )))
+
+    (pass-if (complex-nan?   (* 0   +inf.0-inf.0i )))
+    (pass-if (complex-nan?   (* -inf.0+inf.0i   0 ))))
+
   (with-test-prefix "inum * bignum"
 
     (pass-if "0 * 2^256 = 0"
@@ -2752,13 +2861,13 @@
 
   (with-test-prefix "inum * flonum"
 
-    (pass-if "0 * 1.0 = 0"
-      (eqv? 0 (* 0 1.0))))
+    (pass-if "0 * 1.0 = 0.0"
+      (eqv? 0.0 (* 0 1.0))))
 
   (with-test-prefix "inum * complex"
 
-    (pass-if "0 * 1+1i = 0"
-      (eqv? 0 (* 0 1+1i))))
+    (pass-if "0 * 1+1i = 0.0+0.0i"
+      (eqv? 0.0+0.0i (* 0 1+1i))))
 
   (with-test-prefix "inum * frac"
 
@@ -2771,16 +2880,12 @@
       (eqv? 0 (* (ash 1 256) 0))))
 
   (with-test-prefix "flonum * inum"
-
-    ;; in guile 1.6.8 and 1.8.1 and earlier this returned inexact 0.0
-    (pass-if "1.0 * 0 = 0"
-      (eqv? 0 (* 1.0 0))))
+    (pass-if "1.0 * 0 = 0.0"
+      (eqv? 0.0 (* 1.0 0))))
 
   (with-test-prefix "complex * inum"
-
-    ;; in guile 1.6.8 and 1.8.1 and earlier this returned inexact 0.0
-    (pass-if "1+1i * 0 = 0"
-      (eqv? 0 (* 1+1i 0))))
+    (pass-if "1+1i * 0 = 0.0+0.0i"
+      (eqv? 0.0+0.0i (* 1+1i 0))))
 
   (pass-if "complex * bignum"
     (let ((big (ash 1 90)))
-- 
1.5.6.5


[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #5: Move comment about trig functions back where it belongs --]
[-- Type: text/x-diff, Size: 1612 bytes --]

From c7d7dec54e0c6ff75d3a98cc2e5f4e750e9c5e62 Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Tue, 1 Feb 2011 06:50:48 -0500
Subject: [PATCH] Move comment about trig functions back where it belongs

* libguile/numbers.c: Move a comment about the trigonometric functions
  next to those functions.  At some point they became separated, when
  scm_expt was placed between them.
---
 libguile/numbers.c |   12 ++++++------
 1 files changed, 6 insertions(+), 6 deletions(-)

diff --git a/libguile/numbers.c b/libguile/numbers.c
index 9ba340f..f9e00e6 100644
--- a/libguile/numbers.c
+++ b/libguile/numbers.c
@@ -6692,12 +6692,6 @@ SCM_PRIMITIVE_GENERIC (scm_ceiling, "ceiling", 1, 0, 0,
 }
 #undef FUNC_NAME
 
-/* sin/cos/tan/asin/acos/atan
-   sinh/cosh/tanh/asinh/acosh/atanh
-   Derived from "Transcen.scm", Complex trancendental functions for SCM.
-   Written by Jerry D. Hedden, (C) FSF.
-   See the file `COPYING' for terms applying to this program. */
-
 SCM_PRIMITIVE_GENERIC (scm_expt, "expt", 2, 0, 0,
 		       (SCM x, SCM y),
 		       "Return @var{x} raised to the power of @var{y}.")
@@ -6739,6 +6733,12 @@ SCM_PRIMITIVE_GENERIC (scm_expt, "expt", 2, 0, 0,
 }
 #undef FUNC_NAME
 
+/* sin/cos/tan/asin/acos/atan
+   sinh/cosh/tanh/asinh/acosh/atanh
+   Derived from "Transcen.scm", Complex trancendental functions for SCM.
+   Written by Jerry D. Hedden, (C) FSF.
+   See the file `COPYING' for terms applying to this program. */
+
 SCM_PRIMITIVE_GENERIC (scm_sin, "sin", 1, 0, 0,
                        (SCM z),
                        "Compute the sine of @var{z}.")
-- 
1.5.6.5


[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #6: Trigonometric functions return exact numbers in some cases --]
[-- Type: text/x-diff, Size: 9439 bytes --]

From 702c1210e420a0fcd68b9c62f85633c5401a3a28 Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Tue, 1 Feb 2011 06:56:02 -0500
Subject: [PATCH] Trigonometric functions return exact numbers in some cases

* libguile/numbers.c (scm_sin, scm_cos, scm_tan, scm_asin, scm_acos,
  scm_atan, scm_sinh, scm_cosh, scm_tanh, scm_sys_asinh, scm_sys_acosh,
  scm_sys_atanh): Return an exact result in some cases.

* test-suite/tests/numbers.test: Add test cases.

* NEWS: Add NEWS entry
---
 NEWS                          |    7 +++
 libguile/numbers.c            |   48 ++++++++++++++-----
 test-suite/tests/numbers.test |  102 +++++++++++++++++++++++++++++++++++++++-
 3 files changed, 142 insertions(+), 15 deletions(-)

diff --git a/NEWS b/NEWS
index 63df7db..64d2864 100644
--- a/NEWS
+++ b/NEWS
@@ -187,6 +187,13 @@ was at least 1 or inexact, e.g. (rationalize 4 1) should return 3 per
 R5RS and R6RS, but previously it returned 4.  It also now handles
 cases involving infinities and NaNs properly, per R6RS.
 
+*** Trigonometric functions now return exact numbers in some cases
+
+scm_sin `sin', scm_cos `cos', scm_tan `tan', scm_asin `asin', scm_acos
+`acos', scm_atan `atan', scm_sinh `sinh', scm_cosh `cosh', scm_tanh
+`tanh', scm_sys_asinh `asinh', scm_sys_acosh `acosh', and
+scm_sys_atanh `atanh' now return exact results in some cases.
+
 *** New procedure: `finite?'
 
 Add scm_finite_p `finite?' from R6RS to guile core, which returns #t
diff --git a/libguile/numbers.c b/libguile/numbers.c
index f9e00e6..df95c32 100644
--- a/libguile/numbers.c
+++ b/libguile/numbers.c
@@ -6744,7 +6744,9 @@ SCM_PRIMITIVE_GENERIC (scm_sin, "sin", 1, 0, 0,
                        "Compute the sine of @var{z}.")
 #define FUNC_NAME s_scm_sin
 {
-  if (scm_is_real (z))
+  if (SCM_UNLIKELY (scm_is_eq (z, SCM_INUM0)))
+    return z;  /* sin(exact0) = exact0 */
+  else if (scm_is_real (z))
     return scm_from_double (sin (scm_to_double (z)));
   else if (SCM_COMPLEXP (z))
     { double x, y;
@@ -6763,7 +6765,9 @@ SCM_PRIMITIVE_GENERIC (scm_cos, "cos", 1, 0, 0,
                        "Compute the cosine of @var{z}.")
 #define FUNC_NAME s_scm_cos
 {
-  if (scm_is_real (z))
+  if (SCM_UNLIKELY (scm_is_eq (z, SCM_INUM0)))
+    return SCM_INUM1;  /* cos(exact0) = exact1 */
+  else if (scm_is_real (z))
     return scm_from_double (cos (scm_to_double (z)));
   else if (SCM_COMPLEXP (z))
     { double x, y;
@@ -6782,7 +6786,9 @@ SCM_PRIMITIVE_GENERIC (scm_tan, "tan", 1, 0, 0,
                        "Compute the tangent of @var{z}.")
 #define FUNC_NAME s_scm_tan
 {
-  if (scm_is_real (z))
+  if (SCM_UNLIKELY (scm_is_eq (z, SCM_INUM0)))
+    return z;  /* tan(exact0) = exact0 */
+  else if (scm_is_real (z))
     return scm_from_double (tan (scm_to_double (z)));
   else if (SCM_COMPLEXP (z))
     { double x, y, w;
@@ -6805,7 +6811,9 @@ SCM_PRIMITIVE_GENERIC (scm_sinh, "sinh", 1, 0, 0,
                        "Compute the hyperbolic sine of @var{z}.")
 #define FUNC_NAME s_scm_sinh
 {
-  if (scm_is_real (z))
+  if (SCM_UNLIKELY (scm_is_eq (z, SCM_INUM0)))
+    return z;  /* sinh(exact0) = exact0 */
+  else if (scm_is_real (z))
     return scm_from_double (sinh (scm_to_double (z)));
   else if (SCM_COMPLEXP (z))
     { double x, y;
@@ -6824,7 +6832,9 @@ SCM_PRIMITIVE_GENERIC (scm_cosh, "cosh", 1, 0, 0,
                        "Compute the hyperbolic cosine of @var{z}.")
 #define FUNC_NAME s_scm_cosh
 {
-  if (scm_is_real (z))
+  if (SCM_UNLIKELY (scm_is_eq (z, SCM_INUM0)))
+    return SCM_INUM1;  /* cosh(exact0) = exact1 */
+  else if (scm_is_real (z))
     return scm_from_double (cosh (scm_to_double (z)));
   else if (SCM_COMPLEXP (z))
     { double x, y;
@@ -6843,7 +6853,9 @@ SCM_PRIMITIVE_GENERIC (scm_tanh, "tanh", 1, 0, 0,
                        "Compute the hyperbolic tangent of @var{z}.")
 #define FUNC_NAME s_scm_tanh
 {
-  if (scm_is_real (z))
+  if (SCM_UNLIKELY (scm_is_eq (z, SCM_INUM0)))
+    return z;  /* tanh(exact0) = exact0 */
+  else if (scm_is_real (z))
     return scm_from_double (tanh (scm_to_double (z)));
   else if (SCM_COMPLEXP (z))
     { double x, y, w;
@@ -6866,7 +6878,9 @@ SCM_PRIMITIVE_GENERIC (scm_asin, "asin", 1, 0, 0,
                        "Compute the arc sine of @var{z}.")
 #define FUNC_NAME s_scm_asin
 {
-  if (scm_is_real (z))
+  if (SCM_UNLIKELY (scm_is_eq (z, SCM_INUM0)))
+    return z;  /* asin(exact0) = exact0 */
+  else if (scm_is_real (z))
     {
       double w = scm_to_double (z);
       if (w >= -1.0 && w <= 1.0)
@@ -6892,7 +6906,9 @@ SCM_PRIMITIVE_GENERIC (scm_acos, "acos", 1, 0, 0,
                        "Compute the arc cosine of @var{z}.")
 #define FUNC_NAME s_scm_acos
 {
-  if (scm_is_real (z))
+  if (SCM_UNLIKELY (scm_is_eq (z, SCM_INUM1)))
+    return SCM_INUM0;  /* acos(exact1) = exact0 */
+  else if (scm_is_real (z))
     {
       double w = scm_to_double (z);
       if (w >= -1.0 && w <= 1.0)
@@ -6924,7 +6940,9 @@ SCM_PRIMITIVE_GENERIC (scm_atan, "atan", 1, 1, 0,
 {
   if (SCM_UNBNDP (y))
     {
-      if (scm_is_real (z))
+      if (SCM_UNLIKELY (scm_is_eq (z, SCM_INUM0)))
+	return z;  /* atan(exact0) = exact0 */
+      else if (scm_is_real (z))
         return scm_from_double (atan (scm_to_double (z)));
       else if (SCM_COMPLEXP (z))
         {
@@ -6955,7 +6973,9 @@ SCM_PRIMITIVE_GENERIC (scm_sys_asinh, "asinh", 1, 0, 0,
                        "Compute the inverse hyperbolic sine of @var{z}.")
 #define FUNC_NAME s_scm_sys_asinh
 {
-  if (scm_is_real (z))
+  if (SCM_UNLIKELY (scm_is_eq (z, SCM_INUM0)))
+    return z;  /* asinh(exact0) = exact0 */
+  else if (scm_is_real (z))
     return scm_from_double (asinh (scm_to_double (z)));
   else if (scm_is_number (z))
     return scm_log (scm_sum (z,
@@ -6971,7 +6991,9 @@ SCM_PRIMITIVE_GENERIC (scm_sys_acosh, "acosh", 1, 0, 0,
                        "Compute the inverse hyperbolic cosine of @var{z}.")
 #define FUNC_NAME s_scm_sys_acosh
 {
-  if (scm_is_real (z) && scm_to_double (z) >= 1.0)
+  if (SCM_UNLIKELY (scm_is_eq (z, SCM_INUM1)))
+    return SCM_INUM0;  /* acosh(exact1) = exact0 */
+  else if (scm_is_real (z) && scm_to_double (z) >= 1.0)
     return scm_from_double (acosh (scm_to_double (z)));
   else if (scm_is_number (z))
     return scm_log (scm_sum (z,
@@ -6987,7 +7009,9 @@ SCM_PRIMITIVE_GENERIC (scm_sys_atanh, "atanh", 1, 0, 0,
                        "Compute the inverse hyperbolic tangent of @var{z}.")
 #define FUNC_NAME s_scm_sys_atanh
 {
-  if (scm_is_real (z) && scm_to_double (z) >= -1.0 && scm_to_double (z) <= 1.0)
+  if (SCM_UNLIKELY (scm_is_eq (z, SCM_INUM0)))
+    return z;  /* atanh(exact0) = exact0 */
+  else if (scm_is_real (z) && scm_to_double (z) >= -1.0 && scm_to_double (z) <= 1.0)
     return scm_from_double (atanh (scm_to_double (z)));
   else if (scm_is_number (z))
     return scm_divide (scm_log (scm_divide (scm_sum (SCM_INUM1, z),
diff --git a/test-suite/tests/numbers.test b/test-suite/tests/numbers.test
index 96fb6d9..9c01fa1 100644
--- a/test-suite/tests/numbers.test
+++ b/test-suite/tests/numbers.test
@@ -3316,25 +3316,121 @@
 
 
 ;;;
+;;; sin
+;;;
+
+(with-test-prefix "sin"
+  (pass-if (eqv? 0   (sin 0)))
+  (pass-if (eqv? 0.0 (sin 0.0)))
+  (pass-if (eqv-loosely? 1.0 (sin 1.57)))
+  (pass-if (eqv-loosely? +1.175i (sin +i)))
+  (pass-if (real-nan? (sin +nan.0)))
+  (pass-if (real-nan? (sin +inf.0)))
+  (pass-if (real-nan? (sin -inf.0))))
+
+;;;
+;;; cos
+;;;
+
+(with-test-prefix "cos"
+  (pass-if (eqv? 1   (cos 0)))
+  (pass-if (eqv? 1.0 (cos 0.0)))
+  (pass-if (eqv-loosely? 0.0 (cos 1.57)))
+  (pass-if (eqv-loosely? 1.543 (cos +i)))
+  (pass-if (real-nan? (cos +nan.0)))
+  (pass-if (real-nan? (cos +inf.0)))
+  (pass-if (real-nan? (cos -inf.0))))
+
+;;;
+;;; tan
+;;;
+
+(with-test-prefix "tan"
+  (pass-if (eqv? 0   (tan 0)))
+  (pass-if (eqv? 0.0 (tan 0.0)))
+  (pass-if (eqv-loosely? 1.0 (tan 0.785)))
+  (pass-if (eqv-loosely? +0.76i (tan +i)))
+  (pass-if (real-nan? (tan +nan.0)))
+  (pass-if (real-nan? (tan +inf.0)))
+  (pass-if (real-nan? (tan -inf.0))))
+
+;;;
+;;; asin
+;;;
+
+(with-test-prefix "asin"
+  (pass-if (complex-nan? (asin +nan.0)))
+  (pass-if (eqv? 0 (asin 0)))
+  (pass-if (eqv? 0.0 (asin 0.0))))
+
+;;;
+;;; acos
+;;;
+
+(with-test-prefix "acos"
+  (pass-if (complex-nan? (acos +nan.0)))
+  (pass-if (eqv? 0 (acos 1)))
+  (pass-if (eqv? 0.0 (acos 1.0))))
+
+;;;
+;;; atan
+;;;
+;;; FIXME: add tests for two-argument atan
+;;;
+(with-test-prefix "atan"
+  (pass-if (real-nan? (atan +nan.0)))
+  (pass-if (eqv? 0 (atan 0)))
+  (pass-if (eqv? 0.0 (atan 0.0)))
+  (pass-if (eqv-loosely?  1.57 (atan +inf.0)))
+  (pass-if (eqv-loosely? -1.57 (atan -inf.0))))
+
+;;;
+;;; sinh
+;;;
+
+(with-test-prefix "sinh"
+  (pass-if (= 0 (sinh 0)))
+  (pass-if (= 0.0 (sinh 0.0))))
+
+;;;
+;;; cosh
+;;;
+
+(with-test-prefix "cosh"
+  (pass-if (= 1 (cosh 0)))
+  (pass-if (= 1.0 (cosh 0.0))))
+
+;;;
+;;; tanh
+;;;
+
+(with-test-prefix "tanh"
+  (pass-if (= 0 (tanh 0)))
+  (pass-if (= 0.0 (tanh 0.0))))
+
+;;;
 ;;; asinh
 ;;;
 
 (with-test-prefix "asinh"
-  (pass-if (= 0 (asinh 0))))
+  (pass-if (= 0 (asinh 0)))
+  (pass-if (= 0.0 (asinh 0.0))))
 
 ;;;
 ;;; acosh
 ;;;
 
 (with-test-prefix "acosh"
-  (pass-if (= 0 (acosh 1))))
+  (pass-if (= 0 (acosh 1)))
+  (pass-if (= 0.0 (acosh 1.0))))
 
 ;;;
 ;;; atanh
 ;;;
 
 (with-test-prefix "atanh"
-  (pass-if (= 0 (atanh 0))))
+  (pass-if (= 0 (atanh 0)))
+  (pass-if (= 0.0 (atanh 0.0))))
 
 ;;;
 ;;; make-rectangular
-- 
1.5.6.5


[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #7: Improve discussion of exactness propagation in manual --]
[-- Type: text/x-diff, Size: 1991 bytes --]

From c00d87d077720a895ef8f52732760549e15c3b6d Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Thu, 27 Jan 2011 15:57:38 -0500
Subject: [PATCH] Improve discussion of exactness propagation in manual

* doc/ref/api-data.texi (Exact and Inexact Numbers): Improve the
  discussion of exactness propagation.  Mention that there are
  exceptions to the rule that calculations involving inexact numbers
  must product an inexact result.
---
 doc/ref/api-data.texi |   13 +++++++++----
 1 files changed, 9 insertions(+), 4 deletions(-)

diff --git a/doc/ref/api-data.texi b/doc/ref/api-data.texi
index b819fcb..1ce9e1e 100755
--- a/doc/ref/api-data.texi
+++ b/doc/ref/api-data.texi
@@ -712,14 +712,19 @@ Equivalent to @code{scm_is_true (scm_complex_p (val))}.
 @rnindex exact->inexact
 @rnindex inexact->exact
 
-R5RS requires that a calculation involving inexact numbers always
-produces an inexact result.  To meet this requirement, Guile
-distinguishes between an exact integer value such as @samp{5} and the
-corresponding inexact real value which, to the limited precision
+R5RS requires that, with few exceptions, a calculation involving inexact
+numbers always produces an inexact result.  To meet this requirement,
+Guile distinguishes between an exact integer value such as @samp{5} and
+the corresponding inexact integer value which, to the limited precision
 available, has no fractional part, and is printed as @samp{5.0}.  Guile
 will only convert the latter value to the former when forced to do so by
 an invocation of the @code{inexact->exact} procedure.
 
+The only exception to the above requirement is when the values of the
+inexact numbers do not affect the result.  For example @code{(expt n 0)}
+is @samp{1} for any value of @code{n}, therefore @code{(expt 5.0 0)} is
+permitted to return an exact @samp{1}.
+
 @deffn {Scheme Procedure} exact? z
 @deffnx {C Function} scm_exact_p (z)
 Return @code{#t} if the number @var{z} is exact, @code{#f}
-- 
1.5.6.5

Re: [PATCH] Handle products with exact 0 differently, etc

[Andy Wingo]

Hi Mark,

On Tue 01 Feb 2011 13:09, Mark H Weaver <mhw@netris.org> writes:

> Subject: [PATCH] Trigonometric functions return exact numbers in some
> cases

Why would you want this?

Just wondering.

Andy
-- 
http://wingolog.org/

Re: [PATCH] Handle products with exact 0 differently, etc

[Mark H Weaver]

>> Subject: [PATCH] Trigonometric functions return exact numbers in some
>> cases
>
> Why would you want this?

I'm a developer of Maxima, a free CAS (computer algebra system).  I'm
interested in building a free CAS on top of Guile, probably with large
chunks of functionality adapted from Maxima.  As part of this, I would
like it to return exact results whenever feasible.  If the argument
passed to a trig function is a symbolic expression, that will
automatically be handled by GOOPS methods for the various operators.
However, there are a few exceptions where trig functions applied to
simple numbers will never reach GOOPS.  I'd like to do the right thing
in those cases.

Note that this is also the motivation behind my interest in bigfloats,
and in making sure that the numeric code is fully extensible.

     Best,
      Mark

Re: [PATCH] Handle products with exact 0 differently, etc

[Andy Wingo]

On Tue 01 Feb 2011 13:09, Mark H Weaver <mhw@netris.org> writes:

> Here's another batch of numerics patches.

Thanks, applied all but the trigonometric functions one, pending further
discussion.

Andy
-- 
http://wingolog.org/

Re: [PATCH] Handle products with exact 0 differently, etc

[Neil Jerram]

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

> Here's another batch of numerics patches.  The most important one
> changes the way products involving exact 0 are handled:
>
> * libguile/numbers.c (scm_product): Handle exact 0 differently.  A
>   product containing an exact 0 now returns an exact 0 if and only if
>   the other arguments are all exact.

I don't get this one.  I would expect "... are all finite".  Why would
an exact 0 times a finite number be inexact?

>  If an
>   infinite or NaN value is present, a NaN value is returned.

No problem here - this one makes sense.

>   A note on the rationale for (* 0 0.0) returning 0.0 and not exact 0:
>   The exactness propagation rules allow us to return an exact result in
>   the presence of inexact arguments only if the values of the inexact
>   arguments do not affect the result.

OK.  And in the case where all the inexact arguments are finite, I would
say that "the values of the inexact arguments do not affect the result"
is true.

>  In this case, the value of the
>   inexact argument _does_ affect the result, because an infinite or NaN
>   value causes the result to be a NaN.

But in the (* 0 0.0) case, the inexact argument is not an infinite or
NaN value.  So I don't understand the logic here.

Regards,
        Neil

Re: [PATCH] Handle products with exact 0 differently, etc

[Mark H Weaver]

Neil Jerram <neil@ossau.uklinux.net> writes:
>> * libguile/numbers.c (scm_product): Handle exact 0 differently.  A
>>   product containing an exact 0 now returns an exact 0 if and only if
>>   the other arguments are all exact.
>
> I don't get this one.  I would expect "... are all finite".

I almost wrote "... are all finite and exact", but exact implies finite,
so I simplified it.

> Why would an exact 0 times a finite number be inexact?

The point of the exact/inexact distinction is to ensure that if you see
an exact number, you can trust that it is provably the correct answer,
and that no inaccuracies associated with inexact arithmetic along the
way could possibly have corrupted that answer.

In this case, the inaccuracies associated with inexact arithmetic could
result in an infinity being misrepresented as a finite number, or vice
versa.  For example, if X is inexact, then we cannot claim that the
result of (* 0 (/ X)) is an exact 0, because a small change in X could
affect the final answer.

     Best,
      Mark

Re: [PATCH] Handle products with exact 0 differently, etc

[Neil Jerram]

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

> The point of the exact/inexact distinction is to ensure that if you see
> an exact number, you can trust that it is provably the correct answer,
> and that no inaccuracies associated with inexact arithmetic along the
> way could possibly have corrupted that answer.

Yes, I think I understand that.

> In this case, the inaccuracies associated with inexact arithmetic could
> result in an infinity being misrepresented as a finite number, or vice
> versa.  For example, if X is inexact, then we cannot claim that the
> result of (* 0 (/ X)) is an exact 0, because a small change in X could
> affect the final answer.

In other words, the argument is that any inexact number might actually
be infinite.  (Right?)

That strikes me as stretching the idea of inexactness too far; and also
as not useful, because I'm sure there are many practical ways of
producing inexact numbers that are provably finite (e.g. sqrt(5)).  My
guess is that for the majority of programming situations where I have an
inexact number X, that number is provably finite, and so I'd want (* 0
X) to be an exact 0.

But I'm not an expert, and I'm not actually writing those programs
now...

Are you following specific references here?  FWIW, for cases like this,
the last para of R6RS 11.7.1 allows "either 0 (exact) or 0.0 (inexact)".

        Neil

Re: [PATCH] Handle products with exact 0 differently, etc

[Mark H Weaver]

Neil Jerram <neil@ossau.uklinux.net> writes:
>> In this case, the inaccuracies associated with inexact arithmetic could
>> result in an infinity being misrepresented as a finite number, or vice
>> versa.  For example, if X is inexact, then we cannot claim that the
>> result of (* 0 (/ X)) is an exact 0, because a small change in X could
>> affect the final answer.
>
> In other words, the argument is that any inexact number might actually
> be infinite.  (Right?)

Yes.

> That strikes me as stretching the idea of inexactness too far; and also
> as not useful, because I'm sure there are many practical ways of
> producing inexact numbers that are provably finite (e.g. sqrt(5)).

There are some cases where it can be proven, but anytime you divide by
an inexact number (or apply any other function with singularities), you
must assume that the finiteness of the result is unknown.  Unless you
keep track of the necessary information about how an inexact number was
computed, you cannot know its finiteness, therefore you must assume the
worst.

One approach to proving finiteness in simple cases might be as follows:
Every inexact number could include an additional flag indicating whether
it's provably finite.  Division by an inexact number would have to
return an inexact number with that flag cleared.  Similarly for other
functions with singularities.  Most other operations would AND the
operand flags together to determine the result flag.

However, I don't currently see that the benefits of this would outweigh
the additional implementation complexity and overhead.

> My guess is that for the majority of programming situations where I
> have an inexact number X, that number is provably finite, and so I'd
> want (* 0 X) to be an exact 0.

Can you give a practical example where it is important for (* 0 X) for
inexact X to produce an exact 0?  If I understood why this was important
(besides the aesthetics of an exact result), I might be more motivated
to try to produce an exact 0 when possible.

> Are you following specific references here?  FWIW, for cases like this,
> the last para of R6RS 11.7.1 allows "either 0 (exact) or 0.0 (inexact)".

Yes, the R6RS does allow (* 1.0 0) to be an exact 0, but that's because
it also allows (* 0 +inf.0) and (* 0 +nan.0) to be an exact 0.  Section
11.7.4.3 provides the following examples:

  (* 0 +inf.0)  ==>  0 or +nan.0
  (* 0 +nan.0)  ==>  0 or +nan.0
  (* 1.0 0)     ==>  0 or 0.0

It follows from the rules of exactness propagation that (* 1.0 0) may be
an exact 0 only if (* 0 X) is an exact 0 for _any_ inexact X.  We
_could_ do this, if we made (* 0 +inf.0) and (* 0 +nan.0) return exact 0
as well.  However, doing this would necessarily destroy the equivalence
of (/ X Y) and (* X (/ Y)), as explained in my first message of this
thread.

The following note in section 11.7.4.3, regarding `min' and `max', gives
a hint of how strictly these exactness rules should be enforced:

  Note: If any argument is inexact, then the result is also inexact
  (unless the procedure can prove that the inaccuracy is not large
  enough to affect the result, which is possible only in unusual
  implementations).  If min or max is used to compare number objects of
  mixed exactness, and the numerical value of the result cannot be
  represented as an inexact number object without loss of accuracy, then
  the procedure may raise an exception with condition type
  &implementation-restriction.

     Best,
      Mark

Re: [PATCH] Handle products with exact 0 differently, etc

[Mark H Weaver]

I wrote:
>> Are you following specific references here?  FWIW, for cases like this,
>> the last para of R6RS 11.7.1 allows "either 0 (exact) or 0.0 (inexact)".
>
> Yes, the R6RS does allow (* 1.0 0) to be an exact 0, but that's because
> it also allows (* 0 +inf.0) and (* 0 +nan.0) to be an exact 0.  Section
> 11.7.4.3 provides the following examples:
>
>   (* 0 +inf.0)  ==>  0 or +nan.0
>   (* 0 +nan.0)  ==>  0 or +nan.0
>   (* 1.0 0)     ==>  0 or 0.0
>
> It follows from the rules of exactness propagation that (* 1.0 0) may be
> an exact 0 only if (* 0 X) is an exact 0 for _any_ inexact X.

I apologize, I should have admitted that I am going a bit further than
the R6RS strictly requires, based on my understanding of the rationale
behind exactness in Scheme.

It is true that both the R5RS and the R6RS specifically give (* 0 X),
for inexact X, as an example where it is permissible to return an exact
0, without mentioning any associated caveats, e.g. that one must also
make (* 0 +inf.0) and (* 0 +nan.0) be an exact 0 as well.

I suspect the reason that the R5RS did not mention these caveats is that
it does not specify anything about infinities or NaNs.  For example,
MIT/GNU Scheme evaluates (* 0 X) to exact 0 for inexact X, but that is
justified because its arithmetic operations do not produce infinities or
NaNs; they throw exceptions instead.

As for the R6RS, I confess that a strict reading of the exactness
propagation rule does allow (* 0 X) to be 0 for inexact 0, even
(* 0 +inf.0) yields something else.  That's because their rule requires
only that the same result is returned for all possible substitutions of
_exact_ arguments for the inexact ones.  Since the infinities and NaNs
are not exact, they are not considered as possible substitutions.

  One general exception to the rule above is that an implementation may
  return an exact result despite inexact arguments if that exact result
  would be the correct result for all possible substitutions of exact
  arguments for the inexact ones.

I think that the R6RS is wrong here.  As I understand it, the rationale
for the exactness system is to allow the construction of provably
correct programs, despite the existence of inexact arithmetic.  One
should be able to trust that any exact number is provably correct,
presuming that inexact->exact was not used of course.

Is there a compelling reason why (* 0 X) should yield an exact 0 for
inexact X?  Is there a reason more important than being able to trust
that exact values are provably correct?

    Best,
     Mark

Re: [PATCH] Handle products with exact 0 differently, etc

[Noah Lavine]

Hello Mark,

I haven't read through all of the discussion yet, but it's obvious
that you have good reasons for wanting (* 0 X) to be NaN when X is
inexact. And yet for compatibility reasons it is nice if Guile agrees
with Scheme standards. Therefore I think it would be great if you
would send an email with exactly what you said here to
scheme-reports@scheme-reports.org, which is the public discussion
forum for R7RS.

Of course Guile could just implement this itself, and maybe it will,
but our lives would be much easier if our behavior could be both
useful and standard, and it seems that the best way to do that might
be to change the standard.

Noah

On Tue, Feb 1, 2011 at 11:22 PM, Mark H Weaver <mhw@netris.org> wrote:
> I wrote:
>>> Are you following specific references here?  FWIW, for cases like this,
>>> the last para of R6RS 11.7.1 allows "either 0 (exact) or 0.0 (inexact)".
>>
>> Yes, the R6RS does allow (* 1.0 0) to be an exact 0, but that's because
>> it also allows (* 0 +inf.0) and (* 0 +nan.0) to be an exact 0.  Section
>> 11.7.4.3 provides the following examples:
>>
>>   (* 0 +inf.0)  ==>  0 or +nan.0
>>   (* 0 +nan.0)  ==>  0 or +nan.0
>>   (* 1.0 0)     ==>  0 or 0.0
>>
>> It follows from the rules of exactness propagation that (* 1.0 0) may be
>> an exact 0 only if (* 0 X) is an exact 0 for _any_ inexact X.
>
> I apologize, I should have admitted that I am going a bit further than
> the R6RS strictly requires, based on my understanding of the rationale
> behind exactness in Scheme.
>
> It is true that both the R5RS and the R6RS specifically give (* 0 X),
> for inexact X, as an example where it is permissible to return an exact
> 0, without mentioning any associated caveats, e.g. that one must also
> make (* 0 +inf.0) and (* 0 +nan.0) be an exact 0 as well.
>
> I suspect the reason that the R5RS did not mention these caveats is that
> it does not specify anything about infinities or NaNs.  For example,
> MIT/GNU Scheme evaluates (* 0 X) to exact 0 for inexact X, but that is
> justified because its arithmetic operations do not produce infinities or
> NaNs; they throw exceptions instead.
>
> As for the R6RS, I confess that a strict reading of the exactness
> propagation rule does allow (* 0 X) to be 0 for inexact 0, even
> (* 0 +inf.0) yields something else.  That's because their rule requires
> only that the same result is returned for all possible substitutions of
> _exact_ arguments for the inexact ones.  Since the infinities and NaNs
> are not exact, they are not considered as possible substitutions.
>
>  One general exception to the rule above is that an implementation may
>  return an exact result despite inexact arguments if that exact result
>  would be the correct result for all possible substitutions of exact
>  arguments for the inexact ones.
>
> I think that the R6RS is wrong here.  As I understand it, the rationale
> for the exactness system is to allow the construction of provably
> correct programs, despite the existence of inexact arithmetic.  One
> should be able to trust that any exact number is provably correct,
> presuming that inexact->exact was not used of course.
>
> Is there a compelling reason why (* 0 X) should yield an exact 0 for
> inexact X?  Is there a reason more important than being able to trust
> that exact values are provably correct?
>
>    Best,
>     Mark
>
>

Re: [PATCH] Handle products with exact 0 differently, etc

[Mark H Weaver]

Noah Lavine <noah.b.lavine@gmail.com> writes:
> I haven't read through all of the discussion yet, but it's obvious
> that you have good reasons for wanting (* 0 X) to be NaN when X is
> inexact. And yet for compatibility reasons it is nice if Guile agrees
> with Scheme standards.

In case there is any doubt, the behavior of (* 0 X) that I am advocating
(that it should yield an inexact result when X is inexact) is clearly
permitted by both the R5RS and the R6RS.  That much is beyond any doubt.

My disagreement with the R6RS is that it _permits_ an exact 0 result in
this case, even when infinities and NaNs are supported as numeric
objects.

> Therefore I think it would be great if you would send an email with
> exactly what you said here to scheme-reports@scheme-reports.org, which
> is the public discussion forum for R7RS.

Yes, thank you for suggesting that, I should write something up.

      Mark

Re: [PATCH] Handle products with exact 0 differently, etc

[Neil Jerram]

Noah Lavine <noah.b.lavine@gmail.com> writes:

> Hello Mark,
>
> I haven't read through all of the discussion yet, but it's obvious
> that you have good reasons for wanting (* 0 X) to be NaN when X is

NB (and for the record), "NaN" here should be "inexact".

Re: [PATCH] Handle products with exact 0 differently, etc

[Neil Jerram]

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

> Neil Jerram <neil@ossau.uklinux.net> writes:
>>
>> In other words, the argument is that any inexact number might actually
>> be infinite.  (Right?)
>
> Yes.
>
>> That strikes me as stretching the idea of inexactness too far; and also
>> as not useful, because I'm sure there are many practical ways of
>> producing inexact numbers that are provably finite (e.g. sqrt(5)).
>
> There are some cases where it can be proven, but anytime you divide by
> an inexact number (or apply any other function with singularities), you
> must assume that the finiteness of the result is unknown.

Again, that seems impractical to me.  Just as many inexact number
calculations are provably finite, there are many inexact number
calculations that are provably non-zero.  This approach seems like
letting the tail wag the dog.

However, as you've indicated by your question below, what matters for
our purposes is whether it makes a difference to programming.  So I'll
put my intuition aside and focus on that.

>> My guess is that for the majority of programming situations where I
>> have an inexact number X, that number is provably finite, and so I'd
>> want (* 0 X) to be an exact 0.
>
> Can you give a practical example where it is important for (* 0 X) for
> inexact X to produce an exact 0?  If I understood why this was important
> (besides the aesthetics of an exact result), I might be more motivated
> to try to produce an exact 0 when possible.

No, I can't.

If a program knew for sure that a multiplication included an exact 0,
and that the other arguments were all finite, there would be no point
having the multiplication there.  In other cases, if a program needed to
rely on the result of a multiplication being exact, it would either have
to know that all the arguments were exact, or use inexact->exact.

I think that leaves no cases where it makes a practical difference
whether (* 0 X) is exact or inexact.  So on that basis I'm happy to
concede this point.

Many thanks for your patient explanations, and for all your work on
these numeric patches.

Regards,
        Neil