From: Mark H Weaver <mhw@netris.org>
To: Andy Wingo <wingo@pobox.com>
Cc: guile-devel@gnu.org
Subject: Re: [PATCH] First batch of numerics changes
Date: Sat, 29 Jan 2011 03:20:40 -0500 [thread overview]
Message-ID: <87wrlo2k9z.fsf@yeeloong.netris.org> (raw)
In-Reply-To: <m3ipx9dttf.fsf@unquote.localdomain> (Andy Wingo's message of "Fri, 28 Jan 2011 14:46:52 +0100")
[-- Attachment #1: Type: text/plain, Size: 371 bytes --]
Andy Wingo <wingo@pobox.com> writes:
> I applied up to patch 12, I think, [...]
Great, thanks! Here's another batch which should now apply cleanly.
This includes an efficient implementation of the R6RS division
operations along with documentation and extensive tests. The code has
many paths, so I was careful to include tests for all of them.
Best,
Mark
[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #2: Remove useless test and fix spelling errors --]
[-- Type: text/x-diff, Size: 1887 bytes --]
From 457f9ce87af2e15438662eb4ec4caf7b7a4aa4d1 Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Fri, 28 Jan 2011 19:13:47 -0500
Subject: [PATCH] Remove useless test and fix spelling errors
* test-suite/tests/numbers.test: Remove test for lazy reduction bit of
fractions, which was never implemented. Fix some spelling errors.
---
test-suite/tests/numbers.test | 20 ++++----------------
1 files changed, 4 insertions(+), 16 deletions(-)
diff --git a/test-suite/tests/numbers.test b/test-suite/tests/numbers.test
index f53cb34..4f30f6c 100644
--- a/test-suite/tests/numbers.test
+++ b/test-suite/tests/numbers.test
@@ -318,15 +318,15 @@
(pass-if (not (finite? +inf.0)))
(pass-if (not (finite? -inf.0)))
(pass-if-exception
- "complex numbers not in doman of finite?"
+ "complex numbers not in domain of finite?"
exception:wrong-type-arg
(finite? +inf.0+1i))
(pass-if-exception
- "complex numbers not in doman of finite? (2)"
+ "complex numbers not in domain of finite? (2)"
exception:wrong-type-arg
(finite? +1+inf.0i))
(pass-if-exception
- "complex numbers not in doman of finite? (3)"
+ "complex numbers not in domain of finite? (3)"
exception:wrong-type-arg
(finite? +1+1i))
(pass-if (finite? 3+0i))
@@ -351,7 +351,7 @@
;; (pass-if (inf? (/ 1.0 0.0))
;; (pass-if (inf? (/ 1 0.0))
(pass-if-exception
- "complex numbers not in doman of inf?"
+ "complex numbers not in domain of inf?"
exception:wrong-type-arg
(inf? +1+inf.0i))
(pass-if (inf? +inf.0+0i))
@@ -3386,15 +3386,3 @@
(pass-if "-100i swings back to 45deg down"
(eqv-loosely? +7.071-7.071i (sqrt -100.0i))))
-
-;;
-;; equal?
-;;
-
-
-(with-test-prefix "equal?"
- (pass-if
-
- ;; lazy reduction bit for rationals should not affect equal?
- (equal? 1/2 ((lambda (x) (denominator x) x) 1/2))))
-
--
1.5.6.5
[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #3: `equal?' and `eqv?' are now equivalent for numbers --]
[-- Type: text/x-diff, Size: 12665 bytes --]
From 3afeb53165c69f95120336a5b6cbb83a810be1e9 Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Fri, 28 Jan 2011 19:57:41 -0500
Subject: [PATCH] `equal?' and `eqv?' are now equivalent for numbers
Change `equal?' to work like `eqv?' for numbers.
Previously they worked differently in some cases, e.g.
when comparing signed zeroes or NaNs. For example,
(equal? 0.0 -0.0) returned #t but (eqv? 0.0 -0.0)
returned #f, and (equal? +nan.0 +nan.0) returned #f
but (eqv? +nan.0 +nan.0) returned #t.
* libguile/numbers.c (scm_real_equalp, scm_bigequal,
scm_complex_equalp, scm_i_fraction_equalp): Move to eq.c.
* libguile/eq.c (scm_real_equalp): Compare flonums using
real_eqv instead of ==, so that NaNs are now considered
equal, and to distinguish signed zeroes.
(scm_complex_equalp): Compare real and imaginary
components using real_eqv instead of ==, so that NaNs are
now considered equal, and to distinguish signed zeroes.
(scm_bigequal): Use scm_i_bigcmp instead of duplicating it.
(real_eqv): Test for NaNs using isnan(x) instead of
(x != x), and use SCM_UNLIKELY for optimization.
(scm_eqv_p): Use scm_bigequal, scm_real_equalp,
scm_complex_equalp, and scm_i_fraction_equalp to compare
numbers, instead of inline code. Those predicates now do
what scm_eqv_p formerly did internally. Replace if
statements with switch statements, as is done in
scm_equal_p. Remove useless code to check equality of
fractions with different SCM_CELL_TYPEs; this was for a
tentative "lazy reduction bit" which was never developed.
(scm_eqv_p, scm_equal_p): Remove useless code to check
equality between inexact reals and non-real complex numbers
with zero imaginary part. Such numbers do not exist,
because the current code is careful to never create them.
* test-suite/tests/numbers.test: Add test cases for
`eqv?' and `equal?'. Change existing test case for
`(equal? +nan.0 +nan.0)' to expect #t instead of #f.
* NEWS: Add NEWS entries.
---
NEWS | 15 ++++++
libguile/eq.c | 106 ++++++++++++++++++++---------------------
libguile/numbers.c | 34 -------------
test-suite/tests/numbers.test | 86 +++++++++++++++++++++++++++++++++-
4 files changed, 152 insertions(+), 89 deletions(-)
diff --git a/NEWS b/NEWS
index 9938204..2979849 100644
--- a/NEWS
+++ b/NEWS
@@ -12,6 +12,21 @@ Changes in 1.9.15 (since the 1.9.14 prerelease):
** Changes and bugfixes in numerics code
+*** `eqv?' and `equal?' now compare numbers equivalently
+
+scm_equal_p `equal?' now behaves equivalently to scm_eqv_p `eqv?' for
+numeric values, per R5RS. Previously, equal? worked differently,
+e.g. `(equal? 0.0 -0.0)' returned #t but `(eqv? 0.0 -0.0)' returned #f,
+and `(equal? +nan.0 +nan.0)' returned #f but `(eqv? +nan.0 +nan.0)'
+returned #t.
+
+*** `(equal? +nan.0 +nan.0)' now returns #t
+
+Previously, `(equal? +nan.0 +nan.0)' returned #f, although
+`(let ((x +nan.0)) (equal? x x))' and `(eqv? +nan.0 +nan.0)'
+both returned #t. R5RS requires that `equal?' behave like
+`eqv?' when comparing numbers.
+
*** Infinities are no longer integers.
Following the R6RS, infinities (+inf.0 and -inf.0) are no longer
diff --git a/libguile/eq.c b/libguile/eq.c
index 7502559..00abdd8 100644
--- a/libguile/eq.c
+++ b/libguile/eq.c
@@ -1,4 +1,4 @@
-/* Copyright (C) 1995,1996,1997,1998,2000,2001,2003, 2004, 2006, 2009, 2010 Free Software Foundation, Inc.
+/* Copyright (C) 1995,1996,1997,1998,2000,2001,2003, 2004, 2006, 2009, 2010, 2011 Free Software Foundation, Inc.
*
* This library is free software; you can redistribute it and/or
* modify it under the terms of the GNU Lesser General Public License
@@ -118,7 +118,40 @@ scm_eq_p (SCM x, SCM y)
static int
real_eqv (double x, double y)
{
- return !memcmp (&x, &y, sizeof(double)) || (x != x && y != y);
+ return !memcmp (&x, &y, sizeof(double))
+ || (SCM_UNLIKELY (isnan (x)) && SCM_UNLIKELY (isnan (y)));
+}
+
+SCM
+scm_real_equalp (SCM x, SCM y)
+{
+ return scm_from_bool (real_eqv (SCM_REAL_VALUE (x),
+ SCM_REAL_VALUE (y)));
+}
+
+SCM
+scm_bigequal (SCM x, SCM y)
+{
+ return scm_from_bool (scm_i_bigcmp (x, y) == 0);
+}
+
+SCM
+scm_complex_equalp (SCM x, SCM y)
+{
+ return scm_from_bool (real_eqv (SCM_COMPLEX_REAL (x),
+ SCM_COMPLEX_REAL (y))
+ && real_eqv (SCM_COMPLEX_IMAG (x),
+ SCM_COMPLEX_IMAG (y)));
+}
+
+SCM
+scm_i_fraction_equalp (SCM x, SCM y)
+{
+ return scm_from_bool
+ (scm_is_true (scm_equal_p (SCM_FRACTION_NUMERATOR (x),
+ SCM_FRACTION_NUMERATOR (y)))
+ && scm_is_true (scm_equal_p (SCM_FRACTION_DENOMINATOR (x),
+ SCM_FRACTION_DENOMINATOR (y))));
}
static SCM scm_i_eqv_p (SCM x, SCM y, SCM rest);
@@ -166,48 +199,26 @@ SCM scm_eqv_p (SCM x, SCM y)
return SCM_BOOL_F;
if (SCM_IMP (y))
return SCM_BOOL_F;
- /* this ensures that types and scm_length are the same. */
+ /* this ensures that types and scm_length are the same. */
if (SCM_CELL_TYPE (x) != SCM_CELL_TYPE (y))
+ return SCM_BOOL_F;
+ switch (SCM_TYP7 (x))
{
- /* fractions use 0x10000 as a flag (at the suggestion of Marius Vollmer),
- but this checks the entire type word, so fractions may be accidentally
- flagged here as unequal. Perhaps I should use the 4th double_cell word?
- */
-
- /* treat mixes of real and complex types specially */
- if (SCM_INEXACTP (x))
- {
- if (SCM_REALP (x))
- return scm_from_bool (SCM_COMPLEXP (y)
- && real_eqv (SCM_REAL_VALUE (x),
- SCM_COMPLEX_REAL (y))
- && SCM_COMPLEX_IMAG (y) == 0.0);
- else
- return scm_from_bool (SCM_REALP (y)
- && real_eqv (SCM_COMPLEX_REAL (x),
- SCM_REAL_VALUE (y))
- && SCM_COMPLEX_IMAG (x) == 0.0);
- }
-
- if (SCM_FRACTIONP (x) && SCM_FRACTIONP (y))
- return scm_i_fraction_equalp (x, y);
- return SCM_BOOL_F;
- }
- if (SCM_NUMP (x))
- {
- if (SCM_BIGP (x)) {
- return scm_from_bool (scm_i_bigcmp (x, y) == 0);
- } else if (SCM_REALP (x)) {
- return scm_from_bool (real_eqv (SCM_REAL_VALUE (x), SCM_REAL_VALUE (y)));
- } else if (SCM_FRACTIONP (x)) {
- return scm_i_fraction_equalp (x, y);
- } else { /* complex */
- return scm_from_bool (real_eqv (SCM_COMPLEX_REAL (x),
- SCM_COMPLEX_REAL (y))
- && real_eqv (SCM_COMPLEX_IMAG (x),
- SCM_COMPLEX_IMAG (y)));
- }
+ default:
+ break;
+ case scm_tc7_number:
+ switch SCM_TYP16 (x)
+ {
+ case scm_tc16_big:
+ return scm_bigequal (x, y);
+ case scm_tc16_real:
+ return scm_real_equalp (x, y);
+ case scm_tc16_complex:
+ return scm_complex_equalp (x, y);
+ case scm_tc16_fraction:
+ return scm_i_fraction_equalp (x, y);
+ }
}
return SCM_BOOL_F;
}
@@ -309,19 +320,6 @@ scm_equal_p (SCM x, SCM y)
/* This ensures that types and scm_length are the same. */
if (SCM_CELL_TYPE (x) != SCM_CELL_TYPE (y))
{
- /* treat mixes of real and complex types specially */
- if (SCM_INEXACTP (x) && SCM_INEXACTP (y))
- {
- if (SCM_REALP (x))
- return scm_from_bool (SCM_COMPLEXP (y)
- && SCM_REAL_VALUE (x) == SCM_COMPLEX_REAL (y)
- && SCM_COMPLEX_IMAG (y) == 0.0);
- else
- return scm_from_bool (SCM_REALP (y)
- && SCM_COMPLEX_REAL (x) == SCM_REAL_VALUE (y)
- && SCM_COMPLEX_IMAG (x) == 0.0);
- }
-
/* Vectors can be equal to one-dimensional arrays.
*/
if (scm_is_array (x) && scm_is_array (y))
diff --git a/libguile/numbers.c b/libguile/numbers.c
index 9998ab7..8513fea 100644
--- a/libguile/numbers.c
+++ b/libguile/numbers.c
@@ -3249,40 +3249,6 @@ SCM_DEFINE (scm_string_to_number, "string->number", 1, 1, 0,
/*** END strs->nums ***/
-SCM
-scm_bigequal (SCM x, SCM y)
-{
- int result = mpz_cmp (SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
- scm_remember_upto_here_2 (x, y);
- return scm_from_bool (0 == result);
-}
-
-SCM
-scm_real_equalp (SCM x, SCM y)
-{
- return scm_from_bool (SCM_REAL_VALUE (x) == SCM_REAL_VALUE (y));
-}
-
-SCM
-scm_complex_equalp (SCM x, SCM y)
-{
- return scm_from_bool (SCM_COMPLEX_REAL (x) == SCM_COMPLEX_REAL (y)
- && SCM_COMPLEX_IMAG (x) == SCM_COMPLEX_IMAG (y));
-}
-
-SCM
-scm_i_fraction_equalp (SCM x, SCM y)
-{
- if (scm_is_false (scm_equal_p (SCM_FRACTION_NUMERATOR (x),
- SCM_FRACTION_NUMERATOR (y)))
- || scm_is_false (scm_equal_p (SCM_FRACTION_DENOMINATOR (x),
- SCM_FRACTION_DENOMINATOR (y))))
- return SCM_BOOL_F;
- else
- return SCM_BOOL_T;
-}
-
-
SCM_DEFINE (scm_number_p, "number?", 1, 0, 0,
(SCM x),
"Return @code{#t} if @var{x} is a number, @code{#f}\n"
diff --git a/test-suite/tests/numbers.test b/test-suite/tests/numbers.test
index 4f30f6c..d116b6f 100644
--- a/test-suite/tests/numbers.test
+++ b/test-suite/tests/numbers.test
@@ -1605,12 +1605,24 @@
(with-test-prefix "equal?"
(pass-if (documented? equal?))
+
+ ;; The following test will fail on platforms
+ ;; without distinct signed zeroes 0.0 and -0.0.
+ (pass-if (not (equal? 0.0 -0.0)))
+
(pass-if (equal? 0 0))
(pass-if (equal? 7 7))
(pass-if (equal? -7 -7))
(pass-if (equal? (+ 1 fixnum-max) (+ 1 fixnum-max)))
(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 (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? fixnum-max (+ 1 fixnum-max))))
(pass-if (not (equal? (+ 1 fixnum-max) fixnum-max)))
(pass-if (not (equal? (+ 1 fixnum-max) (+ 2 fixnum-max))))
@@ -1631,7 +1643,10 @@
(pass-if (not (equal? (- (ash 1 1024)) -inf.0)))
(pass-if (not (equal? -inf.0 (- (ash 1 1024)))))
- (pass-if (not (equal? +nan.0 +nan.0)))
+ (pass-if (equal? +nan.0 +nan.0))
+ (pass-if (equal? +nan.0 +nan.0))
+ (pass-if (not (equal? +nan.0 0.0+nan.0i)))
+
(pass-if (not (equal? 0 +nan.0)))
(pass-if (not (equal? +nan.0 0)))
(pass-if (not (equal? 1 +nan.0)))
@@ -1655,6 +1670,75 @@
(pass-if (not (equal? +nan.0 (ash 3 1023)))))
;;;
+;;; eqv?
+;;;
+
+(with-test-prefix "eqv?"
+ (pass-if (documented? eqv?))
+
+ ;; The following test will fail on platforms
+ ;; without distinct signed zeroes 0.0 and -0.0.
+ (pass-if (not (eqv? 0.0 -0.0)))
+
+ (pass-if (eqv? 0 0))
+ (pass-if (eqv? 7 7))
+ (pass-if (eqv? -7 -7))
+ (pass-if (eqv? (+ 1 fixnum-max) (+ 1 fixnum-max)))
+ (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 (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? fixnum-max (+ 1 fixnum-max))))
+ (pass-if (not (eqv? (+ 1 fixnum-max) fixnum-max)))
+ (pass-if (not (eqv? (+ 1 fixnum-max) (+ 2 fixnum-max))))
+ (pass-if (not (eqv? fixnum-min (- fixnum-min 1))))
+ (pass-if (not (eqv? (- fixnum-min 1) fixnum-min)))
+ (pass-if (not (eqv? (- fixnum-min 1) (- fixnum-min 2))))
+ (pass-if (not (eqv? (+ fixnum-max 1) (- fixnum-min 1))))
+
+ (pass-if (not (eqv? (ash 1 256) +inf.0)))
+ (pass-if (not (eqv? +inf.0 (ash 1 256))))
+ (pass-if (not (eqv? (ash 1 256) -inf.0)))
+ (pass-if (not (eqv? -inf.0 (ash 1 256))))
+
+ ;; in gmp prior to 4.2, mpz_cmp_d ended up treating Inf as 2^1024, make
+ ;; sure we've avoided that
+ (pass-if (not (eqv? (ash 1 1024) +inf.0)))
+ (pass-if (not (eqv? +inf.0 (ash 1 1024))))
+ (pass-if (not (eqv? (- (ash 1 1024)) -inf.0)))
+ (pass-if (not (eqv? -inf.0 (- (ash 1 1024)))))
+
+ (pass-if (eqv? +nan.0 +nan.0))
+ (pass-if (not (eqv? +nan.0 0.0+nan.0i)))
+
+ (pass-if (not (eqv? 0 +nan.0)))
+ (pass-if (not (eqv? +nan.0 0)))
+ (pass-if (not (eqv? 1 +nan.0)))
+ (pass-if (not (eqv? +nan.0 1)))
+ (pass-if (not (eqv? -1 +nan.0)))
+ (pass-if (not (eqv? +nan.0 -1)))
+
+ (pass-if (not (eqv? (ash 1 256) +nan.0)))
+ (pass-if (not (eqv? +nan.0 (ash 1 256))))
+ (pass-if (not (eqv? (- (ash 1 256)) +nan.0)))
+ (pass-if (not (eqv? +nan.0 (- (ash 1 256)))))
+
+ (pass-if (not (eqv? (ash 1 8192) +nan.0)))
+ (pass-if (not (eqv? +nan.0 (ash 1 8192))))
+ (pass-if (not (eqv? (- (ash 1 8192)) +nan.0)))
+ (pass-if (not (eqv? +nan.0 (- (ash 1 8192)))))
+
+ ;; 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 (not (eqv? (ash 3 1023) +nan.0)))
+ (pass-if (not (eqv? +nan.0 (ash 3 1023)))))
+
+;;;
;;; =
;;;
--
1.5.6.5
[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #4: Infinities and NaNs are no longer rational --]
[-- Type: text/x-diff, Size: 11536 bytes --]
From ee5315abcbd527613caea504f738136d0e58274e Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Fri, 28 Jan 2011 23:32:20 -0500
Subject: [PATCH] Infinities and NaNs are no longer rational
* libguile/numbers.c (scm_rational_p): Return #f for infinities and
NaNs, per R6RS. Previously it returned #t for real infinities
and NaNs. They are still considered real by scm_real `real?'
however, per R6RS. Also simplify the code.
(scm_real_p): New implementation to reflect the fact that the
rationals and reals are no longer the same set. Previously it just
called scm_rational_p.
(scm_integer_p): Simplify the code.
* test-suite/tests/numbers.test: Add test cases for `rational?'
and `real?' applied to infinities and NaNs.
* doc/ref/api-data.texi (Real and Rational Numbers): Update docs to
reflect the fact that infinities and NaNs are no longer rational, and
that `real?' no longer implies `rational?'. Improve discussion of
infinities and NaNs.
* NEWS: Add NEWS entries, and combine with an earlier entry about
infinities no longer being integers.
---
NEWS | 18 +++++++---
doc/ref/api-data.texi | 73 ++++++++++++++++++++++-------------------
libguile/numbers.c | 40 +++++++---------------
test-suite/tests/numbers.test | 12 ++++++-
4 files changed, 76 insertions(+), 67 deletions(-)
diff --git a/NEWS b/NEWS
index 2979849..5651b17 100644
--- a/NEWS
+++ b/NEWS
@@ -27,11 +27,6 @@ Previously, `(equal? +nan.0 +nan.0)' returned #f, although
both returned #t. R5RS requires that `equal?' behave like
`eqv?' when comparing numbers.
-*** Infinities are no longer integers.
-
-Following the R6RS, infinities (+inf.0 and -inf.0) are no longer
-considered to be integers.
-
*** `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
@@ -40,6 +35,19 @@ integer-expt. This is more correct, and conforming to R6RS, but seems
to be incompatible with R5RS, which would return 0 for all non-zero
values of N.
+*** Infinities are no longer integers, nor rationals
+
+scm_integer_p `integer?' and scm_rational_p `rational?' now return #f
+for infinities, per R6RS. Previously they returned #t for real
+infinities. The real infinities and NaNs are still considered real by
+scm_real `real?' however, per R6RS.
+
+*** NaNs are no longer rationals
+
+scm_rational_p `rational?' now returns #f for NaN values, per R6RS.
+Previously it returned #t for real NaN values. They are still
+considered real by scm_real `real?' however, per R6RS.
+
*** `inf?' and `nan?' now throw exceptions for non-reals
The domain of `inf?' and `nan?' is the real numbers. Guile now signals
diff --git a/doc/ref/api-data.texi b/doc/ref/api-data.texi
index a0ab258..4256e18 100755
--- a/doc/ref/api-data.texi
+++ b/doc/ref/api-data.texi
@@ -492,10 +492,10 @@ are not rational, for example @m{\sqrt2, the square root of 2}, and
@m{\pi,pi}.
Guile can represent both exact and inexact rational numbers, but it
-can not represent irrational numbers. Exact rationals are represented
-by storing the numerator and denominator as two exact integers.
-Inexact rationals are stored as floating point numbers using the C
-type @code{double}.
+cannot represent precise finite irrational numbers. Exact rationals are
+represented by storing the numerator and denominator as two exact
+integers. Inexact rationals are stored as floating point numbers using
+the C type @code{double}.
Exact rationals are written as a fraction of integers. There must be
no whitespace around the slash:
@@ -518,26 +518,41 @@ example:
4.0
@end lisp
-The limited precision of Guile's encoding means that any ``real'' number
-in Guile can be written in a rational form, by multiplying and then dividing
-by sufficient powers of 10 (or in fact, 2). For example,
-@samp{-0.00000142857931198} is the same as @minus{}142857931198 divided by
-100000000000000000. In Guile's current incarnation, therefore, the
-@code{rational?} and @code{real?} predicates are equivalent.
-
-
-Dividing by an exact zero leads to a error message, as one might
-expect. However, dividing by an inexact zero does not produce an
-error. Instead, the result of the division is either plus or minus
-infinity, depending on the sign of the divided number.
+The limited precision of Guile's encoding means that any finite ``real''
+number in Guile can be written in a rational form, by multiplying and
+then dividing by sufficient powers of 10 (or in fact, 2). For example,
+@samp{-0.00000142857931198} is the same as @minus{}142857931198 divided
+by 100000000000000000. In Guile's current incarnation, therefore, the
+@code{rational?} and @code{real?} predicates are equivalent for finite
+numbers.
-The infinities are written @samp{+inf.0} and @samp{-inf.0},
-respectively. This syntax is also recognized by @code{read} as an
-extension to the usual Scheme syntax. The infinities are considered to
-be inexact, non-integer values.
-Dividing zero by zero yields something that is not a number at all:
-@samp{+nan.0}. This is the special `not a number' value.
+Dividing by an exact zero leads to a error message, as one might expect.
+However, dividing by an inexact zero does not produce an error.
+Instead, the result of the division is either plus or minus infinity,
+depending on the sign of the divided number and the sign of the zero
+divisor (some platforms support signed zeroes @samp{-0.0} and
+@samp{+0.0}; @samp{0.0} is the same as @samp{+0.0}).
+
+Dividing zero by an inexact zero yields a @acronym{NaN} (`not a number')
+value, although they are actually considered numbers by Scheme.
+Attempts to compare a @acronym{NaN} value with any number (including
+itself) using @code{=}, @code{<}, @code{>}, @code{<=} or @code{>=}
+always returns @code{#f}. Although a @acronym{NaN} value is not
+@code{=} to itself, it is both @code{eqv?} and @code{equal?} to itself
+and other @acronym{NaN} values. However, the preferred way to test for
+them is by using @code{nan?}.
+
+The real @acronym{NaN} values and infinities are written @samp{+nan.0},
+@samp{+inf.0} and @samp{-inf.0}. This syntax is also recognized by
+@code{read} as an extension to the usual Scheme syntax. These special
+values are considered by Scheme to be inexact real numbers but not
+rational. Note that non-real complex numbers may also contain
+infinities or @acronym{NaN} values in their real or imaginary parts. To
+test a real number to see if it is infinite, a @acronym{NaN} value, or
+neither, use @code{inf?}, @code{nan?}, or @code{finite?}, respectively.
+Every real number in Scheme belongs to precisely one of those three
+classes.
On platforms that follow @acronym{IEEE} 754 for their floating point
arithmetic, the @samp{+inf.0}, @samp{-inf.0}, and @samp{+nan.0} values
@@ -545,13 +560,6 @@ are implemented using the corresponding @acronym{IEEE} 754 values.
They behave in arithmetic operations like @acronym{IEEE} 754 describes
it, i.e., @code{(= +nan.0 +nan.0)} @result{} @code{#f}.
-While @samp{+nan.0} is not @code{=} to itself, it is @code{eqv?} to
-itself.
-
-To test for the special values, use the functions @code{inf?} and
-@code{nan?}. To test for numbers than are neither infinite nor a NaN,
-use @code{finite?}.
-
@deffn {Scheme Procedure} real? obj
@deffnx {C Function} scm_real_p (obj)
Return @code{#t} if @var{obj} is a real number, else @code{#f}. Note
@@ -566,9 +574,6 @@ Return @code{#t} if @var{x} is a rational number, @code{#f} otherwise.
Note that the set of integer values forms a subset of the set of
rational numbers, i. e. the predicate will also be fulfilled if
@var{x} is an integer number.
-
-Since Guile can not represent irrational numbers, every number
-satisfying @code{real?} also satisfies @code{rational?} in Guile.
@end deffn
@deffn {Scheme Procedure} rationalize x eps
@@ -607,12 +612,12 @@ NaN, @code{#f} otherwise.
@deffn {Scheme Procedure} nan
@deffnx {C Function} scm_nan ()
-Return NaN.
+Return @samp{+nan.0}, a @acronym{NaN} value.
@end deffn
@deffn {Scheme Procedure} inf
@deffnx {C Function} scm_inf ()
-Return Inf.
+Return @samp{+inf.0}, positive infinity.
@end deffn
@deffn {Scheme Procedure} numerator x
diff --git a/libguile/numbers.c b/libguile/numbers.c
index 8513fea..608cf7a 100644
--- a/libguile/numbers.c
+++ b/libguile/numbers.c
@@ -3281,8 +3281,8 @@ SCM_DEFINE (scm_real_p, "real?", 1, 0, 0,
"fulfilled if @var{x} is an integer number.")
#define FUNC_NAME s_scm_real_p
{
- /* we can't represent irrational numbers. */
- return scm_rational_p (x);
+ return scm_from_bool
+ (SCM_I_INUMP (x) || SCM_REALP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x));
}
#undef FUNC_NAME
@@ -3294,18 +3294,12 @@ SCM_DEFINE (scm_rational_p, "rational?", 1, 0, 0,
"fulfilled if @var{x} is an integer number.")
#define FUNC_NAME s_scm_rational_p
{
- if (SCM_I_INUMP (x))
- return SCM_BOOL_T;
- else if (SCM_IMP (x))
- return SCM_BOOL_F;
- else if (SCM_BIGP (x))
- return SCM_BOOL_T;
- else if (SCM_FRACTIONP (x))
+ if (SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x))
return SCM_BOOL_T;
else if (SCM_REALP (x))
- /* due to their limited precision, all floating point numbers are
- rational as well. */
- return SCM_BOOL_T;
+ /* due to their limited precision, finite floating point numbers are
+ rational as well. (finite means neither infinity nor a NaN) */
+ return scm_from_bool (DOUBLE_IS_FINITE (SCM_REAL_VALUE (x)));
else
return SCM_BOOL_F;
}
@@ -3317,23 +3311,15 @@ SCM_DEFINE (scm_integer_p, "integer?", 1, 0, 0,
"else.")
#define FUNC_NAME s_scm_integer_p
{
- double r;
- if (SCM_I_INUMP (x))
- return SCM_BOOL_T;
- if (SCM_IMP (x))
- return SCM_BOOL_F;
- if (SCM_BIGP (x))
+ if (SCM_I_INUMP (x) || SCM_BIGP (x))
return SCM_BOOL_T;
- if (!SCM_INEXACTP (x))
- return SCM_BOOL_F;
- if (SCM_COMPLEXP (x))
- return SCM_BOOL_F;
- r = SCM_REAL_VALUE (x);
- if (isinf (r))
+ else if (SCM_REALP (x))
+ {
+ double val = SCM_REAL_VALUE (x);
+ return scm_from_bool (!isinf (val) && (val == floor (val)));
+ }
+ else
return SCM_BOOL_F;
- if (r == floor (r))
- return SCM_BOOL_T;
- return SCM_BOOL_F;
}
#undef FUNC_NAME
diff --git a/test-suite/tests/numbers.test b/test-suite/tests/numbers.test
index d116b6f..36e3128 100644
--- a/test-suite/tests/numbers.test
+++ b/test-suite/tests/numbers.test
@@ -1505,6 +1505,11 @@
(pass-if (real? (+ 1 fixnum-max)))
(pass-if (real? (- 1 fixnum-min)))
(pass-if (real? 1.3))
+ (pass-if (real? +inf.0))
+ (pass-if (real? -inf.0))
+ (pass-if (real? +nan.0))
+ (pass-if (not (real? +inf.0-inf.0i)))
+ (pass-if (not (real? +nan.0+nan.0i)))
(pass-if (not (real? 3+4i)))
(pass-if (not (real? #\a)))
(pass-if (not (real? "a")))
@@ -1515,7 +1520,7 @@
(pass-if (not (real? (current-input-port)))))
;;;
-;;; rational? (same as real? right now)
+;;; rational?
;;;
(with-test-prefix "rational?"
@@ -1526,6 +1531,11 @@
(pass-if (rational? (+ 1 fixnum-max)))
(pass-if (rational? (- 1 fixnum-min)))
(pass-if (rational? 1.3))
+ (pass-if (not (rational? +inf.0)))
+ (pass-if (not (rational? -inf.0)))
+ (pass-if (not (rational? +nan.0)))
+ (pass-if (not (rational? +inf.0-inf.0i)))
+ (pass-if (not (rational? +nan.0+nan.0i)))
(pass-if (not (rational? 3+4i)))
(pass-if (not (rational? #\a)))
(pass-if (not (rational? "a")))
--
1.5.6.5
[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #5: Implement R6RS `real-valued?', `rational-valued?', `integer-valued?' --]
[-- Type: text/x-diff, Size: 6369 bytes --]
From b2d9e082b5740d6f722533d4ce30f3fbda955a9b Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Fri, 28 Jan 2011 23:42:01 -0500
Subject: [PATCH] Implement R6RS `real-valued?', `rational-valued?', `integer-valued?'
* module/rnrs/base.scm (real-valued?, rational-valued?,
integer-valued?): Implement in compliance with R6RS.
* test-suite/tests/r6rs-base.test: Add test cases for
`real-valued?', `rational-valued?', and `integer-valued?'.
* NEWS: Add NEWS entries.
---
NEWS | 4 ++
module/rnrs/base.scm | 19 +++++----
test-suite/tests/r6rs-base.test | 89 ++++++++++++++++++++++++++++++++++++++-
3 files changed, 103 insertions(+), 9 deletions(-)
diff --git a/NEWS b/NEWS
index 5651b17..f45795e 100644
--- a/NEWS
+++ b/NEWS
@@ -76,6 +76,10 @@ by scheme, despite their name).
throws exceptions for non-numbers. (Note that NaNs _are_ considered
numbers by scheme, despite their name).
+**** `real-valued?', `rational-valued?' and `integer-valued?' changes
+
+These predicates are now implemented in accordance with R6RS.
+
** New reader option: `hungry-eol-escapes'
Guile's string syntax is more compatible with R6RS when the
diff --git a/module/rnrs/base.scm b/module/rnrs/base.scm
index c7579c3..04a7e23 100644
--- a/module/rnrs/base.scm
+++ b/module/rnrs/base.scm
@@ -102,14 +102,17 @@
(define (exact-integer-sqrt x)
(let* ((s (exact (floor (sqrt x)))) (e (- x (* s s)))) (values s e)))
- ;; These definitions should be revisited, since the behavior of Guile's
- ;; implementations of `integer?', `rational?', and `real?' (exported from this
- ;; library) is not entirely consistent with R6RS's requirements for those
- ;; functions.
-
- (define integer-valued? integer?)
- (define rational-valued? rational?)
- (define real-valued? real?)
+ (define (real-valued? x)
+ (and (complex? x)
+ (zero? (imag-part x))))
+
+ (define (rational-valued? x)
+ (and (real-valued? x)
+ (rational? (real-part x))))
+
+ (define (integer-valued? x)
+ (and (rational-valued? x)
+ (= x (floor (real-part x)))))
(define (vector-for-each proc . vecs)
(apply for-each (cons proc (map vector->list vecs))))
diff --git a/test-suite/tests/r6rs-base.test b/test-suite/tests/r6rs-base.test
index a3603a1..1509b04 100644
--- a/test-suite/tests/r6rs-base.test
+++ b/test-suite/tests/r6rs-base.test
@@ -1,6 +1,6 @@
;;; r6rs-base.test --- Test suite for R6RS (rnrs base)
-;; Copyright (C) 2010 Free Software Foundation, Inc.
+;; Copyright (C) 2010, 2011 Free Software Foundation, Inc.
;;
;; This library is free software; you can redistribute it and/or
;; modify it under the terms of the GNU Lesser General Public
@@ -85,3 +85,90 @@
(pass-if "vector-map simple"
(equal? '#(3 2 1) (vector-map (lambda (x) (- 4 x)) '#(1 2 3)))))
+(with-test-prefix "real-valued?"
+ (pass-if (real-valued? +nan.0))
+ (pass-if (real-valued? +nan.0+0i))
+ (pass-if (real-valued? +nan.0+0.0i))
+ (pass-if (real-valued? +inf.0))
+ (pass-if (real-valued? -inf.0))
+ (pass-if (real-valued? +inf.0+0.0i))
+ (pass-if (real-valued? -inf.0-0.0i))
+ (pass-if (real-valued? 3))
+ (pass-if (real-valued? -2.5))
+ (pass-if (real-valued? -2.5+0i))
+ (pass-if (real-valued? -2.5+0.0i))
+ (pass-if (real-valued? -2.5-0i))
+ (pass-if (real-valued? #e1e10))
+ (pass-if (real-valued? 1e200))
+ (pass-if (real-valued? 1e200+0.0i))
+ (pass-if (real-valued? 6/10))
+ (pass-if (real-valued? 6/10+0.0i))
+ (pass-if (real-valued? 6/10+0i))
+ (pass-if (real-valued? 6/3))
+ (pass-if (not (real-valued? 3+i)))
+ (pass-if (not (real-valued? -2.5+0.01i)))
+ (pass-if (not (real-valued? +nan.0+0.01i)))
+ (pass-if (not (real-valued? +nan.0+nan.0i)))
+ (pass-if (not (real-valued? +inf.0-0.01i)))
+ (pass-if (not (real-valued? +0.01i)))
+ (pass-if (not (real-valued? -inf.0i))))
+
+(with-test-prefix "rational-valued?"
+ (pass-if (not (rational-valued? +nan.0)))
+ (pass-if (not (rational-valued? +nan.0+0i)))
+ (pass-if (not (rational-valued? +nan.0+0.0i)))
+ (pass-if (not (rational-valued? +inf.0)))
+ (pass-if (not (rational-valued? -inf.0)))
+ (pass-if (not (rational-valued? +inf.0+0.0i)))
+ (pass-if (not (rational-valued? -inf.0-0.0i)))
+ (pass-if (rational-valued? 3))
+ (pass-if (rational-valued? -2.5))
+ (pass-if (rational-valued? -2.5+0i))
+ (pass-if (rational-valued? -2.5+0.0i))
+ (pass-if (rational-valued? -2.5-0i))
+ (pass-if (rational-valued? #e1e10))
+ (pass-if (rational-valued? 1e200))
+ (pass-if (rational-valued? 1e200+0.0i))
+ (pass-if (rational-valued? 6/10))
+ (pass-if (rational-valued? 6/10+0.0i))
+ (pass-if (rational-valued? 6/10+0i))
+ (pass-if (rational-valued? 6/3))
+ (pass-if (not (rational-valued? 3+i)))
+ (pass-if (not (rational-valued? -2.5+0.01i)))
+ (pass-if (not (rational-valued? +nan.0+0.01i)))
+ (pass-if (not (rational-valued? +nan.0+nan.0i)))
+ (pass-if (not (rational-valued? +inf.0-0.01i)))
+ (pass-if (not (rational-valued? +0.01i)))
+ (pass-if (not (rational-valued? -inf.0i))))
+
+(with-test-prefix "integer-valued?"
+ (pass-if (not (integer-valued? +nan.0)))
+ (pass-if (not (integer-valued? +nan.0+0i)))
+ (pass-if (not (integer-valued? +nan.0+0.0i)))
+ (pass-if (not (integer-valued? +inf.0)))
+ (pass-if (not (integer-valued? -inf.0)))
+ (pass-if (not (integer-valued? +inf.0+0.0i)))
+ (pass-if (not (integer-valued? -inf.0-0.0i)))
+ (pass-if (integer-valued? 3))
+ (pass-if (integer-valued? 3.0))
+ (pass-if (integer-valued? 3+0i))
+ (pass-if (integer-valued? 3+0.0i))
+ (pass-if (integer-valued? 8/4))
+ (pass-if (integer-valued? #e1e10))
+ (pass-if (integer-valued? 1e200))
+ (pass-if (integer-valued? 1e200+0.0i))
+ (pass-if (not (integer-valued? -2.5)))
+ (pass-if (not (integer-valued? -2.5+0i)))
+ (pass-if (not (integer-valued? -2.5+0.0i)))
+ (pass-if (not (integer-valued? -2.5-0i)))
+ (pass-if (not (integer-valued? 6/10)))
+ (pass-if (not (integer-valued? 6/10+0.0i)))
+ (pass-if (not (integer-valued? 6/10+0i)))
+ (pass-if (not (integer-valued? 3+i)))
+ (pass-if (not (integer-valued? -2.5+0.01i)))
+ (pass-if (not (integer-valued? +nan.0+0.01i)))
+ (pass-if (not (integer-valued? +nan.0+nan.0i)))
+ (pass-if (not (integer-valued? +inf.0-0.01i)))
+ (pass-if (not (integer-valued? +0.01i)))
+ (pass-if (not (integer-valued? -inf.0i))))
+
--
1.5.6.5
[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #6: Add SCM_LIKELY and SCM_UNLIKELY for optimization --]
[-- Type: text/x-diff, Size: 4974 bytes --]
From 3fdddf143b231b989c78c3f7875d367eb42e72cd Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Fri, 28 Jan 2011 23:58:02 -0500
Subject: [PATCH] Add SCM_LIKELY and SCM_UNLIKELY for optimization
* libguile/numbers.c (scm_abs, scm_quotient, scm_remainder, scm_modulo):
Add SCM_LIKELY and SCM_UNLIKELY in several places for optimization.
(scm_remainder): Add comment about C99 "%" semantics.
Strip away a redundant set of braces.
---
libguile/numbers.c | 67 ++++++++++++++++++++++++++-------------------------
1 files changed, 34 insertions(+), 33 deletions(-)
diff --git a/libguile/numbers.c b/libguile/numbers.c
index 608cf7a..0fae4cb 100644
--- a/libguile/numbers.c
+++ b/libguile/numbers.c
@@ -728,7 +728,7 @@ SCM_PRIMITIVE_GENERIC (scm_abs, "abs", 1, 0, 0,
"Return the absolute value of @var{x}.")
#define FUNC_NAME
{
- if (SCM_I_INUMP (x))
+ if (SCM_LIKELY (SCM_I_INUMP (x)))
{
scm_t_inum xx = SCM_I_INUM (x);
if (xx >= 0)
@@ -774,18 +774,18 @@ SCM_GPROC (s_quotient, "quotient", 2, 0, 0, scm_quotient, g_quotient);
SCM
scm_quotient (SCM x, SCM y)
{
- if (SCM_I_INUMP (x))
+ if (SCM_LIKELY (SCM_I_INUMP (x)))
{
scm_t_inum xx = SCM_I_INUM (x);
- if (SCM_I_INUMP (y))
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
{
scm_t_inum yy = SCM_I_INUM (y);
- if (yy == 0)
+ if (SCM_UNLIKELY (yy == 0))
scm_num_overflow (s_quotient);
else
{
scm_t_inum z = xx / yy;
- if (SCM_FIXABLE (z))
+ if (SCM_LIKELY (SCM_FIXABLE (z)))
return SCM_I_MAKINUM (z);
else
return scm_i_inum2big (z);
@@ -809,12 +809,12 @@ scm_quotient (SCM x, SCM y)
}
else if (SCM_BIGP (x))
{
- if (SCM_I_INUMP (y))
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
{
scm_t_inum yy = SCM_I_INUM (y);
- if (yy == 0)
+ if (SCM_UNLIKELY (yy == 0))
scm_num_overflow (s_quotient);
- else if (yy == 1)
+ else if (SCM_UNLIKELY (yy == 1))
return x;
else
{
@@ -858,15 +858,18 @@ SCM_GPROC (s_remainder, "remainder", 2, 0, 0, scm_remainder, g_remainder);
SCM
scm_remainder (SCM x, SCM y)
{
- if (SCM_I_INUMP (x))
+ if (SCM_LIKELY (SCM_I_INUMP (x)))
{
- if (SCM_I_INUMP (y))
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
{
scm_t_inum yy = SCM_I_INUM (y);
- if (yy == 0)
+ if (SCM_UNLIKELY (yy == 0))
scm_num_overflow (s_remainder);
else
{
+ /* C99 specifies that "%" is the remainder corresponding to a
+ quotient rounded towards zero, and that's also traditional
+ for machine division, so z here should be well defined. */
scm_t_inum z = SCM_I_INUM (x) % yy;
return SCM_I_MAKINUM (z);
}
@@ -889,10 +892,10 @@ scm_remainder (SCM x, SCM y)
}
else if (SCM_BIGP (x))
{
- if (SCM_I_INUMP (y))
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
{
scm_t_inum yy = SCM_I_INUM (y);
- if (yy == 0)
+ if (SCM_UNLIKELY (yy == 0))
scm_num_overflow (s_remainder);
else
{
@@ -931,13 +934,13 @@ SCM_GPROC (s_modulo, "modulo", 2, 0, 0, scm_modulo, g_modulo);
SCM
scm_modulo (SCM x, SCM y)
{
- if (SCM_I_INUMP (x))
+ if (SCM_LIKELY (SCM_I_INUMP (x)))
{
scm_t_inum xx = SCM_I_INUM (x);
- if (SCM_I_INUMP (y))
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
{
scm_t_inum yy = SCM_I_INUM (y);
- if (yy == 0)
+ if (SCM_UNLIKELY (yy == 0))
scm_num_overflow (s_modulo);
else
{
@@ -1008,10 +1011,10 @@ scm_modulo (SCM x, SCM y)
}
else if (SCM_BIGP (x))
{
- if (SCM_I_INUMP (y))
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
{
scm_t_inum yy = SCM_I_INUM (y);
- if (yy == 0)
+ if (SCM_UNLIKELY (yy == 0))
scm_num_overflow (s_modulo);
else
{
@@ -1029,22 +1032,20 @@ scm_modulo (SCM x, SCM y)
}
else if (SCM_BIGP (y))
{
- {
- SCM result = scm_i_mkbig ();
- int y_sgn = mpz_sgn (SCM_I_BIG_MPZ (y));
- SCM pos_y = scm_i_clonebig (y, y_sgn >= 0);
- mpz_mod (SCM_I_BIG_MPZ (result),
- SCM_I_BIG_MPZ (x),
- SCM_I_BIG_MPZ (pos_y));
+ SCM result = scm_i_mkbig ();
+ int y_sgn = mpz_sgn (SCM_I_BIG_MPZ (y));
+ SCM pos_y = scm_i_clonebig (y, y_sgn >= 0);
+ mpz_mod (SCM_I_BIG_MPZ (result),
+ SCM_I_BIG_MPZ (x),
+ SCM_I_BIG_MPZ (pos_y));
- scm_remember_upto_here_1 (x);
- if ((y_sgn < 0) && (mpz_sgn (SCM_I_BIG_MPZ (result)) != 0))
- mpz_add (SCM_I_BIG_MPZ (result),
- SCM_I_BIG_MPZ (y),
- SCM_I_BIG_MPZ (result));
- scm_remember_upto_here_2 (y, pos_y);
- return scm_i_normbig (result);
- }
+ scm_remember_upto_here_1 (x);
+ if ((y_sgn < 0) && (mpz_sgn (SCM_I_BIG_MPZ (result)) != 0))
+ mpz_add (SCM_I_BIG_MPZ (result),
+ SCM_I_BIG_MPZ (y),
+ SCM_I_BIG_MPZ (result));
+ scm_remember_upto_here_2 (y, pos_y);
+ return scm_i_normbig (result);
}
else
SCM_WTA_DISPATCH_2 (g_modulo, x, y, SCM_ARG2, s_modulo);
--
1.5.6.5
[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #7: Implement efficient R6RS `div', `mod', et al --]
[-- Type: text/x-diff, Size: 60957 bytes --]
From a1dda78005b13c4b9dfa97b636f21e62dd3b0f38 Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Sat, 29 Jan 2011 02:36:02 -0500
Subject: [PATCH] Implement efficient R6RS `div', `mod', et al
* libguile/numbers.c (scm_div, scm_mod, scm_div_and_mod, scm_div0,
scm_mod0, scm_div0_and_mod0): New extensible procedures `div', `mod',
`div-and-mod', `div0', `mod0', `div0-and-mod0'.
(scm_i_inexact_div, scm_i_inexact_mod, scm_i_inexact_div_and_mod,
scm_i_inexact_div0, scm_i_inexact_mod0, scm_i_inexact_div0_and_mod0,
scm_i_slow_exact_div, scm_i_slow_exact_mod, scm_i_slow_exact_div_and_mod,
scm_i_slow_exact_div, scm_i_slow_exact_mod, scm_i_slow_exact_div_and_mod,
scm_i_bigint_div0, scm_i_bigint_mod0, scm_i_bigint_div0_and_mod0):
New internal static procedures, not intended to be used except by
scm_div, scm_mod, scm_div_and_mod, scm_div0, scm_mod0,
and scm_div0_and_mod0.
* libguile/numbers.h: Add function prototypes.
* module/rnrs/base.scm: Remove incorrect stub implementations of `div',
`mod', `div-and-mod', `div0', `mod0', and `div0-and-mod0'.
* module/rnrs/arithmetic/fixnums.scm (fxdiv, fxmod, fxdiv-and-mod,
fxdiv0, fxmod0, fxdiv0-and-mod0): Remove redundant checks for division
by zero and unnecessary complexity.
(fx+/carry): Remove unneeded calls to `inexact->exact'.
* module/rnrs/arithmetic/flonums.scm (fldiv, flmod, fldiv-and-mod,
fldiv0, flmod0, fldiv0-and-mod0): Remove redundant checks for division
by zero and unnecessary complexity. Remove unneeded calls to
`inexact->exact' and `exact->inexact'
* test-suite/tests/numbers.test: (test-eqv?): New internal predicate for
comparing numerical outputs with expected values.
Add extensive test code for `div', `mod', `div-and-mod', `div0',
`mod0', and `div0-and-mod0'.
* test-suite/tests/r6rs-arithmetic-fixnums.test: Fix some broken test
cases, and remove `unresolved' test markers for `fxdiv', `fxmod',
`fxdiv-and-mod', `fxdiv0', `fxmod0', and `fxdiv0-and-mod0'.
* test-suite/tests/r6rs-arithmetic-flonums.test: Remove `unresolved'
test markers for `fldiv', `flmod', `fldiv-and-mod', `fldiv0',
`flmod0', and `fldiv0-and-mod0'.
* doc/ref/api-data.texi (Arithmetic): Document `div', `mod',
`div-and-mod', `div0', `mod0', and `div0-and-mod0'.
(Operations on Integer Values): Add cross-references to `div', `mod',
et al, from `quotient', `remainder', and `modulo'.
* doc/ref/r6rs.texi (rnrs base): Remove stub descriptions for `div',
`mod', `div-and-mod', `div0', `mod0', and `div0-and-mod0'. Instead,
cross reference to their descriptions in the core arithmetic section.
* NEWS: Add NEWS entry.
---
NEWS | 16 +
doc/ref/api-data.texi | 67 ++
doc/ref/r6rs.texi | 19 +-
libguile/numbers.c | 1172 ++++++++++++++++++++++++-
libguile/numbers.h | 6 +
module/rnrs/arithmetic/fixnums.scm | 23 +-
module/rnrs/arithmetic/flonums.scm | 31 +-
module/rnrs/base.scm | 17 -
test-suite/tests/numbers.test | 166 ++++-
test-suite/tests/r6rs-arithmetic-fixnums.test | 23 +-
test-suite/tests/r6rs-arithmetic-flonums.test | 9 +-
11 files changed, 1457 insertions(+), 92 deletions(-)
diff --git a/NEWS b/NEWS
index f45795e..085f2b9 100644
--- a/NEWS
+++ b/NEWS
@@ -12,6 +12,22 @@ Changes in 1.9.15 (since the 1.9.14 prerelease):
** Changes and bugfixes in numerics code
+**** New procedures: `div', `mod', `div-and-mod' et al
+
+Added efficient R6RS division operations to Guile core. These
+procedures each accept two real numbers X and Y, where Y must be
+non-zero. `div' returns an integer Q and `mod' returns a real R such
+that X = R + Q * Y and 0 <= R < abs(Y). `div-and-mod' returns both Q
+and R, and is more efficient than calling `div' and `mod' separately.
+`div0', `mod0', and `div0-and-mod0' are similar except that
+-abs(Y/2) <= R < abs(Y/2).
+
+**** `div0', `mod0', and `div0-and-mod0' now implemented correctly
+
+These functions are now implemented correctly (though admittedly
+inefficiently). R6RS states that (div0-and-mod0 123 -10) should
+return -12 and 3, but previously it returned -12 and -7.
+
*** `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 4256e18..41702a9 100755
--- a/doc/ref/api-data.texi
+++ b/doc/ref/api-data.texi
@@ -897,6 +897,9 @@ sign as @var{n}. In all cases quotient and remainder satisfy
(remainder 13 4) @result{} 1
(remainder -13 4) @result{} -1
@end lisp
+
+See also @code{div}, @code{mod} and related operations in
+@ref{Arithmetic}.
@end deffn
@c begin (texi-doc-string "guile" "modulo")
@@ -911,6 +914,9 @@ sign as @var{d}.
(modulo 13 -4) @result{} -3
(modulo -13 -4) @result{} -1
@end lisp
+
+See also @code{div}, @code{mod} and related operations in
+@ref{Arithmetic}.
@end deffn
@c begin (texi-doc-string "guile" "gcd")
@@ -1130,6 +1136,12 @@ Returns the magnitude or angle of @var{z} as a @code{double}.
@rnindex ceiling
@rnindex truncate
@rnindex round
+@rnindex div
+@rnindex mod
+@rnindex div-and-mod
+@rnindex div0
+@rnindex mod0
+@rnindex div0-and-mod0
The C arithmetic functions below always takes two arguments, while the
Scheme functions can take an arbitrary number. When you need to
@@ -1229,6 +1241,61 @@ respectively, but these functions take and return @code{double}
values.
@end deftypefn
+@deffn {Scheme Procedure} div x y
+@deffnx {Scheme Procedure} mod x y
+@deffnx {Scheme Procedure} div-and-mod x y
+@deffnx {C Function} scm_div (x y)
+@deffnx {C Function} scm_mod (x y)
+@deffnx {C Function} scm_div_and_mod (x y)
+These procedures implement number-theoretic division.
+
+Each accepts two real numbers @var{x} and @var{y}, where @var{y} is
+non-zero. @code{div} returns an integer @var{q} and @code{mod} returns
+a real @var{r} such that @math{@var{x} = @var{r} + @var{q}*@var{y}} and
+@math{0 <= @var{r} < abs(@var{y})}. @code{div-and-mod} returns both
+values, and is more efficient than calling @code{div} and @code{mod}
+separately.
+
+@lisp
+(div 123 10) @result{} 12
+(mod 123 10) @result{} 3
+(div-and-mod 123 10) @result{} 12 and 3
+(div-and-mod 123 -10) @result{} -12 and 3
+(div-and-mod -123 10) @result{} -13 and 7
+(div-and-mod -123 -10) @result{} 13 and 7
+(div-and-mod -123.2 -63.5) @result{} 2.0 and 3.8
+(div-and-mod 125/7 -10/7) @result{} -12 and 5/7
+@end lisp
+@end deffn
+
+@deffn {Scheme Procedure} div0 x y
+@deffnx {Scheme Procedure} mod0 x y
+@deffnx {Scheme Procedure} div0-and-mod0 x y
+@deffnx {C Function} scm_div0 (x y)
+@deffnx {C Function} scm_mod0 (x y)
+@deffnx {C Function} scm_div0_and_mod0 (x y)
+These procedures are similar to @code{div}, @code{mod}, and
+@code{div-and-mod}, except that @code{mod0} returns values that lie
+within a half-open interval centered on zero.
+
+Precisely, @code{div0} returns an integer @var{q} and @code{mod0}
+returns a real @var{r} such that @math{@var{x} = @var{r} +
+@var{q}*@var{y}} and @math{-abs(@var{y}/2) <= @var{r} < abs(@var{y}/2)}.
+@code{div0-and-mod0} returns both values, and is more efficient than
+calling @code{div0} and @code{mod0} separately.
+
+@lisp
+(div0 123 10) @result{} 12
+(mod0 123 10) @result{} 3
+(div0-and-mod0 123 10) @result{} 12 and 3
+(div0-and-mod0 123 -10) @result{} -12 and 3
+(div0-and-mod0 -123 10) @result{} -12 and -3
+(div0-and-mod0 -123 -10) @result{} 12 and -3
+(div0-and-mod0 -123.2 -63.5) @result{} 2.0 and 3.8
+(div0-and-mod0 125/7 -10/7) @result{} -13 and -5/7
+@end lisp
+@end deffn
+
@node Scientific
@subsubsection Scientific Functions
diff --git a/doc/ref/r6rs.texi b/doc/ref/r6rs.texi
index 5fee65f..6439478 100644
--- a/doc/ref/r6rs.texi
+++ b/doc/ref/r6rs.texi
@@ -1,6 +1,6 @@
@c -*-texinfo-*-
@c This is part of the GNU Guile Reference Manual.
-@c Copyright (C) 2010
+@c Copyright (C) 2010, 2011
@c Free Software Foundation, Inc.
@c See the file guile.texi for copying conditions.
@@ -461,24 +461,13 @@ grouped below by the existing manual sections to which they correspond.
@deffnx {Scheme Procedure} floor x
@deffnx {Scheme Procedure} ceiling x
@deffnx {Scheme Procedure} round x
-@xref{Arithmetic}, for documentation.
-@end deffn
-
-@deffn {Scheme Procedure} div x1 x2
+@deffnx {Scheme Procedure} div x1 x2
@deffnx {Scheme Procedure} mod x1 x2
@deffnx {Scheme Procedure} div-and-mod x1 x2
-These procedures implement number-theoretic division.
-
-@code{div-and-mod} returns two values, the respective results of
-@code{(div x1 x2)} and @code{(mod x1 x2)}.
-@end deffn
-
-@deffn {Scheme Procedure} div0 x1 x2
+@deffnx {Scheme Procedure} div0 x1 x2
@deffnx {Scheme Procedure} mod0 x1 x2
@deffnx {Scheme Procedure} div0-and-mod0 x1 x2
-These procedures are similar to @code{div}, @code{mod}, and
-@code{div-and-mod}, except that @code{mod0} returns values that lie
-within a half-open interval centered on zero.
+@xref{Arithmetic}, for documentation.
@end deffn
@deffn {Scheme Procedure} exact-integer-sqrt k
diff --git a/libguile/numbers.c b/libguile/numbers.c
index 0fae4cb..5ade135 100644
--- a/libguile/numbers.c
+++ b/libguile/numbers.c
@@ -105,6 +105,7 @@ typedef scm_t_signed_bits scm_t_inum;
static SCM flo0;
+static SCM exactly_one_half;
#define SCM_SWAP(x, y) do { SCM __t = x; x = y; y = __t; } while (0)
@@ -1054,6 +1055,1175 @@ scm_modulo (SCM x, SCM y)
SCM_WTA_DISPATCH_2 (g_modulo, x, y, SCM_ARG1, s_modulo);
}
+static SCM scm_i_inexact_div (double x, double y);
+static SCM scm_i_slow_exact_div (SCM x, SCM y);
+
+SCM_GPROC (s_div, "div", 2, 0, 0, scm_div, g_div);
+/* "Return q = @var{x} div @var{y}, where x = r + q*y,\n"
+ * "q is an integer and 0 <= r < abs(y)."
+ * "@lisp\n"
+ * "(div 123 10) @result{} 12\n"
+ * "(div 123 -10) @result{} -12\n"
+ * "(div -123 10) @result{} -13\n"
+ * "(div -123 -10) @result{} 13\n"
+ * "@end lisp"
+ */
+SCM
+scm_div (SCM x, SCM y)
+{
+ if (SCM_LIKELY (SCM_I_INUMP (x)))
+ {
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
+ {
+ scm_t_inum yy = SCM_I_INUM (y);
+ if (SCM_UNLIKELY (yy == 0))
+ scm_num_overflow (s_div);
+ else
+ {
+ scm_t_inum xx = SCM_I_INUM (x);
+ scm_t_inum qq = xx / yy;
+ if (xx < 0 && xx < qq * yy)
+ {
+ if (yy > 0)
+ qq--;
+ else
+ qq++;
+ }
+ return SCM_I_MAKINUM (qq);
+ }
+ }
+ else if (SCM_BIGP (y))
+ {
+ if (SCM_I_INUM (x) >= 0)
+ return SCM_INUM0;
+ else
+ return SCM_I_MAKINUM (- mpz_sgn (SCM_I_BIG_MPZ (y)));
+ }
+ else if (SCM_REALP (y))
+ return scm_i_inexact_div (SCM_I_INUM (x), SCM_REAL_VALUE (y));
+ else if (SCM_FRACTIONP (y))
+ return scm_i_slow_exact_div (x, y);
+ else
+ SCM_WTA_DISPATCH_2 (g_div, x, y, SCM_ARG2, s_div);
+ }
+ else if (SCM_BIGP (x))
+ {
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
+ {
+ scm_t_inum yy = SCM_I_INUM (y);
+ if (SCM_UNLIKELY (yy == 0))
+ scm_num_overflow (s_div);
+ else
+ {
+ SCM q = scm_i_mkbig ();
+ if (yy > 0)
+ mpz_fdiv_q_ui (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (x), yy);
+ else
+ {
+ mpz_fdiv_q_ui (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (x), -yy);
+ mpz_neg (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (q));
+ }
+ scm_remember_upto_here_1 (x);
+ return scm_i_normbig (q);
+ }
+ }
+ else if (SCM_BIGP (y))
+ {
+ SCM q = scm_i_mkbig ();
+ if (mpz_sgn (SCM_I_BIG_MPZ (y)) > 0)
+ mpz_fdiv_q (SCM_I_BIG_MPZ (q),
+ SCM_I_BIG_MPZ (x),
+ SCM_I_BIG_MPZ (y));
+ else
+ mpz_cdiv_q (SCM_I_BIG_MPZ (q),
+ SCM_I_BIG_MPZ (x),
+ SCM_I_BIG_MPZ (y));
+ scm_remember_upto_here_2 (x, y);
+ return scm_i_normbig (q);
+ }
+ else if (SCM_REALP (y))
+ return scm_i_inexact_div (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+ else if (SCM_FRACTIONP (y))
+ return scm_i_slow_exact_div (x, y);
+ else
+ SCM_WTA_DISPATCH_2 (g_div, x, y, SCM_ARG2, s_div);
+ }
+ else if (SCM_REALP (x))
+ {
+ if (!(SCM_REALP (y) || SCM_I_INUMP (y) ||
+ SCM_BIGP (y) || SCM_FRACTIONP (y)))
+ SCM_WTA_DISPATCH_2 (g_div, x, y, SCM_ARG2, s_div);
+ else
+ return scm_i_inexact_div (SCM_REAL_VALUE (x), scm_to_double (y));
+ }
+ else if (SCM_FRACTIONP (x))
+ {
+ if (SCM_REALP (y))
+ return scm_i_inexact_div (scm_i_fraction2double (x),
+ SCM_REAL_VALUE (y));
+ else
+ return scm_i_slow_exact_div (x, y);
+ }
+ else
+ SCM_WTA_DISPATCH_2 (g_div, x, y, SCM_ARG1, s_div);
+}
+
+static SCM
+scm_i_inexact_div (double x, double y)
+{
+ if (SCM_LIKELY (y > 0))
+ return scm_from_double (floor(x / y));
+ else if (SCM_LIKELY (y < 0))
+ return scm_from_double (ceil(x / y));
+ else if (y == 0)
+ scm_num_overflow (s_div); /* or should we return a NaN? */
+ else
+ return scm_nan ();
+}
+
+/* Compute exact div the slow way.
+ We use this only if both arguments are exact,
+ and at least one of them is a fraction */
+static SCM
+scm_i_slow_exact_div (SCM x, SCM y)
+{
+ if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+ SCM_WTA_DISPATCH_2 (g_div, x, y, SCM_ARG1, s_div);
+ else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+ SCM_WTA_DISPATCH_2 (g_div, x, y, SCM_ARG2, s_div);
+ else if (scm_is_true (scm_positive_p (y)))
+ return scm_floor (scm_divide (x, y));
+ else if (scm_is_true (scm_negative_p (y)))
+ return scm_ceiling (scm_divide (x, y));
+ else
+ scm_num_overflow (s_div);
+}
+
+static SCM scm_i_inexact_mod (double x, double y);
+static SCM scm_i_slow_exact_mod (SCM x, SCM y);
+
+SCM_GPROC (s_mod, "mod", 2, 0, 0, scm_mod, g_mod);
+/* "Return r = @var{x} mod @var{y}, where x = r + q*y,\n"
+ * "q is an integer and 0 <= r < abs(y)."
+ * "@lisp\n"
+ * "(mod 123 10) @result{} 3\n"
+ * "(mod 123 -10) @result{} 3\n"
+ * "(mod -123 10) @result{} 7\n"
+ * "(mod -123 -10) @result{} 7\n"
+ * "@end lisp"
+ */
+SCM
+scm_mod (SCM x, SCM y)
+{
+ if (SCM_LIKELY (SCM_I_INUMP (x)))
+ {
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
+ {
+ scm_t_inum yy = SCM_I_INUM (y);
+ if (SCM_UNLIKELY (yy == 0))
+ scm_num_overflow (s_mod);
+ else
+ {
+ scm_t_inum rr = SCM_I_INUM (x) % yy;
+ if (rr >= 0)
+ return SCM_I_MAKINUM (rr);
+ else if (yy > 0)
+ return SCM_I_MAKINUM (rr + yy);
+ else
+ return SCM_I_MAKINUM (rr - yy);
+ }
+ }
+ else if (SCM_BIGP (y))
+ {
+ scm_t_inum xx = SCM_I_INUM (x);
+ if ((xx == SCM_MOST_NEGATIVE_FIXNUM) &&
+ (0 == mpz_cmp_ui (SCM_I_BIG_MPZ (y),
+ - SCM_MOST_NEGATIVE_FIXNUM)))
+ {
+ /* Special case: x == fixnum-min && y == abs (fixnum-min) */
+ scm_remember_upto_here_1 (y);
+ return SCM_INUM0;
+ }
+ else if (xx >= 0)
+ return x;
+ else if (mpz_sgn (SCM_I_BIG_MPZ (y)) > 0)
+ {
+ SCM r = scm_i_mkbig ();
+ mpz_sub_ui (SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (y), -xx);
+ scm_remember_upto_here_1 (y);
+ return scm_i_normbig (r);
+ }
+ else
+ {
+ SCM r = scm_i_mkbig ();
+ mpz_add_ui (SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (y), -xx);
+ scm_remember_upto_here_1 (y);
+ mpz_neg (SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (r));
+ return scm_i_normbig (r);
+ }
+ }
+ else if (SCM_REALP (y))
+ return scm_i_inexact_mod (SCM_I_INUM (x), SCM_REAL_VALUE (y));
+ else if (SCM_FRACTIONP (y))
+ return scm_i_slow_exact_mod (x, y);
+ else
+ SCM_WTA_DISPATCH_2 (g_mod, x, y, SCM_ARG2, s_mod);
+ }
+ else if (SCM_BIGP (x))
+ {
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
+ {
+ scm_t_inum yy = SCM_I_INUM (y);
+ if (SCM_UNLIKELY (yy == 0))
+ scm_num_overflow (s_mod);
+ else
+ {
+ scm_t_inum rr;
+ if (yy < 0)
+ yy = - yy;
+ rr = mpz_fdiv_ui (SCM_I_BIG_MPZ (x), yy);
+ scm_remember_upto_here_1 (x);
+ return SCM_I_MAKINUM (rr);
+ }
+ }
+ else if (SCM_BIGP (y))
+ {
+ SCM r = scm_i_mkbig ();
+ mpz_mod (SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
+ scm_remember_upto_here_2 (x, y);
+ return scm_i_normbig (r);
+ }
+ else if (SCM_REALP (y))
+ return scm_i_inexact_mod (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+ else if (SCM_FRACTIONP (y))
+ return scm_i_slow_exact_mod (x, y);
+ else
+ SCM_WTA_DISPATCH_2 (g_mod, x, y, SCM_ARG2, s_mod);
+ }
+ else if (SCM_REALP (x))
+ {
+ if (!(SCM_REALP (y) || SCM_I_INUMP (y) ||
+ SCM_BIGP (y) || SCM_FRACTIONP (y)))
+ SCM_WTA_DISPATCH_2 (g_mod, x, y, SCM_ARG2, s_mod);
+ else
+ return scm_i_inexact_mod (SCM_REAL_VALUE (x), scm_to_double (y));
+ }
+ else if (SCM_FRACTIONP (x))
+ {
+ if (SCM_REALP (y))
+ return scm_i_inexact_mod (scm_i_fraction2double (x),
+ SCM_REAL_VALUE (y));
+ else
+ return scm_i_slow_exact_mod (x, y);
+ }
+ else
+ SCM_WTA_DISPATCH_2 (g_mod, x, y, SCM_ARG1, s_mod);
+}
+
+static SCM
+scm_i_inexact_mod (double x, double y)
+{
+ double q;
+
+ /* Although it would be more efficient to use fmod here, we can't
+ because it would in some cases produce results inconsistent with
+ scm_i_inexact_div, such that x != r + q * y (not even close). In
+ particular, when x is very close to a multiple of y, then r might
+ be either 0.0 or abs(y)-epsilon, but those two cases must
+ correspond with different choices of q. If r = 0.0 then q must be
+ x/y, and if r = abs(y) then q must be (x-r)/y. If div chooses one
+ way and mod chooses the other, it would be bad. This problem
+ actually happened with (div 130.0 10/7) and (mod 130.0 10/7) on one
+ platform. */
+ if (SCM_LIKELY (y > 0))
+ q = floor(x / y);
+ else if (SCM_LIKELY (y < 0))
+ q = ceil(x / y);
+ else if (y == 0)
+ scm_num_overflow (s_mod); /* or should we return a NaN? */
+ else
+ return scm_nan ();
+ return scm_from_double (x - q * y);
+}
+
+/* Compute exact mod the slow way: x-y*(x div y)
+ We use this only if both arguments are exact,
+ and at least one of them is a fraction */
+static SCM
+scm_i_slow_exact_mod (SCM x, SCM y)
+{
+ if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+ SCM_WTA_DISPATCH_2 (g_mod, x, y, SCM_ARG1, s_mod);
+ else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+ SCM_WTA_DISPATCH_2 (g_mod, x, y, SCM_ARG2, s_mod);
+ else if (scm_is_true (scm_positive_p (y)))
+ return scm_difference
+ (x, scm_product (y, scm_floor (scm_divide (x, y))));
+ else if (scm_is_true (scm_negative_p (y)))
+ return scm_difference
+ (x, scm_product (y, scm_ceiling (scm_divide (x, y))));
+ else
+ scm_num_overflow (s_mod);
+}
+
+
+static SCM scm_i_inexact_div_and_mod (double x, double y);
+static SCM scm_i_slow_exact_div_and_mod (SCM x, SCM y);
+
+SCM_GPROC (s_div_and_mod, "div-and-mod", 2, 0, 0,
+ scm_div_and_mod, g_div_and_mod);
+/* "Return q and r, where x = r + q*y,"
+ * "q is an integer, and 0 <= r < abs(y)."
+ * "@lisp\n"
+ * "(div-and-mod 123 10) @result{} 12 and 3\n"
+ * "(div-and-mod 123 -10) @result{} -12 and 3\n"
+ * "(div-and-mod -123 10) @result{} -13 and 7\n"
+ * "(div-and-mod -123 -10) @result{} 13 and 7\n"
+ * "@end lisp"
+ */
+SCM
+scm_div_and_mod (SCM x, SCM y)
+{
+ if (SCM_LIKELY (SCM_I_INUMP (x)))
+ {
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
+ {
+ scm_t_inum yy = SCM_I_INUM (y);
+ if (SCM_UNLIKELY (yy == 0))
+ scm_num_overflow (s_div_and_mod);
+ else
+ {
+ scm_t_inum xx = SCM_I_INUM (x);
+ scm_t_inum qq = xx / yy;
+ scm_t_inum rr = xx - qq * yy;
+ if (rr < 0)
+ {
+ if (yy > 0)
+ { rr += yy; qq--; }
+ else
+ { rr -= yy; qq++; }
+ }
+ return scm_values (scm_list_2 (SCM_I_MAKINUM (qq),
+ SCM_I_MAKINUM (rr)));
+ }
+ }
+ else if (SCM_BIGP (y))
+ {
+ scm_t_inum xx = SCM_I_INUM (x);
+ if (xx >= 0)
+ return scm_values (scm_list_2 (SCM_INUM0, x));
+ else if ((xx == SCM_MOST_NEGATIVE_FIXNUM) &&
+ (0 == mpz_cmp_ui (SCM_I_BIG_MPZ (y),
+ - SCM_MOST_NEGATIVE_FIXNUM)))
+ {
+ /* Special case: x == fixnum-min && y == abs (fixnum-min) */
+ scm_remember_upto_here_1 (y);
+ return scm_values
+ (scm_list_2 (SCM_I_MAKINUM (-1), SCM_INUM0));
+ }
+ else if (mpz_sgn (SCM_I_BIG_MPZ (y)) > 0)
+ {
+ SCM r = scm_i_mkbig ();
+ mpz_sub_ui (SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (y), -xx);
+ scm_remember_upto_here_1 (y);
+ return scm_values
+ (scm_list_2 (SCM_I_MAKINUM (-1), scm_i_normbig (r)));
+ }
+ else
+ {
+ SCM r = scm_i_mkbig ();
+ mpz_add_ui (SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (y), -xx);
+ scm_remember_upto_here_1 (y);
+ mpz_neg (SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (r));
+ return scm_values (scm_list_2 (SCM_INUM1, scm_i_normbig (r)));
+ }
+ }
+ else if (SCM_REALP (y))
+ return scm_i_inexact_div_and_mod (SCM_I_INUM (x), SCM_REAL_VALUE (y));
+ else if (SCM_FRACTIONP (y))
+ return scm_i_slow_exact_div_and_mod (x, y);
+ else
+ SCM_WTA_DISPATCH_2 (g_div_and_mod, x, y, SCM_ARG2, s_div_and_mod);
+ }
+ else if (SCM_BIGP (x))
+ {
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
+ {
+ scm_t_inum yy = SCM_I_INUM (y);
+ if (SCM_UNLIKELY (yy == 0))
+ scm_num_overflow (s_div_and_mod);
+ else
+ {
+ SCM q = scm_i_mkbig ();
+ SCM r = scm_i_mkbig ();
+ if (yy > 0)
+ mpz_fdiv_qr_ui (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (x), yy);
+ else
+ {
+ mpz_fdiv_qr_ui (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (x), -yy);
+ mpz_neg (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (q));
+ }
+ scm_remember_upto_here_1 (x);
+ return scm_values (scm_list_2 (scm_i_normbig (q),
+ scm_i_normbig (r)));
+ }
+ }
+ else if (SCM_BIGP (y))
+ {
+ SCM q = scm_i_mkbig ();
+ SCM r = scm_i_mkbig ();
+ if (mpz_sgn (SCM_I_BIG_MPZ (y)) > 0)
+ mpz_fdiv_qr (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
+ else
+ mpz_cdiv_qr (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
+ scm_remember_upto_here_2 (x, y);
+ return scm_values (scm_list_2 (scm_i_normbig (q),
+ scm_i_normbig (r)));
+ }
+ else if (SCM_REALP (y))
+ return scm_i_inexact_div_and_mod (scm_i_big2dbl (x),
+ SCM_REAL_VALUE (y));
+ else if (SCM_FRACTIONP (y))
+ return scm_i_slow_exact_div_and_mod (x, y);
+ else
+ SCM_WTA_DISPATCH_2 (g_div_and_mod, x, y, SCM_ARG2, s_div_and_mod);
+ }
+ else if (SCM_REALP (x))
+ {
+ if (!(SCM_REALP (y) || SCM_I_INUMP (y) ||
+ SCM_BIGP (y) || SCM_FRACTIONP (y)))
+ SCM_WTA_DISPATCH_2 (g_div_and_mod, x, y, SCM_ARG2, s_div_and_mod);
+ else
+ return scm_i_inexact_div_and_mod (SCM_REAL_VALUE (x),
+ scm_to_double (y));
+ }
+ else if (SCM_FRACTIONP (x))
+ {
+ if (SCM_REALP (y))
+ return scm_i_inexact_div_and_mod (scm_i_fraction2double (x),
+ SCM_REAL_VALUE (y));
+ else
+ return scm_i_slow_exact_div_and_mod (x, y);
+ }
+ else
+ SCM_WTA_DISPATCH_2 (g_div_and_mod, x, y, SCM_ARG1, s_div_and_mod);
+}
+
+static SCM
+scm_i_inexact_div_and_mod (double x, double y)
+{
+ double q, r;
+
+ if (SCM_LIKELY (y > 0))
+ q = floor(x / y);
+ else if (SCM_LIKELY (y < 0))
+ q = ceil(x / y);
+ else if (y == 0)
+ scm_num_overflow (s_div_and_mod); /* or should we return a NaN? */
+ else
+ q = guile_NaN;
+ r = x - q * y;
+ return scm_values (scm_list_2 (scm_from_double (q),
+ scm_from_double (r)));
+}
+
+/* Compute exact div and mod the slow way.
+ We use this only if both arguments are exact,
+ and at least one of them is a fraction */
+static SCM
+scm_i_slow_exact_div_and_mod (SCM x, SCM y)
+{
+ SCM q, r;
+
+ if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+ SCM_WTA_DISPATCH_2 (g_div_and_mod, x, y, SCM_ARG1, s_div_and_mod);
+ else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+ SCM_WTA_DISPATCH_2 (g_div_and_mod, x, y, SCM_ARG2, s_div_and_mod);
+ else if (scm_is_true (scm_positive_p (y)))
+ q = scm_floor (scm_divide (x, y));
+ else if (scm_is_true (scm_negative_p (y)))
+ q = scm_ceiling (scm_divide (x, y));
+ else
+ scm_num_overflow (s_div_and_mod);
+ r = scm_difference (x, scm_product (q, y));
+ return scm_values (scm_list_2 (q, r));
+}
+
+static SCM scm_i_inexact_div0 (double x, double y);
+static SCM scm_i_bigint_div0 (SCM x, SCM y);
+static SCM scm_i_slow_exact_div0 (SCM x, SCM y);
+
+SCM_GPROC (s_div0, "div0", 2, 0, 0, scm_div0, g_div0);
+/* "Return q = @var{x} div0 @var{y}, where x = r + q*y,\n"
+ * "q is an integer and -abs(y/2) <= r < abs(y/2)."
+ * "@lisp\n"
+ * "(div0 123 10) @result{} 12\n"
+ * "(div0 123 -10) @result{} -12\n"
+ * "(div0 -123 10) @result{} -12\n"
+ * "(div0 -123 -10) @result{} 12\n"
+ * "@end lisp"
+ */
+SCM
+scm_div0 (SCM x, SCM y)
+{
+ if (SCM_LIKELY (SCM_I_INUMP (x)))
+ {
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
+ {
+ scm_t_inum yy = SCM_I_INUM (y);
+ if (SCM_UNLIKELY (yy == 0))
+ scm_num_overflow (s_div0);
+ else
+ {
+ scm_t_inum xx = SCM_I_INUM (x);
+ scm_t_inum qq = xx / yy;
+ scm_t_inum rr = xx - qq * yy;
+ if (SCM_LIKELY (xx > 0))
+ {
+ if (SCM_LIKELY (yy > 0))
+ {
+ if (rr >= (yy + 1) / 2)
+ qq++;
+ }
+ else
+ {
+ if (rr >= (1 - yy) / 2)
+ qq--;
+ }
+ }
+ else
+ {
+ if (SCM_LIKELY (yy > 0))
+ {
+ if (rr < -yy / 2)
+ qq--;
+ }
+ else
+ {
+ if (rr < yy / 2)
+ qq++;
+ }
+ }
+ return SCM_I_MAKINUM (qq);
+ }
+ }
+ else if (SCM_BIGP (y))
+ {
+ /* Pass a denormalized bignum version of x (even though it
+ can fit in a fixnum) to scm_i_bigint_div0 */
+ return scm_i_bigint_div0
+ (scm_i_long2big (SCM_I_INUM (x)), y);
+ }
+ else if (SCM_REALP (y))
+ return scm_i_inexact_div0 (SCM_I_INUM (x), SCM_REAL_VALUE (y));
+ else if (SCM_FRACTIONP (y))
+ return scm_i_slow_exact_div0 (x, y);
+ else
+ SCM_WTA_DISPATCH_2 (g_div0, x, y, SCM_ARG2, s_div0);
+ }
+ else if (SCM_BIGP (x))
+ {
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
+ {
+ scm_t_inum yy = SCM_I_INUM (y);
+ if (SCM_UNLIKELY (yy == 0))
+ scm_num_overflow (s_div0);
+ else
+ {
+ SCM q = scm_i_mkbig ();
+ scm_t_inum rr;
+ /* Arrange for rr to initially be non-positive,
+ because that simplifies the test to see
+ if it is within the needed bounds. */
+ if (yy > 0)
+ {
+ rr = - mpz_cdiv_q_ui (SCM_I_BIG_MPZ (q),
+ SCM_I_BIG_MPZ (x), yy);
+ scm_remember_upto_here_1 (x);
+ if (rr < -yy / 2)
+ mpz_sub_ui (SCM_I_BIG_MPZ (q),
+ SCM_I_BIG_MPZ (q), 1);
+ }
+ else
+ {
+ rr = - mpz_cdiv_q_ui (SCM_I_BIG_MPZ (q),
+ SCM_I_BIG_MPZ (x), -yy);
+ scm_remember_upto_here_1 (x);
+ mpz_neg (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (q));
+ if (rr < yy / 2)
+ mpz_add_ui (SCM_I_BIG_MPZ (q),
+ SCM_I_BIG_MPZ (q), 1);
+ }
+ return scm_i_normbig (q);
+ }
+ }
+ else if (SCM_BIGP (y))
+ return scm_i_bigint_div0 (x, y);
+ else if (SCM_REALP (y))
+ return scm_i_inexact_div0 (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+ else if (SCM_FRACTIONP (y))
+ return scm_i_slow_exact_div0 (x, y);
+ else
+ SCM_WTA_DISPATCH_2 (g_div0, x, y, SCM_ARG2, s_div0);
+ }
+ else if (SCM_REALP (x))
+ {
+ if (!(SCM_REALP (y) || SCM_I_INUMP (y) ||
+ SCM_BIGP (y) || SCM_FRACTIONP (y)))
+ SCM_WTA_DISPATCH_2 (g_div0, x, y, SCM_ARG2, s_div0);
+ else
+ return scm_i_inexact_div0 (SCM_REAL_VALUE (x), scm_to_double (y));
+ }
+ else if (SCM_FRACTIONP (x))
+ {
+ if (SCM_REALP (y))
+ return scm_i_inexact_div0 (scm_i_fraction2double (x),
+ SCM_REAL_VALUE (y));
+ else
+ return scm_i_slow_exact_div0 (x, y);
+ }
+ else
+ SCM_WTA_DISPATCH_2 (g_div0, x, y, SCM_ARG1, s_div0);
+}
+
+static SCM
+scm_i_inexact_div0 (double x, double y)
+{
+ if (SCM_LIKELY (y > 0))
+ return scm_from_double (floor(x / y + 0.5));
+ else if (SCM_LIKELY (y < 0))
+ return scm_from_double (ceil(x / y - 0.5));
+ else if (y == 0)
+ scm_num_overflow (s_div0); /* or should we return a NaN? */
+ else
+ return scm_nan ();
+}
+
+/* Assumes that both x and y are bigints, though
+ x might be able to fit into a fixnum. */
+static SCM
+scm_i_bigint_div0 (SCM x, SCM y)
+{
+ SCM q, r, min_r;
+
+ /* Note that x might be small enough to fit into a
+ fixnum, so we must not let it escape into the wild */
+ q = scm_i_mkbig ();
+ r = scm_i_mkbig ();
+
+ /* min_r will eventually become -abs(y)/2 */
+ min_r = scm_i_mkbig ();
+ mpz_tdiv_q_2exp (SCM_I_BIG_MPZ (min_r),
+ SCM_I_BIG_MPZ (y), 1);
+
+ /* Arrange for rr to initially be non-positive,
+ because that simplifies the test to see
+ if it is within the needed bounds. */
+ if (mpz_sgn (SCM_I_BIG_MPZ (y)) > 0)
+ {
+ mpz_cdiv_qr (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
+ scm_remember_upto_here_2 (x, y);
+ mpz_neg (SCM_I_BIG_MPZ (min_r), SCM_I_BIG_MPZ (min_r));
+ if (mpz_cmp (SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (min_r)) < 0)
+ mpz_sub_ui (SCM_I_BIG_MPZ (q),
+ SCM_I_BIG_MPZ (q), 1);
+ }
+ else
+ {
+ mpz_fdiv_qr (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
+ scm_remember_upto_here_2 (x, y);
+ if (mpz_cmp (SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (min_r)) < 0)
+ mpz_add_ui (SCM_I_BIG_MPZ (q),
+ SCM_I_BIG_MPZ (q), 1);
+ }
+ scm_remember_upto_here_2 (r, min_r);
+ return scm_i_normbig (q);
+}
+
+/* Compute exact div0 the slow way.
+ We use this only if both arguments are exact,
+ and at least one of them is a fraction */
+static SCM
+scm_i_slow_exact_div0 (SCM x, SCM y)
+{
+ if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+ SCM_WTA_DISPATCH_2 (g_div0, x, y, SCM_ARG1, s_div0);
+ else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+ SCM_WTA_DISPATCH_2 (g_div0, x, y, SCM_ARG2, s_div0);
+ else if (scm_is_true (scm_positive_p (y)))
+ return scm_floor (scm_sum (scm_divide (x, y), exactly_one_half));
+ else if (scm_is_true (scm_negative_p (y)))
+ return scm_ceiling (scm_difference (scm_divide (x, y), exactly_one_half));
+ else
+ scm_num_overflow (s_div0);
+}
+
+static SCM scm_i_inexact_mod0 (double x, double y);
+static SCM scm_i_bigint_mod0 (SCM x, SCM y);
+static SCM scm_i_slow_exact_mod0 (SCM x, SCM y);
+
+SCM_GPROC (s_mod0, "mod0", 2, 0, 0, scm_mod0, g_mod0);
+/* "Return r = @var{x} mod0 @var{y}, where x = r + q*y,\n"
+ * "q is an integer and -abs(y/2) <= r < abs(y/2)."
+ * "@lisp\n"
+ * "(mod0 123 10) @result{} 3\n"
+ * "(mod0 123 -10) @result{} 3\n"
+ * "(mod0 -123 10) @result{} -3\n"
+ * "(mod0 -123 -10) @result{} -3\n"
+ * "@end lisp"
+ */
+SCM
+scm_mod0 (SCM x, SCM y)
+{
+ if (SCM_LIKELY (SCM_I_INUMP (x)))
+ {
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
+ {
+ scm_t_inum yy = SCM_I_INUM (y);
+ if (SCM_UNLIKELY (yy == 0))
+ scm_num_overflow (s_mod0);
+ else
+ {
+ scm_t_inum xx = SCM_I_INUM (x);
+ scm_t_inum rr = xx % yy;
+ if (SCM_LIKELY (xx > 0))
+ {
+ if (SCM_LIKELY (yy > 0))
+ {
+ if (rr >= (yy + 1) / 2)
+ rr -= yy;
+ }
+ else
+ {
+ if (rr >= (1 - yy) / 2)
+ rr += yy;
+ }
+ }
+ else
+ {
+ if (SCM_LIKELY (yy > 0))
+ {
+ if (rr < -yy / 2)
+ rr += yy;
+ }
+ else
+ {
+ if (rr < yy / 2)
+ rr -= yy;
+ }
+ }
+ return SCM_I_MAKINUM (rr);
+ }
+ }
+ else if (SCM_BIGP (y))
+ {
+ /* Pass a denormalized bignum version of x (even though it
+ can fit in a fixnum) to scm_i_bigint_mod0 */
+ return scm_i_bigint_mod0
+ (scm_i_long2big (SCM_I_INUM (x)), y);
+ }
+ else if (SCM_REALP (y))
+ return scm_i_inexact_mod0 (SCM_I_INUM (x), SCM_REAL_VALUE (y));
+ else if (SCM_FRACTIONP (y))
+ return scm_i_slow_exact_mod0 (x, y);
+ else
+ SCM_WTA_DISPATCH_2 (g_mod0, x, y, SCM_ARG2, s_mod0);
+ }
+ else if (SCM_BIGP (x))
+ {
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
+ {
+ scm_t_inum yy = SCM_I_INUM (y);
+ if (SCM_UNLIKELY (yy == 0))
+ scm_num_overflow (s_mod0);
+ else
+ {
+ scm_t_inum rr;
+ /* Arrange for rr to initially be non-positive,
+ because that simplifies the test to see
+ if it is within the needed bounds. */
+ if (yy > 0)
+ {
+ rr = - mpz_cdiv_ui (SCM_I_BIG_MPZ (x), yy);
+ scm_remember_upto_here_1 (x);
+ if (rr < -yy / 2)
+ rr += yy;
+ }
+ else
+ {
+ rr = - mpz_cdiv_ui (SCM_I_BIG_MPZ (x), -yy);
+ scm_remember_upto_here_1 (x);
+ if (rr < yy / 2)
+ rr -= yy;
+ }
+ return SCM_I_MAKINUM (rr);
+ }
+ }
+ else if (SCM_BIGP (y))
+ return scm_i_bigint_mod0 (x, y);
+ else if (SCM_REALP (y))
+ return scm_i_inexact_mod0 (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+ else if (SCM_FRACTIONP (y))
+ return scm_i_slow_exact_mod0 (x, y);
+ else
+ SCM_WTA_DISPATCH_2 (g_mod0, x, y, SCM_ARG2, s_mod0);
+ }
+ else if (SCM_REALP (x))
+ {
+ if (!(SCM_REALP (y) || SCM_I_INUMP (y) ||
+ SCM_BIGP (y) || SCM_FRACTIONP (y)))
+ SCM_WTA_DISPATCH_2 (g_mod0, x, y, SCM_ARG2, s_mod0);
+ else
+ return scm_i_inexact_mod0 (SCM_REAL_VALUE (x), scm_to_double (y));
+ }
+ else if (SCM_FRACTIONP (x))
+ {
+ if (SCM_REALP (y))
+ return scm_i_inexact_mod0 (scm_i_fraction2double (x),
+ SCM_REAL_VALUE (y));
+ else
+ return scm_i_slow_exact_mod0 (x, y);
+ }
+ else
+ SCM_WTA_DISPATCH_2 (g_mod0, x, y, SCM_ARG1, s_mod0);
+}
+
+static SCM
+scm_i_inexact_mod0 (double x, double y)
+{
+ double q;
+
+ /* Although it would be more efficient to use fmod here, we can't
+ because it would in some cases produce results inconsistent with
+ scm_i_inexact_div0, such that x != r + q * y (not even close). In
+ particular, when x-y/2 is very close to a multiple of y, then r
+ might be either -abs(y/2) or abs(y/2)-epsilon, but those two cases
+ must correspond with different choices of q. If div0 chooses one
+ way and mod0 chooses the other, it would be bad. */
+ if (SCM_LIKELY (y > 0))
+ q = floor(x / y + 0.5);
+ else if (SCM_LIKELY (y < 0))
+ q = ceil(x / y - 0.5);
+ else if (y == 0)
+ scm_num_overflow (s_mod0); /* or should we return a NaN? */
+ else
+ return scm_nan ();
+ return scm_from_double (x - q * y);
+}
+
+/* Assumes that both x and y are bigints, though
+ x might be able to fit into a fixnum. */
+static SCM
+scm_i_bigint_mod0 (SCM x, SCM y)
+{
+ SCM r, min_r;
+
+ /* Note that x might be small enough to fit into a
+ fixnum, so we must not let it escape into the wild */
+ r = scm_i_mkbig ();
+
+ /* min_r will eventually become -abs(y)/2 */
+ min_r = scm_i_mkbig ();
+ mpz_tdiv_q_2exp (SCM_I_BIG_MPZ (min_r),
+ SCM_I_BIG_MPZ (y), 1);
+
+ /* Arrange for rr to initially be non-positive,
+ because that simplifies the test to see
+ if it is within the needed bounds. */
+ if (mpz_sgn (SCM_I_BIG_MPZ (y)) > 0)
+ {
+ mpz_cdiv_r (SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
+ mpz_neg (SCM_I_BIG_MPZ (min_r), SCM_I_BIG_MPZ (min_r));
+ if (mpz_cmp (SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (min_r)) < 0)
+ mpz_add (SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (y));
+ }
+ else
+ {
+ mpz_fdiv_r (SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
+ if (mpz_cmp (SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (min_r)) < 0)
+ mpz_sub (SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (y));
+ }
+ scm_remember_upto_here_2 (x, y);
+ return scm_i_normbig (r);
+}
+
+/* Compute exact mod0 the slow way: x-y*(x div0 y)
+ We use this only if both arguments are exact,
+ and at least one of them is a fraction */
+static SCM
+scm_i_slow_exact_mod0 (SCM x, SCM y)
+{
+ if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+ SCM_WTA_DISPATCH_2 (g_mod0, x, y, SCM_ARG1, s_mod0);
+ else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+ SCM_WTA_DISPATCH_2 (g_mod0, x, y, SCM_ARG2, s_mod0);
+ else if (scm_is_true (scm_positive_p (y)))
+ return scm_difference
+ (x, scm_product (y, scm_floor (scm_sum (scm_divide (x, y),
+ exactly_one_half))));
+ else if (scm_is_true (scm_negative_p (y)))
+ return scm_difference
+ (x, scm_product (y, scm_ceiling (scm_difference (scm_divide (x, y),
+ exactly_one_half))));
+ else
+ scm_num_overflow (s_mod0);
+}
+
+
+static SCM scm_i_inexact_div0_and_mod0 (double x, double y);
+static SCM scm_i_bigint_div0_and_mod0 (SCM x, SCM y);
+static SCM scm_i_slow_exact_div0_and_mod0 (SCM x, SCM y);
+
+SCM_GPROC (s_div0_and_mod0, "div0-and-mod0", 2, 0, 0,
+ scm_div0_and_mod0, g_div0_and_mod0);
+/* "Return q and r, where x = r + q*y,"
+ * "q is an integer and -abs(y/2) <= r < abs(y/2)."
+ * "@lisp\n"
+ * "(div0-and-mod0 123 10) @result{} 12 and 3\n"
+ * "(div0-and-mod0 123 -10) @result{} -12 and 3\n"
+ * "(div0-and-mod0 -123 10) @result{} -12 and -3\n"
+ * "(div0-and-mod0 -123 -10) @result{} 12 and -3\n"
+ * "@end lisp"
+ */
+SCM
+scm_div0_and_mod0 (SCM x, SCM y)
+{
+ if (SCM_LIKELY (SCM_I_INUMP (x)))
+ {
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
+ {
+ scm_t_inum yy = SCM_I_INUM (y);
+ if (SCM_UNLIKELY (yy == 0))
+ scm_num_overflow (s_div0_and_mod0);
+ else
+ {
+ scm_t_inum xx = SCM_I_INUM (x);
+ scm_t_inum qq = xx / yy;
+ scm_t_inum rr = xx - qq * yy;
+ if (SCM_LIKELY (xx > 0))
+ {
+ if (SCM_LIKELY (yy > 0))
+ {
+ if (rr >= (yy + 1) / 2)
+ { qq++; rr -= yy; }
+ }
+ else
+ {
+ if (rr >= (1 - yy) / 2)
+ { qq--; rr += yy; }
+ }
+ }
+ else
+ {
+ if (SCM_LIKELY (yy > 0))
+ {
+ if (rr < -yy / 2)
+ { qq--; rr += yy; }
+ }
+ else
+ {
+ if (rr < yy / 2)
+ { qq++; rr -= yy; }
+ }
+ }
+ return scm_values (scm_list_2 (SCM_I_MAKINUM (qq),
+ SCM_I_MAKINUM (rr)));
+ }
+ }
+ else if (SCM_BIGP (y))
+ {
+ /* Pass a denormalized bignum version of x (even though it
+ can fit in a fixnum) to scm_i_bigint_div0_and_mod0 */
+ return scm_i_bigint_div0_and_mod0
+ (scm_i_long2big (SCM_I_INUM (x)), y);
+ }
+ else if (SCM_REALP (y))
+ return scm_i_inexact_div0_and_mod0 (SCM_I_INUM (x),
+ SCM_REAL_VALUE (y));
+ else if (SCM_FRACTIONP (y))
+ return scm_i_slow_exact_div0_and_mod0 (x, y);
+ else
+ SCM_WTA_DISPATCH_2 (g_div0_and_mod0, x, y, SCM_ARG2, s_div0_and_mod0);
+ }
+ else if (SCM_BIGP (x))
+ {
+ if (SCM_LIKELY (SCM_I_INUMP (y)))
+ {
+ scm_t_inum yy = SCM_I_INUM (y);
+ if (SCM_UNLIKELY (yy == 0))
+ scm_num_overflow (s_div0_and_mod0);
+ else
+ {
+ SCM q = scm_i_mkbig ();
+ scm_t_inum rr;
+ /* Arrange for rr to initially be non-positive,
+ because that simplifies the test to see
+ if it is within the needed bounds. */
+ if (yy > 0)
+ {
+ rr = - mpz_cdiv_q_ui (SCM_I_BIG_MPZ (q),
+ SCM_I_BIG_MPZ (x), yy);
+ scm_remember_upto_here_1 (x);
+ if (rr < -yy / 2)
+ {
+ mpz_sub_ui (SCM_I_BIG_MPZ (q),
+ SCM_I_BIG_MPZ (q), 1);
+ rr += yy;
+ }
+ }
+ else
+ {
+ rr = - mpz_cdiv_q_ui (SCM_I_BIG_MPZ (q),
+ SCM_I_BIG_MPZ (x), -yy);
+ scm_remember_upto_here_1 (x);
+ mpz_neg (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (q));
+ if (rr < yy / 2)
+ {
+ mpz_add_ui (SCM_I_BIG_MPZ (q),
+ SCM_I_BIG_MPZ (q), 1);
+ rr -= yy;
+ }
+ }
+ return scm_values (scm_list_2 (scm_i_normbig (q),
+ SCM_I_MAKINUM (rr)));
+ }
+ }
+ else if (SCM_BIGP (y))
+ return scm_i_bigint_div0_and_mod0 (x, y);
+ else if (SCM_REALP (y))
+ return scm_i_inexact_div0_and_mod0 (scm_i_big2dbl (x),
+ SCM_REAL_VALUE (y));
+ else if (SCM_FRACTIONP (y))
+ return scm_i_slow_exact_div0_and_mod0 (x, y);
+ else
+ SCM_WTA_DISPATCH_2 (g_div0_and_mod0, x, y, SCM_ARG2, s_div0_and_mod0);
+ }
+ else if (SCM_REALP (x))
+ {
+ if (!(SCM_REALP (y) || SCM_I_INUMP (y) ||
+ SCM_BIGP (y) || SCM_FRACTIONP (y)))
+ SCM_WTA_DISPATCH_2 (g_div0_and_mod0, x, y, SCM_ARG2, s_div0_and_mod0);
+ else
+ return scm_i_inexact_div0_and_mod0 (SCM_REAL_VALUE (x),
+ scm_to_double (y));
+ }
+ else if (SCM_FRACTIONP (x))
+ {
+ if (SCM_REALP (y))
+ return scm_i_inexact_div0_and_mod0 (scm_i_fraction2double (x),
+ SCM_REAL_VALUE (y));
+ else
+ return scm_i_slow_exact_div0_and_mod0 (x, y);
+ }
+ else
+ SCM_WTA_DISPATCH_2 (g_div0_and_mod0, x, y, SCM_ARG1, s_div0_and_mod0);
+}
+
+static SCM
+scm_i_inexact_div0_and_mod0 (double x, double y)
+{
+ double q, r;
+
+ if (SCM_LIKELY (y > 0))
+ q = floor(x / y + 0.5);
+ else if (SCM_LIKELY (y < 0))
+ q = ceil(x / y - 0.5);
+ else if (y == 0)
+ scm_num_overflow (s_div0_and_mod0); /* or should we return a NaN? */
+ else
+ q = guile_NaN;
+ r = x - q * y;
+ return scm_values (scm_list_2 (scm_from_double (q),
+ scm_from_double (r)));
+}
+
+/* Assumes that both x and y are bigints, though
+ x might be able to fit into a fixnum. */
+static SCM
+scm_i_bigint_div0_and_mod0 (SCM x, SCM y)
+{
+ SCM q, r, min_r;
+
+ /* Note that x might be small enough to fit into a
+ fixnum, so we must not let it escape into the wild */
+ q = scm_i_mkbig ();
+ r = scm_i_mkbig ();
+
+ /* min_r will eventually become -abs(y/2) */
+ min_r = scm_i_mkbig ();
+ mpz_tdiv_q_2exp (SCM_I_BIG_MPZ (min_r),
+ SCM_I_BIG_MPZ (y), 1);
+
+ /* Arrange for rr to initially be non-positive,
+ because that simplifies the test to see
+ if it is within the needed bounds. */
+ if (mpz_sgn (SCM_I_BIG_MPZ (y)) > 0)
+ {
+ mpz_cdiv_qr (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
+ mpz_neg (SCM_I_BIG_MPZ (min_r), SCM_I_BIG_MPZ (min_r));
+ if (mpz_cmp (SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (min_r)) < 0)
+ {
+ mpz_sub_ui (SCM_I_BIG_MPZ (q),
+ SCM_I_BIG_MPZ (q), 1);
+ mpz_add (SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (y));
+ }
+ }
+ else
+ {
+ mpz_fdiv_qr (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
+ if (mpz_cmp (SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (min_r)) < 0)
+ {
+ mpz_add_ui (SCM_I_BIG_MPZ (q),
+ SCM_I_BIG_MPZ (q), 1);
+ mpz_sub (SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (r),
+ SCM_I_BIG_MPZ (y));
+ }
+ }
+ scm_remember_upto_here_2 (x, y);
+ return scm_values (scm_list_2 (scm_i_normbig (q),
+ scm_i_normbig (r)));
+}
+
+/* Compute exact div0 and mod0 the slow way.
+ We use this only if both arguments are exact,
+ and at least one of them is a fraction */
+static SCM
+scm_i_slow_exact_div0_and_mod0 (SCM x, SCM y)
+{
+ SCM q, r;
+
+ if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+ SCM_WTA_DISPATCH_2 (g_div0_and_mod0, x, y, SCM_ARG1, s_div0_and_mod0);
+ else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+ SCM_WTA_DISPATCH_2 (g_div0_and_mod0, x, y, SCM_ARG2, s_div0_and_mod0);
+ else if (scm_is_true (scm_positive_p (y)))
+ q = scm_floor (scm_sum (scm_divide (x, y), exactly_one_half));
+ else if (scm_is_true (scm_negative_p (y)))
+ q = scm_ceiling (scm_difference (scm_divide (x, y), exactly_one_half));
+ else
+ scm_num_overflow (s_div0_and_mod0);
+ r = scm_difference (x, scm_product (q, y));
+ return scm_values (scm_list_2 (q, r));
+}
+
+
SCM_PRIMITIVE_GENERIC (scm_i_gcd, "gcd", 0, 2, 1,
(SCM x, SCM y, SCM rest),
"Return the greatest common divisor of all parameter values.\n"
@@ -5356,8 +6526,6 @@ SCM_DEFINE (scm_truncate_number, "truncate", 1, 0, 0,
}
#undef FUNC_NAME
-static SCM exactly_one_half;
-
SCM_DEFINE (scm_round_number, "round", 1, 0, 0,
(SCM x),
"Round the number @var{x} towards the nearest integer. "
diff --git a/libguile/numbers.h b/libguile/numbers.h
index 740dc80..4cc6095 100644
--- a/libguile/numbers.h
+++ b/libguile/numbers.h
@@ -177,6 +177,12 @@ SCM_API SCM scm_abs (SCM x);
SCM_API SCM scm_quotient (SCM x, SCM y);
SCM_API SCM scm_remainder (SCM x, SCM y);
SCM_API SCM scm_modulo (SCM x, SCM y);
+SCM_API SCM scm_div (SCM x, SCM y);
+SCM_API SCM scm_mod (SCM x, SCM y);
+SCM_API SCM scm_div_and_mod (SCM x, SCM y);
+SCM_API SCM scm_div0 (SCM x, SCM y);
+SCM_API SCM scm_mod0 (SCM x, SCM y);
+SCM_API SCM scm_div0_and_mod0 (SCM x, SCM y);
SCM_API SCM scm_gcd (SCM x, SCM y);
SCM_API SCM scm_lcm (SCM n1, SCM n2);
SCM_API SCM scm_logand (SCM n1, SCM n2);
diff --git a/module/rnrs/arithmetic/fixnums.scm b/module/rnrs/arithmetic/fixnums.scm
index c1f3571..befbe9d 100644
--- a/module/rnrs/arithmetic/fixnums.scm
+++ b/module/rnrs/arithmetic/fixnums.scm
@@ -1,6 +1,6 @@
;;; fixnums.scm --- The R6RS fixnums arithmetic library
-;; Copyright (C) 2010 Free Software Foundation, Inc.
+;; Copyright (C) 2010, 2011 Free Software Foundation, Inc.
;;
;; This library is free software; you can redistribute it and/or
;; modify it under the terms of the GNU Lesser General Public
@@ -175,40 +175,33 @@
(define (fxdiv fx1 fx2)
(assert-fixnum fx1 fx2)
- (if (zero? fx2) (raise (make-assertion-violation)))
- (let ((r (div fx1 fx2))) r))
+ (div fx1 fx2))
(define (fxmod fx1 fx2)
(assert-fixnum fx1 fx2)
- (if (zero? fx2) (raise (make-assertion-violation)))
- (let ((r (mod fx1 fx2))) r))
+ (mod fx1 fx2))
(define (fxdiv-and-mod fx1 fx2)
(assert-fixnum fx1 fx2)
- (if (zero? fx2) (raise (make-assertion-violation)))
(div-and-mod fx1 fx2))
(define (fxdiv0 fx1 fx2)
(assert-fixnum fx1 fx2)
- (if (zero? fx2) (raise (make-assertion-violation)))
- (let ((r (div0 fx1 fx2))) r))
+ (div0 fx1 fx2))
(define (fxmod0 fx1 fx2)
(assert-fixnum fx1 fx2)
- (if (zero? fx2) (raise (make-assertion-violation)))
- (let ((r (mod0 fx1 fx2))) r))
+ (mod0 fx1 fx2))
(define (fxdiv0-and-mod0 fx1 fx2)
(assert-fixnum fx1 fx2)
- (if (zero? fx2) (raise (make-assertion-violation)))
- (call-with-values (lambda () (div0-and-mod0 fx1 fx2))
- (lambda (q r) (values q r))))
+ (div0-and-mod0 fx1 fx2))
(define (fx+/carry fx1 fx2 fx3)
(assert-fixnum fx1 fx2 fx3)
(let* ((s (+ fx1 fx2 fx3))
- (s0 (mod0 s (inexact->exact (expt 2 (fixnum-width)))))
- (s1 (div0 s (inexact->exact (expt 2 (fixnum-width))))))
+ (s0 (mod0 s (expt 2 (fixnum-width))))
+ (s1 (div0 s (expt 2 (fixnum-width)))))
(values s0 s1)))
(define (fx-/carry fx1 fx2 fx3)
diff --git a/module/rnrs/arithmetic/flonums.scm b/module/rnrs/arithmetic/flonums.scm
index 4fadbd0..b65c294 100644
--- a/module/rnrs/arithmetic/flonums.scm
+++ b/module/rnrs/arithmetic/flonums.scm
@@ -1,6 +1,6 @@
;;; flonums.scm --- The R6RS flonums arithmetic library
-;; Copyright (C) 2010 Free Software Foundation, Inc.
+;; Copyright (C) 2010, 2011 Free Software Foundation, Inc.
;;
;; This library is free software; you can redistribute it and/or
;; modify it under the terms of the GNU Lesser General Public
@@ -127,40 +127,27 @@
(define (fldiv-and-mod fl1 fl2)
(assert-iflonum fl1 fl2)
- (if (zero? fl2) (raise (make-assertion-violation)))
- (let ((fx1 (inexact->exact fl1))
- (fx2 (inexact->exact fl2)))
- (call-with-values (lambda () (div-and-mod fx1 fx2))
- (lambda (div mod) (values (exact->inexact div)
- (exact->inexact mod))))))
+ (div-and-mod fl1 fl2))
(define (fldiv fl1 fl2)
(assert-iflonum fl1 fl2)
- (if (zero? fl2) (raise (make-assertion-violation)))
- (let ((fx1 (inexact->exact fl1))
- (fx2 (inexact->exact fl2)))
- (exact->inexact (quotient fx1 fx2))))
+ (div fl1 fl2))
(define (flmod fl1 fl2)
(assert-iflonum fl1 fl2)
- (if (zero? fl2) (raise (make-assertion-violation)))
- (let ((fx1 (inexact->exact fl1))
- (fx2 (inexact->exact fl2)))
- (exact->inexact (modulo fx1 fx2))))
+ (mod fl1 fl2))
(define (fldiv0-and-mod0 fl1 fl2)
(assert-iflonum fl1 fl2)
- (if (zero? fl2) (raise (make-assertion-violation)))
- (let* ((fx1 (inexact->exact fl1))
- (fx2 (inexact->exact fl2)))
- (call-with-values (lambda () (div0-and-mod0 fx1 fx2))
- (lambda (q r) (values (real->flonum q) (real->flonum r))))))
+ (div0-and-mod0 fl1 fl2))
(define (fldiv0 fl1 fl2)
- (call-with-values (lambda () (fldiv0-and-mod0 fl1 fl2)) (lambda (q r) q)))
+ (assert-iflonum fl1 fl2)
+ (div0 fl1 fl2))
(define (flmod0 fl1 fl2)
- (call-with-values (lambda () (fldiv0-and-mod0 fl1 fl2)) (lambda (q r) r)))
+ (assert-iflonum fl1 fl2)
+ (mod0 fl1 fl2))
(define (flnumerator fl)
(assert-flonum fl)
diff --git a/module/rnrs/base.scm b/module/rnrs/base.scm
index 04a7e23..37c574a 100644
--- a/module/rnrs/base.scm
+++ b/module/rnrs/base.scm
@@ -74,8 +74,6 @@
syntax-rules identifier-syntax)
(import (rename (except (guile) error raise)
- (quotient div)
- (modulo mod)
(inf? infinite?)
(exact->inexact inexact)
(inexact->exact exact))
@@ -119,21 +117,6 @@
(define (vector-map proc . vecs)
(list->vector (apply map (cons proc (map vector->list vecs)))))
- (define (div-and-mod x y) (let ((q (div x y)) (r (mod x y))) (values q r)))
-
- (define (div0 x y)
- (call-with-values (lambda () (div0-and-mod0 x y)) (lambda (q r) q)))
-
- (define (mod0 x y)
- (call-with-values (lambda () (div0-and-mod0 x y)) (lambda (q r) r)))
-
- (define (div0-and-mod0 x y)
- (call-with-values (lambda () (div-and-mod x y))
- (lambda (q r)
- (cond ((< r (abs (/ y 2))) (values q r))
- ((negative? y) (values (- q 1) (+ r y)))
- (else (values (+ q 1) (+ r y)))))))
-
(define raise
(@ (rnrs exceptions) raise))
(define condition
diff --git a/test-suite/tests/numbers.test b/test-suite/tests/numbers.test
index 36e3128..c89b98a 100644
--- a/test-suite/tests/numbers.test
+++ b/test-suite/tests/numbers.test
@@ -17,7 +17,8 @@
(define-module (test-suite test-numbers)
#:use-module (test-suite lib)
- #:use-module (ice-9 documentation))
+ #:use-module (ice-9 documentation)
+ #:use-module (srfi srfi-11)) ; let-values
;;;
;;; miscellaneous
@@ -92,6 +93,35 @@
(negative? obj)
(inf? obj)))
+;;
+;; Tolerance used by test-eqv? for inexact numbers.
+;;
+(define test-epsilon 1e-10)
+
+;;
+;; Like eqv?, except that inexact finite numbers need only be within
+;; test-epsilon (1e-10) to be considered equal. An exception is made
+;; for zeroes, however. If X is zero, then it is tested using eqv?
+;; without any allowance for imprecision. In particular, 0.0 is
+;; considered distinct from -0.0. For non-real complex numbers,
+;; each component is tested according to these rules. The intent
+;; is that the known-correct value will be the first parameter.
+;;
+(define (test-eqv? x y)
+ (cond ((real? x)
+ (and (real? y) (test-real-eqv? x y)))
+ ((complex? x)
+ (and (not (real? y))
+ (test-real-eqv? (real-part x) (real-part y))
+ (test-real-eqv? (imag-part x) (imag-part y))))
+ (else (eqv? x y))))
+
+;; Auxiliary predicate used by test-eqv?
+(define (test-real-eqv? x y)
+ (cond ((or (exact? x) (zero? x) (nan? x) (inf? x))
+ (eqv? x y))
+ (else (and (inexact? y) (> test-epsilon (abs (- x y)))))))
+
(define const-e 2.7182818284590452354)
(define const-e^2 7.3890560989306502274)
(define const-1/e 0.3678794411714423215)
@@ -3480,3 +3510,137 @@
(pass-if "-100i swings back to 45deg down"
(eqv-loosely? +7.071-7.071i (sqrt -100.0i))))
+;;;
+;;; div
+;;; mod
+;;; div-and-mod
+;;; div0
+;;; mod0
+;;; div0-and-mod0
+;;;
+
+(with-test-prefix "Number-theoretic division"
+
+ ;; Tests that (lo <= x < hi),
+ ;; but allowing for imprecision
+ ;; if x is inexact.
+ (define (test-within-range? lo hi x)
+ (if (exact? x)
+ (and (<= lo x) (< x hi))
+ (let ((lo (- lo test-epsilon))
+ (hi (+ hi test-epsilon)))
+ (<= lo x hi))))
+
+ (define (safe-div x y)
+ (cond ((not (and (real? x) (real? y))) (throw 'wrong-type-arg))
+ ((zero? y) (throw 'divide-by-zero))
+ ((nan? y) (nan))
+ ((positive? y) (floor (/ x y)))
+ ((negative? y) (ceiling (/ x y)))
+ (else (throw 'unknown-problem))))
+
+ (define (safe-mod x y)
+ (- x (* y (safe-div x y))))
+
+ (define (safe-div-and-mod x y)
+ (let ((q (safe-div x y))
+ (r (safe-mod x y)))
+ (if (not (and (eq? (exact? q) (exact? r))
+ (eq? (exact? q) (and (exact? x) (exact? y)))
+ (test-real-eqv? r (- x (* q y)))
+ (or (and (integer? q)
+ (test-within-range? 0 (abs y) r))
+ (not (finite? x))
+ (not (finite? y)))))
+ (throw 'safe-div-and-mod-is-broken (list x y q r))
+ (values q r))))
+
+ (define (safe-div0 x y)
+ (cond ((not (and (real? x) (real? y))) (throw 'wrong-type-arg))
+ ((zero? y) (throw 'divide-by-zero))
+ ((nan? y) (nan))
+ ((positive? y) (floor (+ 1/2 (/ x y))))
+ ((negative? y) (ceiling (+ -1/2 (/ x y))))
+ (else (throw 'unknown-problem))))
+
+ (define (safe-mod0 x y)
+ (- x (* y (safe-div0 x y))))
+
+ (define (safe-div0-and-mod0 x y)
+ (let ((q (safe-div0 x y))
+ (r (safe-mod0 x y)))
+ (if (not (and (eq? (exact? q) (exact? r))
+ (eq? (exact? q) (and (exact? x) (exact? y)))
+ (test-real-eqv? r (- x (* q y)))
+ (or (and (integer? q)
+ (test-within-range? (* -1/2 (abs y))
+ (* +1/2 (abs y))
+ r))
+ (not (finite? x))
+ (not (finite? y)))))
+ (throw 'safe-div0-and-mod0-is-broken (list x y q r))
+ (values q r))))
+
+ (define test-numerators
+ (list 123 125 127 130 3 5 10 123.2 125.0
+ -123 -125 -127 -130 -3 -5 -10 -123.2 -125.0
+ 127.2 130.0 123/7 125/7 127/7 130/7
+ -127.2 -130.0 -123/7 -125/7 -127/7 -130/7
+ 0 +0.0 -0.0 +inf.0 -inf.0 +nan.0
+ most-negative-fixnum (1+ most-positive-fixnum) (1- most-negative-fixnum)
+ (* 123 (+ 1 most-positive-fixnum)) (* 125 (+ 1 most-positive-fixnum)) (* 127 (+ 1 most-positive-fixnum))
+ (* 130 (+ 1 most-positive-fixnum)) (* 3 (+ 1 most-positive-fixnum)) (* 5 (+ 1 most-positive-fixnum))
+ (* 10 (+ 1 most-positive-fixnum))
+ (* -123 (+ 1 most-positive-fixnum)) (* -125 (+ 1 most-positive-fixnum)) (* -127 (+ 1 most-positive-fixnum))
+ (* -130 (+ 1 most-positive-fixnum)) (* -3 (+ 1 most-positive-fixnum)) (* -5 (+ 1 most-positive-fixnum))
+ (* -10 (+ 1 most-positive-fixnum))
+ (* 123 (+ 2 most-positive-fixnum)) (* 125 (+ 2 most-positive-fixnum)) (* 127 (+ 2 most-positive-fixnum))
+ (* 130 (+ 2 most-positive-fixnum)) (* 3 (+ 2 most-positive-fixnum)) (* 5 (+ 2 most-positive-fixnum))
+ (* 10 (+ 2 most-positive-fixnum))
+ (* -123 (+ 2 most-positive-fixnum)) (* -125 (+ 2 most-positive-fixnum)) (* -127 (+ 2 most-positive-fixnum))
+ (* -130 (+ 2 most-positive-fixnum)) (* -3 (+ 2 most-positive-fixnum)) (* -5 (+ 2 most-positive-fixnum))
+ (* -10 (+ 2 most-positive-fixnum))))
+
+ (define test-denominators
+ (list 10 5 10/7 127/2 10.0 63.5
+ -10 -5 -10/7 -127/2 -10.0 -63.5
+ +inf.0 -inf.0 +nan.0 most-negative-fixnum
+ (+ 1 most-positive-fixnum) (+ -1 most-negative-fixnum)
+ (+ 2 most-positive-fixnum) (+ -2 most-negative-fixnum)))
+
+ (define (do-tests-1 op-name real-op safe-op)
+ (for-each (lambda (d)
+ (for-each (lambda (n)
+ (run-test (list op-name n d) #t
+ (lambda ()
+ (test-eqv? (real-op n d)
+ (safe-op n d)))))
+ test-numerators))
+ test-denominators))
+
+ (define (do-tests-2 op-name real-op safe-op)
+ (for-each (lambda (d)
+ (for-each (lambda (n)
+ (run-test (list op-name n d) #t
+ (lambda ()
+ (let-values
+ (((q r) (safe-op n d))
+ ((q1 r1) (real-op n d)))
+ (and (test-eqv? q q1)
+ (test-eqv? r r1))))))
+ test-numerators))
+ test-denominators))
+
+ (with-test-prefix "div" (do-tests-1 'div div safe-div))
+ (with-test-prefix "mod" (do-tests-1 'mod mod safe-mod))
+ (with-test-prefix "div-and-mod"
+ (do-tests-2 'div-and-mod
+ div-and-mod
+ safe-div-and-mod))
+
+ (with-test-prefix "div0" (do-tests-1 'div0 div0 safe-div0))
+ (with-test-prefix "mod0" (do-tests-1 'mod0 mod0 safe-mod0))
+ (with-test-prefix "div0-and-mod0"
+ (do-tests-2 'div0-and-mod0
+ div0-and-mod0
+ safe-div0-and-mod0)))
diff --git a/test-suite/tests/r6rs-arithmetic-fixnums.test b/test-suite/tests/r6rs-arithmetic-fixnums.test
index fed72eb..d39d544 100644
--- a/test-suite/tests/r6rs-arithmetic-fixnums.test
+++ b/test-suite/tests/r6rs-arithmetic-fixnums.test
@@ -1,6 +1,6 @@
;;; arithmetic-fixnums.test --- Test suite for R6RS (rnrs arithmetic bitwise)
-;; Copyright (C) 2010 Free Software Foundation, Inc.
+;; Copyright (C) 2010, 2011 Free Software Foundation, Inc.
;;
;; This library is free software; you can redistribute it and/or
;; modify it under the terms of the GNU Lesser General Public
@@ -121,32 +121,25 @@
(pass-if "simple"
(call-with-values (lambda () (fxdiv-and-mod 123 10))
(lambda (d m)
- (or (and (fx=? d 12) (fx=? m 3))
- (throw 'unresolved))))))
+ (and (fx=? d 12) (fx=? m 3))))))
-(with-test-prefix "fxdiv"
- (pass-if "simple" (or (fx=? (fxdiv -123 10) -13) (throw 'unresolved))))
-
-(with-test-prefix "fxmod"
- (pass-if "simple" (or (fx=? (fxmod -123 10) 7) (throw 'unresolved))))
+(with-test-prefix "fxdiv" (pass-if "simple" (fx=? (fxdiv -123 10) -13)))
+(with-test-prefix "fxmod" (pass-if "simple" (fx=? (fxmod -123 10) 7)))
(with-test-prefix "fxdiv0-and-mod0"
(pass-if "simple"
(call-with-values (lambda () (fxdiv0-and-mod0 -123 10))
(lambda (d m)
- (or (and (fx=? d 12) (fx=? m -3))
- (throw 'unresolved))))))
-
-(with-test-prefix "fxdiv0"
- (pass-if "simple" (or (fx=? (fxdiv0 -123 10) 12) (throw 'unresolved))))
+ (and (fx=? d -12) (fx=? m -3))))))
-(with-test-prefix "fxmod0"
- (pass-if "simple" (or (fx=? (fxmod0 -123 10) -3) (throw 'unresolved))))
+(with-test-prefix "fxdiv0" (pass-if "simple" (fx=? (fxdiv0 -123 10) -12)))
+(with-test-prefix "fxmod0" (pass-if "simple" (fx=? (fxmod0 -123 10) -3)))
;; Without working div and mod implementations and without any example results
;; from the spec, I have no idea what the results of these functions should
;; be. -juliang
+;; UPDATE: div and mod implementations are now working properly -mhw
(with-test-prefix "fx+/carry" (pass-if "simple" (throw 'unresolved)))
diff --git a/test-suite/tests/r6rs-arithmetic-flonums.test b/test-suite/tests/r6rs-arithmetic-flonums.test
index 873447b..af9dbbf 100644
--- a/test-suite/tests/r6rs-arithmetic-flonums.test
+++ b/test-suite/tests/r6rs-arithmetic-flonums.test
@@ -1,6 +1,6 @@
;;; arithmetic-flonums.test --- Test suite for R6RS (rnrs arithmetic flonums)
-;; Copyright (C) 2010 Free Software Foundation, Inc.
+;; Copyright (C) 2010, 2011 Free Software Foundation, Inc.
;;
;; This library is free software; you can redistribute it and/or
;; modify it under the terms of the GNU Lesser General Public
@@ -195,14 +195,13 @@
(pass-if "simple"
(call-with-values (lambda () (fldiv0-and-mod0 -123.0 10.0))
(lambda (div mod)
- (or (and (fl=? div -12.0) (fl=? mod -3.0))
- (throw 'unresolved))))))
+ (and (fl=? div -12.0) (fl=? mod -3.0))))))
(with-test-prefix "fldiv0"
- (pass-if "simple" (or (fl=? (fldiv0 -123.0 10.0) -12.0) (throw 'unresolved))))
+ (pass-if "simple" (fl=? (fldiv0 -123.0 10.0) -12.0)))
(with-test-prefix "flmod0"
- (pass-if "simple" (or (fl=? (flmod0 -123.0 10.0) -3.0) (throw 'unresolved))))
+ (pass-if "simple" (fl=? (flmod0 -123.0 10.0) -3.0)))
(with-test-prefix "flnumerator"
(pass-if "simple" (fl=? (flnumerator 0.5) 1.0))
--
1.5.6.5
next prev parent reply other threads:[~2011-01-29 8:20 UTC|newest]
Thread overview: 24+ messages / expand[flat|nested] mbox.gz Atom feed top
2011-01-26 16:32 [PATCH] First batch of numerics changes Mark H Weaver
2011-01-26 18:07 ` Mark H Weaver
2011-01-26 22:46 ` Mark H Weaver
2011-01-27 22:06 ` Mark H Weaver
2011-01-28 12:19 ` Andy Wingo
2011-01-29 0:05 ` Mark H Weaver
2011-01-29 11:29 ` Andy Wingo
2011-01-27 22:32 ` Mark H Weaver
2011-01-28 13:46 ` Andy Wingo
2011-01-28 14:44 ` Noah Lavine
2011-01-28 15:55 ` Andy Wingo
2011-01-29 8:20 ` Mark H Weaver [this message]
2011-01-29 17:42 ` Andy Wingo
2011-01-29 20:20 ` Mark H Weaver
2011-01-30 11:48 ` Andy Wingo
2011-01-29 17:50 ` Andy Wingo
2011-01-29 20:36 ` Mark H Weaver
2011-01-29 22:24 ` Mark H Weaver
2011-01-30 6:02 ` Commentary: R6RS div0-and-mod0 vs Taylor's `round/' Mark H Weaver
2011-01-30 11:50 ` [PATCH] First batch of numerics changes Andy Wingo
2011-01-30 12:12 ` Andy Wingo
2011-01-30 16:33 ` Mark H Weaver
2011-01-28 11:41 ` Andy Wingo
2011-01-28 23:36 ` Mark H Weaver
Reply instructions:
You may reply publicly to this message via plain-text email
using any one of the following methods:
* Save the following mbox file, import it into your mail client,
and reply-to-all from there: mbox
Avoid top-posting and favor interleaved quoting:
https://en.wikipedia.org/wiki/Posting_style#Interleaved_style
List information: https://www.gnu.org/software/guile/
* Reply using the --to, --cc, and --in-reply-to
switches of git-send-email(1):
git send-email \
--in-reply-to=87wrlo2k9z.fsf@yeeloong.netris.org \
--to=mhw@netris.org \
--cc=guile-devel@gnu.org \
--cc=wingo@pobox.com \
/path/to/YOUR_REPLY
https://kernel.org/pub/software/scm/git/docs/git-send-email.html
* If your mail client supports setting the In-Reply-To header
via mailto: links, try the mailto: link
Be sure your reply has a Subject: header at the top and a blank line
before the message body.
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).