unofficial mirror of guile-devel@gnu.org 
 help / color / mirror / Atom feed
[parent not found: <595618.30923.qm@web114107.mail.gq1.yahoo.com>]
* [PATCH] New division operators, and optimization for fractions
@ 2011-02-10 23:42 Mark H Weaver
  2011-02-12 11:55 ` Andy Wingo
  2011-02-15 10:38 ` Ludovic Courtès
  0 siblings, 2 replies; 12+ messages in thread
From: Mark H Weaver @ 2011-02-10 23:42 UTC (permalink / raw)
  To: guile-devel

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

Hello all,

Here are three more patches.  The first adds fast implementations of the
remaining number-theoretic division operators as described in Taylor
Campbell's proposal: floor/, ceiling/, truncate/, round/, as well as the
single-valued variants.  The third patch optimizes the case where both
arguments are exact and at least one is a fraction.

The second patch is needed by the third.  It adds an internal C function
that extracts two values from a values object.  I needed it because my
optimization for fractions works by calling the same division operator
recursively on a subproblem involving only integers, and I need a way to
extract the two values.  Initially, I began to make internal versions of
the divide functions that returned their two values via parameters of
type (SCM *), but then I realized that in the cases where I dispatch
into GOOPS, I'd still need to extract from a values object.

If I've overlooked some existing mechanism for extracting multiple
values from C, or if you have a better idea, please let me know and I'll
rework the patch.

Note that these patches should be applied after the ones in my recent
post entitled "Miscellaneous fixes and improvements".

    Best,
     Mark



[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #2: Add four new sets of fast quotient and remainder operators --]
[-- Type: text/x-diff, Size: 87977 bytes --]

From c5cb51d82ed23c557721a6a68a38a1deb060c0ca Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Thu, 10 Feb 2011 17:04:50 -0500
Subject: [PATCH] Add four new sets of fast quotient and remainder operators

* libguile/numbers.c (scm_floor_divide, scm_floor_quotient,
  scm_floor_remainder, scm_ceiling_divide, scm_ceiling_quotient,
  scm_ceiling_remainder, scm_truncate_divide, scm_truncate_quotient,
  scm_truncate_remainder, scm_round_divide, scm_round_quotient,
  scm_round_remainder): New extensible procedures `floor/',
  `floor-quotient', `floor-remainder', `ceiling/', `ceiling-quotient',
  `ceiling-remainder', `truncate/', `truncate-quotient',
  `truncate-remainder', `round/', `round-quotient', and
  `round-remainder'.

  (scm_round_number, scm_floor, scm_ceiling): When applied to a
  fraction, use scm_round_quotient, scm_floor_quotient and
  scm_ceiling_quotient, respectively, for improved efficiency.

* libguile/numbers.h: Add function prototypes.

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

* doc/ref/api-data.texi (Arithmetic): Add documentation.

* NEWS: Add NEWS entries for new division operators.
---
 NEWS                          |   28 +
 doc/ref/api-data.texi         |  155 +++-
 libguile/numbers.c            | 2222 ++++++++++++++++++++++++++++++++++++++++-
 libguile/numbers.h            |   12 +
 test-suite/tests/numbers.test |   75 ++-
 5 files changed, 2451 insertions(+), 41 deletions(-)

diff --git a/NEWS b/NEWS
index 7912259..0b47121 100644
--- a/NEWS
+++ b/NEWS
@@ -15,6 +15,34 @@ Changes since the 1.9.15 prerelease:
 `define-once' is like Lisp's `defvar': it creates a toplevel binding,
 but only if one does not exist already.
 
+** Added four new sets of fast quotient and remainder operators
+
+Added four new sets of fast quotient and remainder operators with
+different semantics than the R5RS operators.  They support not only
+integers, but all reals, including exact rationals and inexact
+floating point numbers.
+
+These procedures accept two real numbers N and D, where the divisor D
+must be non-zero.  Each set of operators computes an integer quotient
+Q and a real remainder R such that N = Q*D + R and |R| < |D|.  They
+differ only in how N/D is rounded to produce Q.
+
+`floor-quotient' and `floor-remainder' compute Q and R, respectively,
+where Q has been rounded toward negative infinity.  `floor/' returns
+both Q and R, and is more efficient than computing each separately.
+Note that when applied to integers, `floor-remainder' is equivalent to
+the R5RS integer-only `modulo' operator.  `ceiling-quotient',
+`ceiling-remainder', and `ceiling/' are similar except that Q is
+rounded toward positive infinity.
+
+For `truncate-quotient', `truncate-remainder', and `truncate/', Q is
+rounded toward zero.  Note that when applied to integers,
+`truncate-quotient' and `truncate-remainder' are equivalent to the
+R5RS integer-only operators `quotient' and `remainder'.
+
+For `round-quotient', `round-remainder', and `round/', Q is rounded to
+the nearest integer, with ties going to the nearest even integer.
+
 ** Improved exactness handling for complex number parsing
 
 When parsing non-real complex numbers, exactness specifiers are now
diff --git a/doc/ref/api-data.texi b/doc/ref/api-data.texi
index fd2e7ee..7604c2c 100644
--- a/doc/ref/api-data.texi
+++ b/doc/ref/api-data.texi
@@ -907,7 +907,7 @@ sign as @var{n}.  In all cases quotient and remainder satisfy
 (remainder -13 4) @result{} -1
 @end lisp
 
-See also @code{euclidean-quotient}, @code{euclidean-remainder} and
+See also @code{truncate-quotient}, @code{truncate-remainder} and
 related operations in @ref{Arithmetic}.
 @end deffn
 
@@ -924,7 +924,7 @@ sign as @var{d}.
 (modulo -13 -4) @result{} -1
 @end lisp
 
-See also @code{euclidean-quotient}, @code{euclidean-remainder} and
+See also @code{floor-quotient}, @code{floor-remainder} and
 related operations in @ref{Arithmetic}.
 @end deffn
 
@@ -1148,9 +1148,21 @@ Returns the magnitude or angle of @var{z} as a @code{double}.
 @rnindex euclidean/
 @rnindex euclidean-quotient
 @rnindex euclidean-remainder
+@rnindex floor/
+@rnindex floor-quotient
+@rnindex floor-remainder
+@rnindex ceiling/
+@rnindex ceiling-quotient
+@rnindex ceiling-remainder
+@rnindex truncate/
+@rnindex truncate-quotient
+@rnindex truncate-remainder
 @rnindex centered/
 @rnindex centered-quotient
 @rnindex centered-remainder
+@rnindex round/
+@rnindex round-quotient
+@rnindex round-remainder
 
 The C arithmetic functions below always takes two arguments, while the
 Scheme functions can take an arbitrary number.  When you need to
@@ -1260,7 +1272,7 @@ These procedures accept two real numbers @var{x} and @var{y}, where the
 divisor @var{y} must be non-zero.  @code{euclidean-quotient} returns the
 integer @var{q} and @code{euclidean-remainder} returns the real number
 @var{r} such that @math{@var{x} = @var{q}*@var{y} + @var{r}} and
-@math{0 <= @var{r} < abs(@var{y})}.  @code{euclidean/} returns both @var{q} and
+@math{0 <= @var{r} < |@var{y}|}.  @code{euclidean/} returns both @var{q} and
 @var{r}, and is more efficient than computing each separately.  Note
 that when @math{@var{y} > 0}, @code{euclidean-quotient} returns
 @math{floor(@var{x}/@var{y})}, otherwise it returns
@@ -1281,6 +1293,93 @@ Note that these operators are equivalent to the R6RS operators
 @end lisp
 @end deffn
 
+@deffn {Scheme Procedure} floor/ x y
+@deffnx {Scheme Procedure} floor-quotient x y
+@deffnx {Scheme Procedure} floor-remainder x y
+@deffnx {C Function} scm_floor_divide (x y)
+@deffnx {C Function} scm_floor_quotient (x y)
+@deffnx {C Function} scm_floor_remainder (x y)
+These procedures accept two real numbers @var{x} and @var{y}, where the
+divisor @var{y} must be non-zero.  @code{floor-quotient} returns the
+integer @var{q} and @code{floor-remainder} returns the real number
+@var{r} such that @math{@var{q} = floor(@var{x}/@var{y})} and
+@math{@var{x} = @var{q}*@var{y} + @var{r}}.  @code{floor/} returns
+both @var{q} and @var{r}, and is more efficient than computing each
+separately.  Note that @var{r}, if non-zero, will have the same sign
+as @var{y}.
+
+When @var{x} and @var{y} are integers, @code{floor-quotient} is
+equivalent to the R5RS integer-only operator @code{modulo}.
+
+@lisp
+(floor-quotient 123 10) @result{} 12
+(floor-remainder 123 10) @result{} 3
+(floor/ 123 10) @result{} 12 and 3
+(floor/ 123 -10) @result{} -13 and -7
+(floor/ -123 10) @result{} -13 and 7
+(floor/ -123 -10) @result{} 12 and -3
+(floor/ -123.2 -63.5) @result{} 1.0 and -59.7
+(floor/ 16/3 -10/7) @result{} -4 and -8/21
+@end lisp
+@end deffn
+
+@deffn {Scheme Procedure} ceiling/ x y
+@deffnx {Scheme Procedure} ceiling-quotient x y
+@deffnx {Scheme Procedure} ceiling-remainder x y
+@deffnx {C Function} scm_ceiling_divide (x y)
+@deffnx {C Function} scm_ceiling_quotient (x y)
+@deffnx {C Function} scm_ceiling_remainder (x y)
+These procedures accept two real numbers @var{x} and @var{y}, where the
+divisor @var{y} must be non-zero.  @code{ceiling-quotient} returns the
+integer @var{q} and @code{ceiling-remainder} returns the real number
+@var{r} such that @math{@var{q} = ceiling(@var{x}/@var{y})} and
+@math{@var{x} = @var{q}*@var{y} + @var{r}}.  @code{ceiling/} returns
+both @var{q} and @var{r}, and is more efficient than computing each
+separately.  Note that @var{r}, if non-zero, will have the opposite sign
+of @var{y}.
+
+@lisp
+(ceiling-quotient 123 10) @result{} 13
+(ceiling-remainder 123 10) @result{} -7
+(ceiling/ 123 10) @result{} 13 and -7
+(ceiling/ 123 -10) @result{} -12 and 3
+(ceiling/ -123 10) @result{} -12 and -3
+(ceiling/ -123 -10) @result{} 13 and 7
+(ceiling/ -123.2 -63.5) @result{} 2.0 and 3.8
+(ceiling/ 16/3 -10/7) @result{} -3 and 22/21
+@end lisp
+@end deffn
+
+@deffn {Scheme Procedure} truncate/ x y
+@deffnx {Scheme Procedure} truncate-quotient x y
+@deffnx {Scheme Procedure} truncate-remainder x y
+@deffnx {C Function} scm_truncate_divide (x y)
+@deffnx {C Function} scm_truncate_quotient (x y)
+@deffnx {C Function} scm_truncate_remainder (x y)
+These procedures accept two real numbers @var{x} and @var{y}, where the
+divisor @var{y} must be non-zero.  @code{truncate-quotient} returns the
+integer @var{q} and @code{truncate-remainder} returns the real number
+@var{r} such that @var{q} is @math{@var{x}/@var{y}} rounded toward zero,
+and @math{@var{x} = @var{q}*@var{y} + @var{r}}.  @code{truncate/} returns
+both @var{q} and @var{r}, and is more efficient than computing each
+separately.  Note that @var{r}, if non-zero, will have the same sign
+as @var{x}.
+
+When @var{x} and @var{y} are integers, these operators are equivalent to
+the R5RS integer-only operators @code{quotient} and @code{remainder}.
+
+@lisp
+(truncate-quotient 123 10) @result{} 12
+(truncate-remainder 123 10) @result{} 3
+(truncate/ 123 10) @result{} 12 and 3
+(truncate/ 123 -10) @result{} -12 and 3
+(truncate/ -123 10) @result{} -12 and -3
+(truncate/ -123 -10) @result{} 12 and -3
+(truncate/ -123.2 -63.5) @result{} 1.0 and -59.7
+(truncate/ 16/3 -10/7) @result{} -3 and 22/21
+@end lisp
+@end deffn
+
 @deffn {Scheme Procedure} centered/ x y
 @deffnx {Scheme Procedure} centered-quotient x y
 @deffnx {Scheme Procedure} centered-remainder x y
@@ -1291,7 +1390,7 @@ These procedures accept two real numbers @var{x} and @var{y}, where the
 divisor @var{y} must be non-zero.  @code{centered-quotient} returns the
 integer @var{q} and @code{centered-remainder} returns the real number
 @var{r} such that @math{@var{x} = @var{q}*@var{y} + @var{r}} and
-@math{-abs(@var{y}/2) <= @var{r} < abs(@var{y}/2)}.  @code{centered/}
+@math{-|@var{y}/2| <= @var{r} < |@var{y}/2|}.  @code{centered/}
 returns both @var{q} and @var{r}, and is more efficient than computing
 each separately.
 
@@ -1300,7 +1399,8 @@ rounded to the nearest integer.  When @math{@var{x}/@var{y}} lies
 exactly half-way between two integers, the tie is broken according to
 the sign of @var{y}.  If @math{@var{y} > 0}, ties are rounded toward
 positive infinity, otherwise they are rounded toward negative infinity.
-This is a consequence of the requirement that @math{-abs(@var{y}/2) <= @var{r} < abs(@var{y}/2)}.
+This is a consequence of the requirement that
+@math{-|@var{y}/2| <= @var{r} < |@var{y}/2|}.
 
 Note that these operators are equivalent to the R6RS operators
 @code{div0}, @code{mod0}, and @code{div0-and-mod0}.
@@ -1312,11 +1412,56 @@ Note that these operators are equivalent to the R6RS operators
 (centered/ 123 -10) @result{} -12 and 3
 (centered/ -123 10) @result{} -12 and -3
 (centered/ -123 -10) @result{} 12 and -3
+(centered/ 125 10) @result{} 13 and -5
+(centered/ 127 10) @result{} 13 and -3
+(centered/ 135 10) @result{} 14 and -5
 (centered/ -123.2 -63.5) @result{} 2.0 and 3.8
 (centered/ 16/3 -10/7) @result{} -4 and -8/21
 @end lisp
 @end deffn
 
+@deffn {Scheme Procedure} round/ x y
+@deffnx {Scheme Procedure} round-quotient x y
+@deffnx {Scheme Procedure} round-remainder x y
+@deffnx {C Function} scm_round_divide (x y)
+@deffnx {C Function} scm_round_quotient (x y)
+@deffnx {C Function} scm_round_remainder (x y)
+These procedures accept two real numbers @var{x} and @var{y}, where the
+divisor @var{y} must be non-zero.  @code{round-quotient} returns the
+integer @var{q} and @code{round-remainder} returns the real number
+@var{r} such that @math{@var{x} = @var{q}*@var{y} + @var{r}} and
+@var{q} is @math{@var{x}/@var{y}} rounded to the nearest integer,
+with ties going to the nearest even integer.  @code{round/}
+returns both @var{q} and @var{r}, and is more efficient than computing
+each separately.
+
+Note that @code{round/} and @code{centered/} are almost equivalent, but
+their behavior differs when @math{@var{x}/@var{y}} lies exactly half-way
+between two integers.  In this case, @code{round/} chooses the nearest
+even integer, whereas @code{centered/} chooses in such a way to satisfy
+the constraint @math{-|@var{y}/2| <= @var{r} < |@var{y}/2|}, which
+is stronger than the corresponding constraint for @code{round/},
+@math{-|@var{y}/2| <= @var{r} <= |@var{y}/2|}.  In particular,
+when @var{x} and @var{y} are integers, the number of possible remainders
+returned by @code{centered/} is @math{|@var{y}|}, whereas the number of
+possible remainders returned by @code{round/} is @math{|@var{y}|+1} when
+@var{y} is even.
+
+@lisp
+(round-quotient 123 10) @result{} 12
+(round-remainder 123 10) @result{} 3
+(round/ 123 10) @result{} 12 and 3
+(round/ 123 -10) @result{} -12 and 3
+(round/ -123 10) @result{} -12 and -3
+(round/ -123 -10) @result{} 12 and -3
+(round/ 125 10) @result{} 12 and 5
+(round/ 127 10) @result{} 13 and -3
+(round/ 135 10) @result{} 14 and -5
+(round/ -123.2 -63.5) @result{} 2.0 and 3.8
+(round/ 16/3 -10/7) @result{} -4 and -8/21
+@end lisp
+@end deffn
+
 @node Scientific
 @subsubsection Scientific Functions
 
diff --git a/libguile/numbers.c b/libguile/numbers.c
index 05840ef..90e449f 100644
--- a/libguile/numbers.c
+++ b/libguile/numbers.c
@@ -1596,6 +1596,1522 @@ scm_i_slow_exact_euclidean_divide (SCM x, SCM y)
   return scm_values (scm_list_2 (q, r));
 }
 
+static SCM scm_i_inexact_floor_quotient (double x, double y);
+static SCM scm_i_slow_exact_floor_quotient (SCM x, SCM y);
+
+SCM_PRIMITIVE_GENERIC (scm_floor_quotient, "floor-quotient", 2, 0, 0,
+		       (SCM x, SCM y),
+		       "Return the floor of @math{@var{x} / @var{y}}.\n"
+		       "@lisp\n"
+		       "(floor-quotient 123 10) @result{} 12\n"
+		       "(floor-quotient 123 -10) @result{} -13\n"
+		       "(floor-quotient -123 10) @result{} -13\n"
+		       "(floor-quotient -123 -10) @result{} 12\n"
+		       "(floor-quotient -123.2 -63.5) @result{} 1.0\n"
+		       "(floor-quotient 16/3 -10/7) @result{} -4\n"
+		       "@end lisp")
+#define FUNC_NAME s_scm_floor_quotient
+{
+  if (SCM_LIKELY (SCM_I_INUMP (x)))
+    {
+      scm_t_inum xx = SCM_I_INUM (x);
+      if (SCM_LIKELY (SCM_I_INUMP (y)))
+	{
+	  scm_t_inum yy = SCM_I_INUM (y);
+	  scm_t_inum xx1 = xx;
+	  scm_t_inum qq;
+	  if (SCM_LIKELY (yy > 0))
+	    {
+	      if (SCM_UNLIKELY (xx < 0))
+		xx1 = xx - yy + 1;
+	    }
+	  else if (SCM_UNLIKELY (yy == 0))
+	    scm_num_overflow (s_scm_floor_quotient);
+	  else if (xx > 0)
+	    xx1 = xx - yy - 1;
+	  qq = xx1 / yy;
+	  if (SCM_LIKELY (SCM_FIXABLE (qq)))
+	    return SCM_I_MAKINUM (qq);
+	  else
+	    return scm_i_inum2big (qq);
+	}
+      else if (SCM_BIGP (y))
+	{
+	  int sign = mpz_sgn (SCM_I_BIG_MPZ (y));
+	  scm_remember_upto_here_1 (y);
+	  if (sign > 0)
+	    return SCM_I_MAKINUM ((xx < 0) ? -1 : 0);
+	  else
+	    return SCM_I_MAKINUM ((xx > 0) ? -1 : 0);
+	}
+      else if (SCM_REALP (y))
+	return scm_i_inexact_floor_quotient (xx, SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_floor_quotient (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_floor_quotient, x, y, SCM_ARG2,
+			    s_scm_floor_quotient);
+    }
+  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_scm_floor_quotient);
+	  else if (SCM_UNLIKELY (yy == 1))
+	    return x;
+	  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_cdiv_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 ();
+	  mpz_fdiv_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_floor_quotient
+	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_floor_quotient (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_floor_quotient, x, y, SCM_ARG2,
+			    s_scm_floor_quotient);
+    }
+  else if (SCM_REALP (x))
+    {
+      if (SCM_REALP (y) || SCM_I_INUMP (y) ||
+	  SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_inexact_floor_quotient
+	  (SCM_REAL_VALUE (x), scm_to_double (y));
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_floor_quotient, x, y, SCM_ARG2,
+			    s_scm_floor_quotient);
+    }
+  else if (SCM_FRACTIONP (x))
+    {
+      if (SCM_REALP (y))
+	return scm_i_inexact_floor_quotient
+	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else
+	return scm_i_slow_exact_floor_quotient (x, y);
+    }
+  else
+    SCM_WTA_DISPATCH_2 (g_scm_floor_quotient, x, y, SCM_ARG1,
+			s_scm_floor_quotient);
+}
+#undef FUNC_NAME
+
+static SCM
+scm_i_inexact_floor_quotient (double x, double y)
+{
+  if (SCM_UNLIKELY (y == 0))
+    scm_num_overflow (s_scm_floor_quotient);  /* or return a NaN? */
+  else
+    return scm_from_double (floor (x / y));
+}
+
+/* Compute exact floor_quotient 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_floor_quotient (SCM x, SCM y)
+{
+  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+    SCM_WTA_DISPATCH_2 (g_scm_floor_quotient, x, y, SCM_ARG1,
+			s_scm_floor_quotient);
+  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+    SCM_WTA_DISPATCH_2 (g_scm_floor_quotient, x, y, SCM_ARG2,
+			s_scm_floor_quotient);
+  else if (scm_is_true (scm_zero_p (y)))
+    scm_num_overflow (s_scm_floor_quotient);
+  else
+    return scm_floor (scm_divide (x, y));
+}
+
+static SCM scm_i_inexact_floor_remainder (double x, double y);
+static SCM scm_i_slow_exact_floor_remainder (SCM x, SCM y);
+
+SCM_PRIMITIVE_GENERIC (scm_floor_remainder, "floor-remainder", 2, 0, 0,
+		       (SCM x, SCM y),
+		       "Return the real number @var{r} such that\n"
+		       "@math{@var{x} = @var{q}*@var{y} + @var{r}}\n"
+		       "where @math{@var{q} = floor(@var{x} / @var{y})}.\n"
+		       "@lisp\n"
+		       "(floor-remainder 123 10) @result{} 3\n"
+		       "(floor-remainder 123 -10) @result{} -7\n"
+		       "(floor-remainder -123 10) @result{} 7\n"
+		       "(floor-remainder -123 -10) @result{} -3\n"
+		       "(floor-remainder -123.2 -63.5) @result{} -59.7\n"
+		       "(floor-remainder 16/3 -10/7) @result{} -8/21\n"
+		       "@end lisp")
+#define FUNC_NAME s_scm_floor_remainder
+{
+  if (SCM_LIKELY (SCM_I_INUMP (x)))
+    {
+      scm_t_inum xx = SCM_I_INUM (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_scm_floor_remainder);
+	  else
+	    {
+	      scm_t_inum rr = xx % yy;
+	      int needs_adjustment;
+
+	      if (SCM_LIKELY (yy > 0))
+		needs_adjustment = (rr < 0);
+	      else
+		needs_adjustment = (rr > 0);
+
+	      if (needs_adjustment)
+		rr += yy;
+	      return SCM_I_MAKINUM (rr);
+	    }
+	}
+      else if (SCM_BIGP (y))
+	{
+	  int sign = mpz_sgn (SCM_I_BIG_MPZ (y));
+	  scm_remember_upto_here_1 (y);
+	  if (sign > 0)
+	    {
+	      if (xx < 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
+		return x;
+	    }
+	  else if (xx <= 0)
+	    return x;
+	  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);
+	      return scm_i_normbig (r);
+	    }
+	}
+      else if (SCM_REALP (y))
+	return scm_i_inexact_floor_remainder (xx, SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_floor_remainder (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_floor_remainder, x, y, SCM_ARG2,
+			    s_scm_floor_remainder);
+    }
+  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_scm_floor_remainder);
+	  else
+	    {
+	      scm_t_inum rr;
+	      if (yy > 0)
+		rr = mpz_fdiv_ui (SCM_I_BIG_MPZ (x), yy);
+	      else
+		rr = -mpz_cdiv_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_fdiv_r (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_floor_remainder
+	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_floor_remainder (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_floor_remainder, x, y, SCM_ARG2,
+			    s_scm_floor_remainder);
+    }
+  else if (SCM_REALP (x))
+    {
+      if (SCM_REALP (y) || SCM_I_INUMP (y) ||
+	  SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_inexact_floor_remainder
+	  (SCM_REAL_VALUE (x), scm_to_double (y));
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_floor_remainder, x, y, SCM_ARG2,
+			    s_scm_floor_remainder);
+    }
+  else if (SCM_FRACTIONP (x))
+    {
+      if (SCM_REALP (y))
+	return scm_i_inexact_floor_remainder
+	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else
+	return scm_i_slow_exact_floor_remainder (x, y);
+    }
+  else
+    SCM_WTA_DISPATCH_2 (g_scm_floor_remainder, x, y, SCM_ARG1,
+			s_scm_floor_remainder);
+}
+#undef FUNC_NAME
+
+static SCM
+scm_i_inexact_floor_remainder (double x, double y)
+{
+  /* 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_floor_quotient, such that x != q * y + r (not even
+     close).  In particular, when x is very close to a multiple of y,
+     then r might be either 0.0 or y, but those two cases must
+     correspond to different choices of q.  If r = 0.0 then q must be
+     x/y, and if r = y then q must be x/y-1.  If quotient chooses one
+     and remainder chooses the other, it would be bad.  */
+  if (SCM_UNLIKELY (y == 0))
+    scm_num_overflow (s_scm_floor_remainder);  /* or return a NaN? */
+  else
+    return scm_from_double (x - y * floor (x / y));
+}
+
+/* Compute exact floor_remainder 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_floor_remainder (SCM x, SCM y)
+{
+  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+    SCM_WTA_DISPATCH_2 (g_scm_floor_remainder, x, y, SCM_ARG1,
+			s_scm_floor_remainder);
+  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+    SCM_WTA_DISPATCH_2 (g_scm_floor_remainder, x, y, SCM_ARG2,
+			s_scm_floor_remainder);
+  else if (scm_is_true (scm_zero_p (y)))
+    scm_num_overflow (s_scm_floor_remainder);
+  else
+    return scm_difference
+      (x, scm_product (y, scm_floor (scm_divide (x, y))));
+}
+
+
+static SCM scm_i_inexact_floor_divide (double x, double y);
+static SCM scm_i_slow_exact_floor_divide (SCM x, SCM y);
+
+SCM_PRIMITIVE_GENERIC (scm_floor_divide, "floor/", 2, 0, 0,
+		       (SCM x, SCM y),
+		       "Return the integer @var{q} and the real number @var{r}\n"
+		       "such that @math{@var{x} = @var{q}*@var{y} + @var{r}}\n"
+		       "and @math{@var{q} = floor(@var{x} / @var{y})}.\n"
+		       "@lisp\n"
+		       "(floor/ 123 10) @result{} 12 and 3\n"
+		       "(floor/ 123 -10) @result{} -13 and -7\n"
+		       "(floor/ -123 10) @result{} -13 and 7\n"
+		       "(floor/ -123 -10) @result{} 12 and -3\n"
+		       "(floor/ -123.2 -63.5) @result{} 1.0 and -59.7\n"
+		       "(floor/ 16/3 -10/7) @result{} -4 and -8/21\n"
+		       "@end lisp")
+#define FUNC_NAME s_scm_floor_divide
+{
+  if (SCM_LIKELY (SCM_I_INUMP (x)))
+    {
+      scm_t_inum xx = SCM_I_INUM (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_scm_floor_divide);
+	  else
+	    {
+	      scm_t_inum qq = xx / yy;
+	      scm_t_inum rr = xx % yy;
+	      int needs_adjustment;
+	      SCM q;
+
+	      if (SCM_LIKELY (yy > 0))
+		needs_adjustment = (rr < 0);
+	      else
+		needs_adjustment = (rr > 0);
+
+	      if (needs_adjustment)
+		{
+		  rr += yy;
+		  qq--;
+		}
+
+	      if (SCM_LIKELY (SCM_FIXABLE (qq)))
+		q = SCM_I_MAKINUM (qq);
+	      else
+		q = scm_i_inum2big (qq);
+	      return scm_values (scm_list_2 (q, SCM_I_MAKINUM (rr)));
+	    }
+	}
+      else if (SCM_BIGP (y))
+	{
+	  int sign = mpz_sgn (SCM_I_BIG_MPZ (y));
+	  scm_remember_upto_here_1 (y);
+	  if (sign > 0)
+	    {
+	      if (xx < 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
+		return scm_values (scm_list_2 (SCM_INUM0, x));
+	    }
+	  else if (xx <= 0)
+	    return scm_values (scm_list_2 (SCM_INUM0, x));
+	  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);
+	      return scm_values (scm_list_2 (SCM_I_MAKINUM (-1),
+					     scm_i_normbig (r)));
+	    }
+	}
+      else if (SCM_REALP (y))
+	return scm_i_inexact_floor_divide (xx, SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_floor_divide (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_floor_divide, x, y, SCM_ARG2,
+			    s_scm_floor_divide);
+    }
+  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_scm_floor_divide);
+	  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_cdiv_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 ();
+	  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);
+	  return scm_values (scm_list_2 (scm_i_normbig (q),
+					 scm_i_normbig (r)));
+	}
+      else if (SCM_REALP (y))
+	return scm_i_inexact_floor_divide
+	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_floor_divide (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_floor_divide, x, y, SCM_ARG2,
+			    s_scm_floor_divide);
+    }
+  else if (SCM_REALP (x))
+    {
+      if (SCM_REALP (y) || SCM_I_INUMP (y) ||
+	  SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_inexact_floor_divide
+	  (SCM_REAL_VALUE (x), scm_to_double (y));
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_floor_divide, x, y, SCM_ARG2,
+			    s_scm_floor_divide);
+    }
+  else if (SCM_FRACTIONP (x))
+    {
+      if (SCM_REALP (y))
+	return scm_i_inexact_floor_divide
+	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else
+	return scm_i_slow_exact_floor_divide (x, y);
+    }
+  else
+    SCM_WTA_DISPATCH_2 (g_scm_floor_divide, x, y, SCM_ARG1,
+			s_scm_floor_divide);
+}
+#undef FUNC_NAME
+
+static SCM
+scm_i_inexact_floor_divide (double x, double y)
+{
+  double q, r;
+
+  if (SCM_UNLIKELY (y == 0))
+    scm_num_overflow (s_scm_floor_divide);  /* or return a NaN? */
+  else
+    q = floor (x / y);
+  r = x - q * y;
+  return scm_values (scm_list_2 (scm_from_double (q),
+				 scm_from_double (r)));
+}
+
+/* Compute exact floor quotient and remainder 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_floor_divide (SCM x, SCM y)
+{
+  SCM q, r;
+
+  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+    SCM_WTA_DISPATCH_2 (g_scm_floor_divide, x, y, SCM_ARG1,
+			s_scm_floor_divide);
+  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+    SCM_WTA_DISPATCH_2 (g_scm_floor_divide, x, y, SCM_ARG2,
+			s_scm_floor_divide);
+  else if (scm_is_true (scm_zero_p (y)))
+    scm_num_overflow (s_scm_floor_divide);
+  else
+    q = scm_floor (scm_divide (x, y));
+  r = scm_difference (x, scm_product (q, y));
+  return scm_values (scm_list_2 (q, r));
+}
+
+static SCM scm_i_inexact_ceiling_quotient (double x, double y);
+static SCM scm_i_slow_exact_ceiling_quotient (SCM x, SCM y);
+
+SCM_PRIMITIVE_GENERIC (scm_ceiling_quotient, "ceiling-quotient", 2, 0, 0,
+		       (SCM x, SCM y),
+		       "Return the ceiling of @math{@var{x} / @var{y}}.\n"
+		       "@lisp\n"
+		       "(ceiling-quotient 123 10) @result{} 13\n"
+		       "(ceiling-quotient 123 -10) @result{} -12\n"
+		       "(ceiling-quotient -123 10) @result{} -12\n"
+		       "(ceiling-quotient -123 -10) @result{} 13\n"
+		       "(ceiling-quotient -123.2 -63.5) @result{} 2.0\n"
+		       "(ceiling-quotient 16/3 -10/7) @result{} -3\n"
+		       "@end lisp")
+#define FUNC_NAME s_scm_ceiling_quotient
+{
+  if (SCM_LIKELY (SCM_I_INUMP (x)))
+    {
+      scm_t_inum xx = SCM_I_INUM (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_scm_ceiling_quotient);
+	  else
+	    {
+	      scm_t_inum xx1 = xx;
+	      scm_t_inum qq;
+	      if (SCM_LIKELY (yy > 0))
+		{
+		  if (SCM_LIKELY (xx >= 0))
+		    xx1 = xx + yy - 1;
+		}
+	      else if (SCM_UNLIKELY (yy == 0))
+		scm_num_overflow (s_scm_ceiling_quotient);
+	      else if (xx < 0)
+		xx1 = xx + yy + 1;
+	      qq = xx1 / yy;
+	      if (SCM_LIKELY (SCM_FIXABLE (qq)))
+		return SCM_I_MAKINUM (qq);
+	      else
+		return scm_i_inum2big (qq);
+	    }
+	}
+      else if (SCM_BIGP (y))
+	{
+	  int sign = mpz_sgn (SCM_I_BIG_MPZ (y));
+	  scm_remember_upto_here_1 (y);
+	  if (SCM_LIKELY (sign > 0))
+	    {
+	      if (SCM_LIKELY (xx > 0))
+		return SCM_INUM1;
+	      else if (SCM_UNLIKELY (xx == SCM_MOST_NEGATIVE_FIXNUM)
+		       && SCM_UNLIKELY (mpz_cmp_ui (SCM_I_BIG_MPZ (y),
+				       - SCM_MOST_NEGATIVE_FIXNUM) == 0))
+		{
+		  /* Special case: x == fixnum-min && y == abs (fixnum-min) */
+		  scm_remember_upto_here_1 (y);
+		  return SCM_I_MAKINUM (-1);
+		}
+	      else
+		return SCM_INUM0;
+	    }
+	  else if (xx >= 0)
+	    return SCM_INUM0;
+	  else
+	    return SCM_INUM1;
+	}
+      else if (SCM_REALP (y))
+	return scm_i_inexact_ceiling_quotient (xx, SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_ceiling_quotient (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_ceiling_quotient, x, y, SCM_ARG2,
+			    s_scm_ceiling_quotient);
+    }
+  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_scm_ceiling_quotient);
+	  else if (SCM_UNLIKELY (yy == 1))
+	    return x;
+	  else
+	    {
+	      SCM q = scm_i_mkbig ();
+	      if (yy > 0)
+		mpz_cdiv_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 ();
+	  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_ceiling_quotient
+	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_ceiling_quotient (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_ceiling_quotient, x, y, SCM_ARG2,
+			    s_scm_ceiling_quotient);
+    }
+  else if (SCM_REALP (x))
+    {
+      if (SCM_REALP (y) || SCM_I_INUMP (y) ||
+	  SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_inexact_ceiling_quotient
+	  (SCM_REAL_VALUE (x), scm_to_double (y));
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_ceiling_quotient, x, y, SCM_ARG2,
+			    s_scm_ceiling_quotient);
+    }
+  else if (SCM_FRACTIONP (x))
+    {
+      if (SCM_REALP (y))
+	return scm_i_inexact_ceiling_quotient
+	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else
+	return scm_i_slow_exact_ceiling_quotient (x, y);
+    }
+  else
+    SCM_WTA_DISPATCH_2 (g_scm_ceiling_quotient, x, y, SCM_ARG1,
+			s_scm_ceiling_quotient);
+}
+#undef FUNC_NAME
+
+static SCM
+scm_i_inexact_ceiling_quotient (double x, double y)
+{
+  if (SCM_UNLIKELY (y == 0))
+    scm_num_overflow (s_scm_ceiling_quotient);  /* or return a NaN? */
+  else
+    return scm_from_double (ceil (x / y));
+}
+
+/* Compute exact ceiling_quotient 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_ceiling_quotient (SCM x, SCM y)
+{
+  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+    SCM_WTA_DISPATCH_2 (g_scm_ceiling_quotient, x, y, SCM_ARG1,
+			s_scm_ceiling_quotient);
+  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+    SCM_WTA_DISPATCH_2 (g_scm_ceiling_quotient, x, y, SCM_ARG2,
+			s_scm_ceiling_quotient);
+  else if (scm_is_true (scm_zero_p (y)))
+    scm_num_overflow (s_scm_ceiling_quotient);
+  else
+    return scm_ceiling (scm_divide (x, y));
+}
+
+static SCM scm_i_inexact_ceiling_remainder (double x, double y);
+static SCM scm_i_slow_exact_ceiling_remainder (SCM x, SCM y);
+
+SCM_PRIMITIVE_GENERIC (scm_ceiling_remainder, "ceiling-remainder", 2, 0, 0,
+		       (SCM x, SCM y),
+		       "Return the real number @var{r} such that\n"
+		       "@math{@var{x} = @var{q}*@var{y} + @var{r}}\n"
+		       "where @math{@var{q} = ceiling(@var{x} / @var{y})}.\n"
+		       "@lisp\n"
+		       "(ceiling-remainder 123 10) @result{} -7\n"
+		       "(ceiling-remainder 123 -10) @result{} 3\n"
+		       "(ceiling-remainder -123 10) @result{} -3\n"
+		       "(ceiling-remainder -123 -10) @result{} 7\n"
+		       "(ceiling-remainder -123.2 -63.5) @result{} 3.8\n"
+		       "(ceiling-remainder 16/3 -10/7) @result{} 22/21\n"
+		       "@end lisp")
+#define FUNC_NAME s_scm_ceiling_remainder
+{
+  if (SCM_LIKELY (SCM_I_INUMP (x)))
+    {
+      scm_t_inum xx = SCM_I_INUM (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_scm_ceiling_remainder);
+	  else
+	    {
+	      scm_t_inum rr = xx % yy;
+	      int needs_adjustment;
+
+	      if (SCM_LIKELY (yy > 0))
+		needs_adjustment = (rr > 0);
+	      else
+		needs_adjustment = (rr < 0);
+
+	      if (needs_adjustment)
+		rr -= yy;
+	      return SCM_I_MAKINUM (rr);
+	    }
+	}
+      else if (SCM_BIGP (y))
+	{
+	  int sign = mpz_sgn (SCM_I_BIG_MPZ (y));
+	  scm_remember_upto_here_1 (y);
+	  if (SCM_LIKELY (sign > 0))
+	    {
+	      if (SCM_LIKELY (xx > 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);
+		  mpz_neg (SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (r));
+		  return scm_i_normbig (r);
+		}
+	      else if (SCM_UNLIKELY (xx == SCM_MOST_NEGATIVE_FIXNUM)
+		       && SCM_UNLIKELY (mpz_cmp_ui (SCM_I_BIG_MPZ (y),
+				       - SCM_MOST_NEGATIVE_FIXNUM) == 0))
+		{
+		  /* Special case: x == fixnum-min && y == abs (fixnum-min) */
+		  scm_remember_upto_here_1 (y);
+		  return SCM_INUM0;
+		}
+	      else
+		return x;
+	    }
+	  else if (xx >= 0)
+	    return x;
+	  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_ceiling_remainder (xx, SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_ceiling_remainder (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_ceiling_remainder, x, y, SCM_ARG2,
+			    s_scm_ceiling_remainder);
+    }
+  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_scm_ceiling_remainder);
+	  else
+	    {
+	      scm_t_inum rr;
+	      if (yy > 0)
+		rr = -mpz_cdiv_ui (SCM_I_BIG_MPZ (x), yy);
+	      else
+		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_cdiv_r (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_ceiling_remainder
+	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_ceiling_remainder (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_ceiling_remainder, x, y, SCM_ARG2,
+			    s_scm_ceiling_remainder);
+    }
+  else if (SCM_REALP (x))
+    {
+      if (SCM_REALP (y) || SCM_I_INUMP (y) ||
+	  SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_inexact_ceiling_remainder
+	  (SCM_REAL_VALUE (x), scm_to_double (y));
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_ceiling_remainder, x, y, SCM_ARG2,
+			    s_scm_ceiling_remainder);
+    }
+  else if (SCM_FRACTIONP (x))
+    {
+      if (SCM_REALP (y))
+	return scm_i_inexact_ceiling_remainder
+	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else
+	return scm_i_slow_exact_ceiling_remainder (x, y);
+    }
+  else
+    SCM_WTA_DISPATCH_2 (g_scm_ceiling_remainder, x, y, SCM_ARG1,
+			s_scm_ceiling_remainder);
+}
+#undef FUNC_NAME
+
+static SCM
+scm_i_inexact_ceiling_remainder (double x, double y)
+{
+  /* 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_ceiling_quotient, such that x != q * y + r (not even
+     close).  In particular, when x is very close to a multiple of y,
+     then r might be either 0.0 or -y, but those two cases must
+     correspond to different choices of q.  If r = 0.0 then q must be
+     x/y, and if r = -y then q must be x/y+1.  If quotient chooses one
+     and remainder chooses the other, it would be bad.  */
+  if (SCM_UNLIKELY (y == 0))
+    scm_num_overflow (s_scm_ceiling_remainder);  /* or return a NaN? */
+  else
+    return scm_from_double (x - y * ceil (x / y));
+}
+
+/* Compute exact ceiling_remainder 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_ceiling_remainder (SCM x, SCM y)
+{
+  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+    SCM_WTA_DISPATCH_2 (g_scm_ceiling_remainder, x, y, SCM_ARG1,
+			s_scm_ceiling_remainder);
+  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+    SCM_WTA_DISPATCH_2 (g_scm_ceiling_remainder, x, y, SCM_ARG2,
+			s_scm_ceiling_remainder);
+  else if (scm_is_true (scm_zero_p (y)))
+    scm_num_overflow (s_scm_ceiling_remainder);
+  else
+    return scm_difference
+      (x, scm_product (y, scm_ceiling (scm_divide (x, y))));
+}
+
+static SCM scm_i_inexact_ceiling_divide (double x, double y);
+static SCM scm_i_slow_exact_ceiling_divide (SCM x, SCM y);
+
+SCM_PRIMITIVE_GENERIC (scm_ceiling_divide, "ceiling/", 2, 0, 0,
+		       (SCM x, SCM y),
+		       "Return the integer @var{q} and the real number @var{r}\n"
+		       "such that @math{@var{x} = @var{q}*@var{y} + @var{r}}\n"
+		       "and @math{@var{q} = ceiling(@var{x} / @var{y})}.\n"
+		       "@lisp\n"
+		       "(ceiling/ 123 10) @result{} 13 and -7\n"
+		       "(ceiling/ 123 -10) @result{} -12 and 3\n"
+		       "(ceiling/ -123 10) @result{} -12 and -3\n"
+		       "(ceiling/ -123 -10) @result{} 13 and 7\n"
+		       "(ceiling/ -123.2 -63.5) @result{} 2.0 and 3.8\n"
+		       "(ceiling/ 16/3 -10/7) @result{} -3 and 22/21\n"
+		       "@end lisp")
+#define FUNC_NAME s_scm_ceiling_divide
+{
+  if (SCM_LIKELY (SCM_I_INUMP (x)))
+    {
+      scm_t_inum xx = SCM_I_INUM (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_scm_ceiling_divide);
+	  else
+	    {
+	      scm_t_inum qq = xx / yy;
+	      scm_t_inum rr = xx % yy;
+	      int needs_adjustment;
+	      SCM q;
+
+	      if (SCM_LIKELY (yy > 0))
+		needs_adjustment = (rr > 0);
+	      else
+		needs_adjustment = (rr < 0);
+
+	      if (needs_adjustment)
+		{
+		  rr -= yy;
+		  qq++;
+		}
+	      if (SCM_LIKELY (SCM_FIXABLE (qq)))
+		q = SCM_I_MAKINUM (qq);
+	      else
+		q = scm_i_inum2big (qq);
+	      return scm_values (scm_list_2 (q, SCM_I_MAKINUM (rr)));
+	    }
+	}
+      else if (SCM_BIGP (y))
+	{
+	  int sign = mpz_sgn (SCM_I_BIG_MPZ (y));
+	  scm_remember_upto_here_1 (y);
+	  if (SCM_LIKELY (sign > 0))
+	    {
+	      if (SCM_LIKELY (xx > 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);
+		  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_UNLIKELY (xx == SCM_MOST_NEGATIVE_FIXNUM)
+		       && SCM_UNLIKELY (mpz_cmp_ui (SCM_I_BIG_MPZ (y),
+				       - SCM_MOST_NEGATIVE_FIXNUM) == 0))
+		{
+		  /* 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
+		return scm_values (scm_list_2 (SCM_INUM0, x));
+	    }
+	  else if (xx >= 0)
+	    return scm_values (scm_list_2 (SCM_INUM0, x));
+	  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_ceiling_divide (xx, SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_ceiling_divide (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_ceiling_divide, x, y, SCM_ARG2,
+			    s_scm_ceiling_divide);
+    }
+  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_scm_ceiling_divide);
+	  else
+	    {
+	      SCM q = scm_i_mkbig ();
+	      SCM r = scm_i_mkbig ();
+	      if (yy > 0)
+		mpz_cdiv_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 ();
+	  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_ceiling_divide
+	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_ceiling_divide (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_ceiling_divide, x, y, SCM_ARG2,
+			    s_scm_ceiling_divide);
+    }
+  else if (SCM_REALP (x))
+    {
+      if (SCM_REALP (y) || SCM_I_INUMP (y) ||
+	  SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_inexact_ceiling_divide
+	  (SCM_REAL_VALUE (x), scm_to_double (y));
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_ceiling_divide, x, y, SCM_ARG2,
+			    s_scm_ceiling_divide);
+    }
+  else if (SCM_FRACTIONP (x))
+    {
+      if (SCM_REALP (y))
+	return scm_i_inexact_ceiling_divide
+	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else
+	return scm_i_slow_exact_ceiling_divide (x, y);
+    }
+  else
+    SCM_WTA_DISPATCH_2 (g_scm_ceiling_divide, x, y, SCM_ARG1,
+			s_scm_ceiling_divide);
+}
+#undef FUNC_NAME
+
+static SCM
+scm_i_inexact_ceiling_divide (double x, double y)
+{
+  double q, r;
+
+  if (SCM_UNLIKELY (y == 0))
+    scm_num_overflow (s_scm_ceiling_divide);  /* or return a NaN? */
+  else
+    q = ceil (x / y);
+  r = x - q * y;
+  return scm_values (scm_list_2 (scm_from_double (q),
+				 scm_from_double (r)));
+}
+
+/* Compute exact ceiling quotient and remainder 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_ceiling_divide (SCM x, SCM y)
+{
+  SCM q, r;
+
+  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+    SCM_WTA_DISPATCH_2 (g_scm_ceiling_divide, x, y, SCM_ARG1,
+			s_scm_ceiling_divide);
+  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+    SCM_WTA_DISPATCH_2 (g_scm_ceiling_divide, x, y, SCM_ARG2,
+			s_scm_ceiling_divide);
+  else if (scm_is_true (scm_zero_p (y)))
+    scm_num_overflow (s_scm_ceiling_divide);
+  else
+    q = scm_ceiling (scm_divide (x, y));
+  r = scm_difference (x, scm_product (q, y));
+  return scm_values (scm_list_2 (q, r));
+}
+
+static SCM scm_i_inexact_truncate_quotient (double x, double y);
+static SCM scm_i_slow_exact_truncate_quotient (SCM x, SCM y);
+
+SCM_PRIMITIVE_GENERIC (scm_truncate_quotient, "truncate-quotient", 2, 0, 0,
+		       (SCM x, SCM y),
+		       "Return @math{@var{x} / @var{y}} rounded toward zero.\n"
+		       "@lisp\n"
+		       "(truncate-quotient 123 10) @result{} 12\n"
+		       "(truncate-quotient 123 -10) @result{} -12\n"
+		       "(truncate-quotient -123 10) @result{} -12\n"
+		       "(truncate-quotient -123 -10) @result{} 12\n"
+		       "(truncate-quotient -123.2 -63.5) @result{} 1.0\n"
+		       "(truncate-quotient 16/3 -10/7) @result{} -3\n"
+		       "@end lisp")
+#define FUNC_NAME s_scm_truncate_quotient
+{
+  if (SCM_LIKELY (SCM_I_INUMP (x)))
+    {
+      scm_t_inum xx = SCM_I_INUM (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_scm_truncate_quotient);
+	  else
+	    {
+	      scm_t_inum qq = xx / yy;
+	      if (SCM_LIKELY (SCM_FIXABLE (qq)))
+		return SCM_I_MAKINUM (qq);
+	      else
+		return scm_i_inum2big (qq);
+	    }
+	}
+      else if (SCM_BIGP (y))
+	{
+	  if (SCM_UNLIKELY (xx == SCM_MOST_NEGATIVE_FIXNUM)
+	      && SCM_UNLIKELY (mpz_cmp_ui (SCM_I_BIG_MPZ (y),
+					   - SCM_MOST_NEGATIVE_FIXNUM) == 0))
+	    {
+	      /* Special case: x == fixnum-min && y == abs (fixnum-min) */
+	      scm_remember_upto_here_1 (y);
+	      return SCM_I_MAKINUM (-1);
+	    }
+	  else
+	    return SCM_INUM0;
+	}
+      else if (SCM_REALP (y))
+	return scm_i_inexact_truncate_quotient (xx, SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_truncate_quotient (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_truncate_quotient, x, y, SCM_ARG2,
+			    s_scm_truncate_quotient);
+    }
+  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_scm_truncate_quotient);
+	  else if (SCM_UNLIKELY (yy == 1))
+	    return x;
+	  else
+	    {
+	      SCM q = scm_i_mkbig ();
+	      if (yy > 0)
+		mpz_tdiv_q_ui (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (x), yy);
+	      else
+		{
+		  mpz_tdiv_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 ();
+	  mpz_tdiv_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_truncate_quotient
+	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_truncate_quotient (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_truncate_quotient, x, y, SCM_ARG2,
+			    s_scm_truncate_quotient);
+    }
+  else if (SCM_REALP (x))
+    {
+      if (SCM_REALP (y) || SCM_I_INUMP (y) ||
+	  SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_inexact_truncate_quotient
+	  (SCM_REAL_VALUE (x), scm_to_double (y));
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_truncate_quotient, x, y, SCM_ARG2,
+			    s_scm_truncate_quotient);
+    }
+  else if (SCM_FRACTIONP (x))
+    {
+      if (SCM_REALP (y))
+	return scm_i_inexact_truncate_quotient
+	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else
+	return scm_i_slow_exact_truncate_quotient (x, y);
+    }
+  else
+    SCM_WTA_DISPATCH_2 (g_scm_truncate_quotient, x, y, SCM_ARG1,
+			s_scm_truncate_quotient);
+}
+#undef FUNC_NAME
+
+static SCM
+scm_i_inexact_truncate_quotient (double x, double y)
+{
+  if (SCM_UNLIKELY (y == 0))
+    scm_num_overflow (s_scm_truncate_quotient);  /* or return a NaN? */
+  else
+    return scm_from_double (trunc (x / y));
+}
+
+/* Compute exact truncate_quotient 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_truncate_quotient (SCM x, SCM y)
+{
+  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+    SCM_WTA_DISPATCH_2 (g_scm_truncate_quotient, x, y, SCM_ARG1,
+			s_scm_truncate_quotient);
+  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+    SCM_WTA_DISPATCH_2 (g_scm_truncate_quotient, x, y, SCM_ARG2,
+			s_scm_truncate_quotient);
+  else if (scm_is_true (scm_zero_p (y)))
+    scm_num_overflow (s_scm_truncate_quotient);
+  else
+    return scm_truncate_number (scm_divide (x, y));
+}
+
+static SCM scm_i_inexact_truncate_remainder (double x, double y);
+static SCM scm_i_slow_exact_truncate_remainder (SCM x, SCM y);
+
+SCM_PRIMITIVE_GENERIC (scm_truncate_remainder, "truncate-remainder", 2, 0, 0,
+		       (SCM x, SCM y),
+		       "Return the real number @var{r} such that\n"
+		       "@math{@var{x} = @var{q}*@var{y} + @var{r}}\n"
+		       "where @math{@var{q} = truncate(@var{x} / @var{y})}.\n"
+		       "@lisp\n"
+		       "(truncate-remainder 123 10) @result{} 3\n"
+		       "(truncate-remainder 123 -10) @result{} 3\n"
+		       "(truncate-remainder -123 10) @result{} -3\n"
+		       "(truncate-remainder -123 -10) @result{} -3\n"
+		       "(truncate-remainder -123.2 -63.5) @result{} -59.7\n"
+		       "(truncate-remainder 16/3 -10/7) @result{} 22/21\n"
+		       "@end lisp")
+#define FUNC_NAME s_scm_truncate_remainder
+{
+  if (SCM_LIKELY (SCM_I_INUMP (x)))
+    {
+      scm_t_inum xx = SCM_I_INUM (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_scm_truncate_remainder);
+	  else
+	    return SCM_I_MAKINUM (xx % yy);
+	}
+      else if (SCM_BIGP (y))
+	{
+	  if (SCM_UNLIKELY (xx == SCM_MOST_NEGATIVE_FIXNUM)
+	      && SCM_UNLIKELY (mpz_cmp_ui (SCM_I_BIG_MPZ (y),
+					   - SCM_MOST_NEGATIVE_FIXNUM) == 0))
+	    {
+	      /* Special case: x == fixnum-min && y == abs (fixnum-min) */
+	      scm_remember_upto_here_1 (y);
+	      return SCM_INUM0;
+	    }
+	  else
+	    return x;
+	}
+      else if (SCM_REALP (y))
+	return scm_i_inexact_truncate_remainder (xx, SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_truncate_remainder (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_truncate_remainder, x, y, SCM_ARG2,
+			    s_scm_truncate_remainder);
+    }
+  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_scm_truncate_remainder);
+	  else
+	    {
+	      scm_t_inum rr = (mpz_tdiv_ui (SCM_I_BIG_MPZ (x),
+					    (yy > 0) ? yy : -yy)
+			       * mpz_sgn (SCM_I_BIG_MPZ (x)));
+	      scm_remember_upto_here_1 (x);
+	      return SCM_I_MAKINUM (rr);
+	    }
+	}
+      else if (SCM_BIGP (y))
+	{
+	  SCM r = scm_i_mkbig ();
+	  mpz_tdiv_r (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_truncate_remainder
+	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_truncate_remainder (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_truncate_remainder, x, y, SCM_ARG2,
+			    s_scm_truncate_remainder);
+    }
+  else if (SCM_REALP (x))
+    {
+      if (SCM_REALP (y) || SCM_I_INUMP (y) ||
+	  SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_inexact_truncate_remainder
+	  (SCM_REAL_VALUE (x), scm_to_double (y));
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_truncate_remainder, x, y, SCM_ARG2,
+			    s_scm_truncate_remainder);
+    }
+  else if (SCM_FRACTIONP (x))
+    {
+      if (SCM_REALP (y))
+	return scm_i_inexact_truncate_remainder
+	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else
+	return scm_i_slow_exact_truncate_remainder (x, y);
+    }
+  else
+    SCM_WTA_DISPATCH_2 (g_scm_truncate_remainder, x, y, SCM_ARG1,
+			s_scm_truncate_remainder);
+}
+#undef FUNC_NAME
+
+static SCM
+scm_i_inexact_truncate_remainder (double x, double y)
+{
+  /* 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_truncate_quotient, such that x != q * y + r (not even
+     close).  In particular, when x is very close to a multiple of y,
+     then r might be either 0.0 or sgn(x)*|y|, but those two cases must
+     correspond to different choices of q.  If quotient chooses one and
+     remainder chooses the other, it would be bad.  */
+  if (SCM_UNLIKELY (y == 0))
+    scm_num_overflow (s_scm_truncate_remainder);  /* or return a NaN? */
+  else
+    return scm_from_double (x - y * trunc (x / y));
+}
+
+/* Compute exact truncate_remainder 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_truncate_remainder (SCM x, SCM y)
+{
+  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+    SCM_WTA_DISPATCH_2 (g_scm_truncate_remainder, x, y, SCM_ARG1,
+			s_scm_truncate_remainder);
+  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+    SCM_WTA_DISPATCH_2 (g_scm_truncate_remainder, x, y, SCM_ARG2,
+			s_scm_truncate_remainder);
+  else if (scm_is_true (scm_zero_p (y)))
+    scm_num_overflow (s_scm_truncate_remainder);
+  else
+    return scm_difference
+      (x, scm_product (y, scm_truncate_number (scm_divide (x, y))));
+}
+
+
+static SCM scm_i_inexact_truncate_divide (double x, double y);
+static SCM scm_i_slow_exact_truncate_divide (SCM x, SCM y);
+
+SCM_PRIMITIVE_GENERIC (scm_truncate_divide, "truncate/", 2, 0, 0,
+		       (SCM x, SCM y),
+		       "Return the integer @var{q} and the real number @var{r}\n"
+		       "such that @math{@var{x} = @var{q}*@var{y} + @var{r}}\n"
+		       "and @math{@var{q} = truncate(@var{x} / @var{y})}.\n"
+		       "@lisp\n"
+		       "(truncate/ 123 10) @result{} 12 and 3\n"
+		       "(truncate/ 123 -10) @result{} -12 and 3\n"
+		       "(truncate/ -123 10) @result{} -12 and -3\n"
+		       "(truncate/ -123 -10) @result{} 12 and -3\n"
+		       "(truncate/ -123.2 -63.5) @result{} 1.0 and -59.7\n"
+		       "(truncate/ 16/3 -10/7) @result{} -3 and 22/21\n"
+		       "@end lisp")
+#define FUNC_NAME s_scm_truncate_divide
+{
+  if (SCM_LIKELY (SCM_I_INUMP (x)))
+    {
+      scm_t_inum xx = SCM_I_INUM (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_scm_truncate_divide);
+	  else
+	    {
+	      scm_t_inum qq = xx / yy;
+	      scm_t_inum rr = xx % yy;
+	      SCM q;
+
+	      if (SCM_LIKELY (SCM_FIXABLE (qq)))
+		q = SCM_I_MAKINUM (qq);
+	      else
+		q = scm_i_inum2big (qq);
+	      return scm_values (scm_list_2 (q, SCM_I_MAKINUM (rr)));
+	    }
+	}
+      else if (SCM_BIGP (y))
+	{
+	  if (SCM_UNLIKELY (xx == SCM_MOST_NEGATIVE_FIXNUM)
+	      && SCM_UNLIKELY (mpz_cmp_ui (SCM_I_BIG_MPZ (y),
+					   - SCM_MOST_NEGATIVE_FIXNUM) == 0))
+	    {
+	      /* 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
+	    return scm_values (scm_list_2 (SCM_INUM0, x));
+	}
+      else if (SCM_REALP (y))
+	return scm_i_inexact_truncate_divide (xx, SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_truncate_divide (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_truncate_divide, x, y, SCM_ARG2,
+			    s_scm_truncate_divide);
+    }
+  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_scm_truncate_divide);
+	  else
+	    {
+	      SCM q = scm_i_mkbig ();
+	      scm_t_inum rr;
+	      if (yy > 0)
+		rr = mpz_tdiv_q_ui (SCM_I_BIG_MPZ (q),
+				    SCM_I_BIG_MPZ (x), yy);
+	      else
+		{
+		  rr = mpz_tdiv_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));
+		}
+	      rr *= mpz_sgn (SCM_I_BIG_MPZ (x));
+	      scm_remember_upto_here_1 (x);
+	      return scm_values (scm_list_2 (scm_i_normbig (q),
+					     SCM_I_MAKINUM (rr)));
+	    }
+	}
+      else if (SCM_BIGP (y))
+	{
+	  SCM q = scm_i_mkbig ();
+	  SCM r = scm_i_mkbig ();
+	  mpz_tdiv_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_truncate_divide
+	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_truncate_divide (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_truncate_divide, x, y, SCM_ARG2,
+			    s_scm_truncate_divide);
+    }
+  else if (SCM_REALP (x))
+    {
+      if (SCM_REALP (y) || SCM_I_INUMP (y) ||
+	  SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_inexact_truncate_divide
+	  (SCM_REAL_VALUE (x), scm_to_double (y));
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_truncate_divide, x, y, SCM_ARG2,
+			    s_scm_truncate_divide);
+    }
+  else if (SCM_FRACTIONP (x))
+    {
+      if (SCM_REALP (y))
+	return scm_i_inexact_truncate_divide
+	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else
+	return scm_i_slow_exact_truncate_divide (x, y);
+    }
+  else
+    SCM_WTA_DISPATCH_2 (g_scm_truncate_divide, x, y, SCM_ARG1,
+			s_scm_truncate_divide);
+}
+#undef FUNC_NAME
+
+static SCM
+scm_i_inexact_truncate_divide (double x, double y)
+{
+  double q, r;
+
+  if (SCM_UNLIKELY (y == 0))
+    scm_num_overflow (s_scm_truncate_divide);  /* or return a NaN? */
+  else
+    q = trunc (x / y);
+  r = x - q * y;
+  return scm_values (scm_list_2 (scm_from_double (q),
+				 scm_from_double (r)));
+}
+
+/* Compute exact truncate quotient and remainder 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_truncate_divide (SCM x, SCM y)
+{
+  SCM q, r;
+
+  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+    SCM_WTA_DISPATCH_2 (g_scm_truncate_divide, x, y, SCM_ARG1,
+			s_scm_truncate_divide);
+  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+    SCM_WTA_DISPATCH_2 (g_scm_truncate_divide, x, y, SCM_ARG2,
+			s_scm_truncate_divide);
+  else if (scm_is_true (scm_zero_p (y)))
+    scm_num_overflow (s_scm_truncate_divide);
+  else
+    q = scm_truncate_number (scm_divide (x, y));
+  r = scm_difference (x, scm_product (q, y));
+  return scm_values (scm_list_2 (q, r));
+}
+
 static SCM scm_i_inexact_centered_quotient (double x, double y);
 static SCM scm_i_bigint_centered_quotient (SCM x, SCM y);
 static SCM scm_i_slow_exact_centered_quotient (SCM x, SCM y);
@@ -2306,6 +3822,670 @@ scm_i_slow_exact_centered_divide (SCM x, SCM y)
   return scm_values (scm_list_2 (q, r));
 }
 
+static SCM scm_i_inexact_round_quotient (double x, double y);
+static SCM scm_i_bigint_round_quotient (SCM x, SCM y);
+static SCM scm_i_slow_exact_round_quotient (SCM x, SCM y);
+
+SCM_PRIMITIVE_GENERIC (scm_round_quotient, "round-quotient", 2, 0, 0,
+		       (SCM x, SCM y),
+		       "Return @math{@var{x} / @var{y}} to the nearest integer,\n"
+		       "with ties going to the nearest even integer.\n"
+		       "@lisp\n"
+		       "(round-quotient 123 10) @result{} 12\n"
+		       "(round-quotient 123 -10) @result{} -12\n"
+		       "(round-quotient -123 10) @result{} -12\n"
+		       "(round-quotient -123 -10) @result{} 12\n"
+		       "(round-quotient 125 10) @result{} 12\n"
+		       "(round-quotient 127 10) @result{} 13\n"
+		       "(round-quotient 135 10) @result{} 14\n"
+		       "(round-quotient -123.2 -63.5) @result{} 2.0\n"
+		       "(round-quotient 16/3 -10/7) @result{} -4\n"
+		       "@end lisp")
+#define FUNC_NAME s_scm_round_quotient
+{
+  if (SCM_LIKELY (SCM_I_INUMP (x)))
+    {
+      scm_t_inum xx = SCM_I_INUM (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_scm_round_quotient);
+	  else
+	    {
+	      scm_t_inum qq = xx / yy;
+	      scm_t_inum rr = xx % yy;
+	      scm_t_inum ay = yy;
+	      scm_t_inum r2 = 2 * rr;
+
+	      if (SCM_LIKELY (yy < 0))
+		{
+		  ay = -ay;
+		  r2 = -r2;
+		}
+
+	      if (qq & 1L)
+		{
+		  if (r2 >= ay)
+		    qq++;
+		  else if (r2 <= -ay)
+		    qq--;
+		}
+	      else
+		{
+		  if (r2 > ay)
+		    qq++;
+		  else if (r2 < -ay)
+		    qq--;
+		}
+	      if (SCM_LIKELY (SCM_FIXABLE (qq)))
+		return SCM_I_MAKINUM (qq);
+	      else
+		return scm_i_inum2big (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_round_quotient */
+	  return scm_i_bigint_round_quotient (scm_i_long2big (xx), y);
+	}
+      else if (SCM_REALP (y))
+	return scm_i_inexact_round_quotient (xx, SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_round_quotient (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_round_quotient, x, y, SCM_ARG2,
+			    s_scm_round_quotient);
+    }
+  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_scm_round_quotient);
+	  else if (SCM_UNLIKELY (yy == 1))
+	    return x;
+	  else
+	    {
+	      SCM q = scm_i_mkbig ();
+	      scm_t_inum rr;
+	      int needs_adjustment;
+
+	      if (yy > 0)
+		{
+		  rr = mpz_fdiv_q_ui (SCM_I_BIG_MPZ (q),
+				      SCM_I_BIG_MPZ (x), yy);
+		  if (mpz_odd_p (SCM_I_BIG_MPZ (q)))
+		    needs_adjustment = (2*rr >= yy);
+		  else
+		    needs_adjustment = (2*rr > yy);
+		}
+	      else
+		{
+		  rr = - mpz_cdiv_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));
+		  if (mpz_odd_p (SCM_I_BIG_MPZ (q)))
+		    needs_adjustment = (2*rr <= yy);
+		  else
+		    needs_adjustment = (2*rr < yy);
+		}
+	      scm_remember_upto_here_1 (x);
+	      if (needs_adjustment)
+		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_round_quotient (x, y);
+      else if (SCM_REALP (y))
+	return scm_i_inexact_round_quotient
+	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_round_quotient (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_round_quotient, x, y, SCM_ARG2,
+			    s_scm_round_quotient);
+    }
+  else if (SCM_REALP (x))
+    {
+      if (SCM_REALP (y) || SCM_I_INUMP (y) ||
+	  SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_inexact_round_quotient
+	  (SCM_REAL_VALUE (x), scm_to_double (y));
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_round_quotient, x, y, SCM_ARG2,
+			    s_scm_round_quotient);
+    }
+  else if (SCM_FRACTIONP (x))
+    {
+      if (SCM_REALP (y))
+	return scm_i_inexact_round_quotient
+	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else
+	return scm_i_slow_exact_round_quotient (x, y);
+    }
+  else
+    SCM_WTA_DISPATCH_2 (g_scm_round_quotient, x, y, SCM_ARG1,
+			s_scm_round_quotient);
+}
+#undef FUNC_NAME
+
+static SCM
+scm_i_inexact_round_quotient (double x, double y)
+{
+  if (SCM_UNLIKELY (y == 0))
+    scm_num_overflow (s_scm_round_quotient);  /* or return a NaN? */
+  else
+    return scm_from_double (scm_c_round (x / y));
+}
+
+/* Assumes that both x and y are bigints, though
+   x might be able to fit into a fixnum. */
+static SCM
+scm_i_bigint_round_quotient (SCM x, SCM y)
+{
+  SCM q, r, r2;
+  int cmp, needs_adjustment;
+
+  /* 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 ();
+  r2 = scm_i_mkbig ();
+
+  mpz_fdiv_qr (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (r),
+	       SCM_I_BIG_MPZ (x), SCM_I_BIG_MPZ (y));
+  mpz_mul_2exp (SCM_I_BIG_MPZ (r2), SCM_I_BIG_MPZ (r), 1);  /* r2 = 2*r */
+  scm_remember_upto_here_2 (x, r);
+
+  cmp = mpz_cmpabs (SCM_I_BIG_MPZ (r2), SCM_I_BIG_MPZ (y));
+  if (mpz_odd_p (SCM_I_BIG_MPZ (q)))
+    needs_adjustment = (cmp >= 0);
+  else
+    needs_adjustment = (cmp > 0);
+  scm_remember_upto_here_2 (r2, y);
+
+  if (needs_adjustment)
+    mpz_add_ui (SCM_I_BIG_MPZ (q), SCM_I_BIG_MPZ (q), 1);
+
+  return scm_i_normbig (q);
+}
+
+/* Compute exact round quotient 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_round_quotient (SCM x, SCM y)
+{
+  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+    SCM_WTA_DISPATCH_2 (g_scm_round_quotient, x, y, SCM_ARG1,
+			s_scm_round_quotient);
+  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+    SCM_WTA_DISPATCH_2 (g_scm_round_quotient, x, y, SCM_ARG2,
+			s_scm_round_quotient);
+  else if (scm_is_true (scm_zero_p (y)))
+    scm_num_overflow (s_scm_round_quotient);
+  else
+    return scm_round_number (scm_divide (x, y));
+}
+
+static SCM scm_i_inexact_round_remainder (double x, double y);
+static SCM scm_i_bigint_round_remainder (SCM x, SCM y);
+static SCM scm_i_slow_exact_round_remainder (SCM x, SCM y);
+
+SCM_PRIMITIVE_GENERIC (scm_round_remainder, "round-remainder", 2, 0, 0,
+		       (SCM x, SCM y),
+		       "Return the real number @var{r} such that\n"
+		       "@math{@var{x} = @var{q}*@var{y} + @var{r}}, where\n"
+		       "@var{q} is @math{@var{x} / @var{y}} rounded to the\n"
+		       "nearest integer, with ties going to the nearest\n"
+		       "even integer.\n"
+		       "@lisp\n"
+		       "(round-remainder 123 10) @result{} 3\n"
+		       "(round-remainder 123 -10) @result{} 3\n"
+		       "(round-remainder -123 10) @result{} -3\n"
+		       "(round-remainder -123 -10) @result{} -3\n"
+		       "(round-remainder 125 10) @result{} 5\n"
+		       "(round-remainder 127 10) @result{} -3\n"
+		       "(round-remainder 135 10) @result{} -5\n"
+		       "(round-remainder -123.2 -63.5) @result{} 3.8\n"
+		       "(round-remainder 16/3 -10/7) @result{} -8/21\n"
+		       "@end lisp")
+#define FUNC_NAME s_scm_round_remainder
+{
+  if (SCM_LIKELY (SCM_I_INUMP (x)))
+    {
+      scm_t_inum xx = SCM_I_INUM (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_scm_round_remainder);
+	  else
+	    {
+	      scm_t_inum qq = xx / yy;
+	      scm_t_inum rr = xx % yy;
+	      scm_t_inum ay = yy;
+	      scm_t_inum r2 = 2 * rr;
+
+	      if (SCM_LIKELY (yy < 0))
+		{
+		  ay = -ay;
+		  r2 = -r2;
+		}
+
+	      if (qq & 1L)
+		{
+		  if (r2 >= ay)
+		    rr -= yy;
+		  else if (r2 <= -ay)
+		    rr += yy;
+		}
+	      else
+		{
+		  if (r2 > ay)
+		    rr -= yy;
+		  else if (r2 < -ay)
+		    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_round_remainder */
+	  return scm_i_bigint_round_remainder
+	    (scm_i_long2big (xx), y);
+	}
+      else if (SCM_REALP (y))
+	return scm_i_inexact_round_remainder (xx, SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_round_remainder (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_round_remainder, x, y, SCM_ARG2,
+			    s_scm_round_remainder);
+    }
+  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_scm_round_remainder);
+	  else
+	    {
+	      SCM q = scm_i_mkbig ();
+	      scm_t_inum rr;
+	      int needs_adjustment;
+
+	      if (yy > 0)
+		{
+		  rr = mpz_fdiv_q_ui (SCM_I_BIG_MPZ (q),
+				      SCM_I_BIG_MPZ (x), yy);
+		  if (mpz_odd_p (SCM_I_BIG_MPZ (q)))
+		    needs_adjustment = (2*rr >= yy);
+		  else
+		    needs_adjustment = (2*rr > yy);
+		}
+	      else
+		{
+		  rr = - mpz_cdiv_q_ui (SCM_I_BIG_MPZ (q),
+					SCM_I_BIG_MPZ (x), -yy);
+		  if (mpz_odd_p (SCM_I_BIG_MPZ (q)))
+		    needs_adjustment = (2*rr <= yy);
+		  else
+		    needs_adjustment = (2*rr < yy);
+		}
+	      scm_remember_upto_here_2 (x, q);
+	      if (needs_adjustment)
+		rr -= yy;
+	      return SCM_I_MAKINUM (rr);
+	    }
+	}
+      else if (SCM_BIGP (y))
+	return scm_i_bigint_round_remainder (x, y);
+      else if (SCM_REALP (y))
+	return scm_i_inexact_round_remainder
+	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_round_remainder (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_round_remainder, x, y, SCM_ARG2,
+			    s_scm_round_remainder);
+    }
+  else if (SCM_REALP (x))
+    {
+      if (SCM_REALP (y) || SCM_I_INUMP (y) ||
+	  SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_inexact_round_remainder
+	  (SCM_REAL_VALUE (x), scm_to_double (y));
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_round_remainder, x, y, SCM_ARG2,
+			    s_scm_round_remainder);
+    }
+  else if (SCM_FRACTIONP (x))
+    {
+      if (SCM_REALP (y))
+	return scm_i_inexact_round_remainder
+	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else
+	return scm_i_slow_exact_round_remainder (x, y);
+    }
+  else
+    SCM_WTA_DISPATCH_2 (g_scm_round_remainder, x, y, SCM_ARG1,
+			s_scm_round_remainder);
+}
+#undef FUNC_NAME
+
+static SCM
+scm_i_inexact_round_remainder (double x, double y)
+{
+  /* 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_round_quotient, 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), but those two
+     cases must correspond to different choices of q.  If quotient
+     chooses one and remainder chooses the other, it would be bad. */
+
+  if (SCM_UNLIKELY (y == 0))
+    scm_num_overflow (s_scm_round_remainder);  /* or return a NaN? */
+  else
+    {
+      double q = scm_c_round (x / y);
+      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_round_remainder (SCM x, SCM y)
+{
+  SCM q, r, r2;
+  int cmp, needs_adjustment;
+
+  /* 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 ();
+  r2 = scm_i_mkbig ();
+
+  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_1 (x);
+  mpz_mul_2exp (SCM_I_BIG_MPZ (r2), SCM_I_BIG_MPZ (r), 1);  /* r2 = 2*r */
+
+  cmp = mpz_cmpabs (SCM_I_BIG_MPZ (r2), SCM_I_BIG_MPZ (y));
+  if (mpz_odd_p (SCM_I_BIG_MPZ (q)))
+    needs_adjustment = (cmp >= 0);
+  else
+    needs_adjustment = (cmp > 0);
+  scm_remember_upto_here_2 (q, r2);
+
+  if (needs_adjustment)
+    mpz_sub (SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (r), SCM_I_BIG_MPZ (y));
+
+  scm_remember_upto_here_1 (y);
+  return scm_i_normbig (r);
+}
+
+/* Compute exact round_remainder 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_round_remainder (SCM x, SCM y)
+{
+  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+    SCM_WTA_DISPATCH_2 (g_scm_round_remainder, x, y, SCM_ARG1,
+			s_scm_round_remainder);
+  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+    SCM_WTA_DISPATCH_2 (g_scm_round_remainder, x, y, SCM_ARG2,
+			s_scm_round_remainder);
+  else if (scm_is_true (scm_zero_p (y)))
+    scm_num_overflow (s_scm_round_remainder);
+  else
+    {
+      SCM q = scm_round_number (scm_divide (x, y));
+      return scm_difference (x, scm_product (q, y));
+    }
+}
+
+
+static SCM scm_i_inexact_round_divide (double x, double y);
+static SCM scm_i_bigint_round_divide (SCM x, SCM y);
+static SCM scm_i_slow_exact_round_divide (SCM x, SCM y);
+
+SCM_PRIMITIVE_GENERIC (scm_round_divide, "round/", 2, 0, 0,
+		       (SCM x, SCM y),
+		       "Return the integer @var{q} and the real number @var{r}\n"
+		       "such that @math{@var{x} = @var{q}*@var{y} + @var{r}}\n"
+		       "and @var{q} is @math{@var{x} / @var{y}} rounded to the\n"
+		       "nearest integer, with ties going to the nearest even integer.\n"
+		       "@lisp\n"
+		       "(round/ 123 10) @result{} 12 and 3\n"
+		       "(round/ 123 -10) @result{} -12 and 3\n"
+		       "(round/ -123 10) @result{} -12 and -3\n"
+		       "(round/ -123 -10) @result{} 12 and -3\n"
+		       "(round/ 125 10) @result{} 12 and 5\n"
+		       "(round/ 127 10) @result{} 13 and -3\n"
+		       "(round/ 135 10) @result{} 14 and -5\n"
+		       "(round/ -123.2 -63.5) @result{} 2.0 and 3.8\n"
+		       "(round/ 16/3 -10/7) @result{} -4 and -8/21\n"
+		       "@end lisp")
+#define FUNC_NAME s_scm_round_divide
+{
+  if (SCM_LIKELY (SCM_I_INUMP (x)))
+    {
+      scm_t_inum xx = SCM_I_INUM (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_scm_round_divide);
+	  else
+	    {
+	      scm_t_inum qq = xx / yy;
+	      scm_t_inum rr = xx % yy;
+	      scm_t_inum ay = yy;
+	      scm_t_inum r2 = 2 * rr;
+	      SCM q;
+
+	      if (SCM_LIKELY (yy < 0))
+		{
+		  ay = -ay;
+		  r2 = -r2;
+		}
+
+	      if (qq & 1L)
+		{
+		  if (r2 >= ay)
+		    { qq++; rr -= yy; }
+		  else if (r2 <= -ay)
+		    { qq--; rr += yy; }
+		}
+	      else
+		{
+		  if (r2 > ay)
+		    { qq++; rr -= yy; }
+		  else if (r2 < -ay)
+		    { qq--; rr += yy; }
+		}
+	      if (SCM_LIKELY (SCM_FIXABLE (qq)))
+		q = SCM_I_MAKINUM (qq);
+	      else
+		q = scm_i_inum2big (qq);
+	      return scm_values (scm_list_2 (q, 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_round_divide */
+	  return scm_i_bigint_round_divide
+	    (scm_i_long2big (SCM_I_INUM (x)), y);
+	}
+      else if (SCM_REALP (y))
+	return scm_i_inexact_round_divide (xx, SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_round_divide (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_round_divide, x, y, SCM_ARG2,
+			    s_scm_round_divide);
+    }
+  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_scm_round_divide);
+	  else
+	    {
+	      SCM q = scm_i_mkbig ();
+	      scm_t_inum rr;
+	      int needs_adjustment;
+
+	      if (yy > 0)
+		{
+		  rr = mpz_fdiv_q_ui (SCM_I_BIG_MPZ (q),
+				      SCM_I_BIG_MPZ (x), yy);
+		  if (mpz_odd_p (SCM_I_BIG_MPZ (q)))
+		    needs_adjustment = (2*rr >= yy);
+		  else
+		    needs_adjustment = (2*rr > yy);
+		}
+	      else
+		{
+		  rr = - mpz_cdiv_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));
+		  if (mpz_odd_p (SCM_I_BIG_MPZ (q)))
+		    needs_adjustment = (2*rr <= yy);
+		  else
+		    needs_adjustment = (2*rr < yy);
+		}
+	      scm_remember_upto_here_1 (x);
+	      if (needs_adjustment)
+		{
+		  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_round_divide (x, y);
+      else if (SCM_REALP (y))
+	return scm_i_inexact_round_divide
+	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
+      else if (SCM_FRACTIONP (y))
+	return scm_i_slow_exact_round_divide (x, y);
+      else
+	SCM_WTA_DISPATCH_2 (g_scm_round_divide, x, y, SCM_ARG2,
+			    s_scm_round_divide);
+    }
+  else if (SCM_REALP (x))
+    {
+      if (SCM_REALP (y) || SCM_I_INUMP (y) ||
+	  SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_inexact_round_divide
+	  (SCM_REAL_VALUE (x), scm_to_double (y));
+     else
+	SCM_WTA_DISPATCH_2 (g_scm_round_divide, x, y, SCM_ARG2,
+			    s_scm_round_divide);
+    }
+  else if (SCM_FRACTIONP (x))
+    {
+      if (SCM_REALP (y))
+	return scm_i_inexact_round_divide
+	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else
+	return scm_i_slow_exact_round_divide (x, y);
+    }
+  else
+    SCM_WTA_DISPATCH_2 (g_scm_round_divide, x, y, SCM_ARG1,
+			s_scm_round_divide);
+}
+#undef FUNC_NAME
+
+static SCM
+scm_i_inexact_round_divide (double x, double y)
+{
+  if (SCM_UNLIKELY (y == 0))
+    scm_num_overflow (s_scm_round_divide);  /* or return a NaN? */
+  else
+    {
+      double q = scm_c_round (x / y);
+      double 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_round_divide (SCM x, SCM y)
+{
+  SCM q, r, r2;
+  int cmp, needs_adjustment;
+
+  /* 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 ();
+  r2 = scm_i_mkbig ();
+
+  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_1 (x);
+  mpz_mul_2exp (SCM_I_BIG_MPZ (r2), SCM_I_BIG_MPZ (r), 1);  /* r2 = 2*r */
+
+  cmp = mpz_cmpabs (SCM_I_BIG_MPZ (r2), SCM_I_BIG_MPZ (y));
+  if (mpz_odd_p (SCM_I_BIG_MPZ (q)))
+    needs_adjustment = (cmp >= 0);
+  else
+    needs_adjustment = (cmp > 0);
+
+  if (needs_adjustment)
+    {
+      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 (r2, y);
+  return scm_values (scm_list_2 (scm_i_normbig (q),
+				 scm_i_normbig (r)));
+}
+
+/* Compute exact round quotient and remainder 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_round_divide (SCM x, SCM y)
+{
+  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
+    SCM_WTA_DISPATCH_2 (g_scm_round_divide, x, y, SCM_ARG1,
+			s_scm_round_divide);
+  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
+    SCM_WTA_DISPATCH_2 (g_scm_round_divide, x, y, SCM_ARG2,
+			s_scm_round_divide);
+  else if (scm_is_true (scm_zero_p (y)))
+    scm_num_overflow (s_scm_round_divide);
+  else
+    {
+      SCM q = scm_round_number (scm_divide (x, y));
+      SCM 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),
@@ -6727,11 +8907,11 @@ SCM_DEFINE (scm_round_number, "round", 1, 0, 0,
     return x;
   else if (SCM_REALP (x))
     return scm_from_double (scm_c_round (SCM_REAL_VALUE (x)));
+  else if (SCM_FRACTIONP (x))
+    return scm_round_quotient (SCM_FRACTION_NUMERATOR (x),
+			       SCM_FRACTION_DENOMINATOR (x));
   else
     {
-      /* OPTIMIZE-ME: Fraction case could be done more efficiently by a
-         single quotient+remainder division then examining to see which way
-         the rounding should go.  */
       SCM plus_half = scm_sum (x, exactly_one_half);
       SCM result = scm_floor (plus_half);
       /* Adjust so that the rounding is towards even.  */
@@ -6754,22 +8934,8 @@ SCM_PRIMITIVE_GENERIC (scm_floor, "floor", 1, 0, 0,
   else if (SCM_REALP (x))
     return scm_from_double (floor (SCM_REAL_VALUE (x)));
   else if (SCM_FRACTIONP (x))
-    {
-      SCM q = scm_quotient (SCM_FRACTION_NUMERATOR (x),
-			    SCM_FRACTION_DENOMINATOR (x));
-      if (scm_is_false (scm_negative_p (x)))
-	{
-	  /* For positive x, rounding towards zero is correct. */
-	  return q;
-	}
-      else
-	{
-	  /* For negative x, we need to return q-1 unless x is an
-	     integer.  But fractions are never integer, per our
-	     assumptions. */
-	  return scm_difference (q, SCM_INUM1);
-	}
-    }
+    return scm_floor_quotient (SCM_FRACTION_NUMERATOR (x),
+			       SCM_FRACTION_DENOMINATOR (x));
   else
     SCM_WTA_DISPATCH_1 (g_scm_floor, x, 1, s_scm_floor);
 }  
@@ -6785,22 +8951,8 @@ SCM_PRIMITIVE_GENERIC (scm_ceiling, "ceiling", 1, 0, 0,
   else if (SCM_REALP (x))
     return scm_from_double (ceil (SCM_REAL_VALUE (x)));
   else if (SCM_FRACTIONP (x))
-    {
-      SCM q = scm_quotient (SCM_FRACTION_NUMERATOR (x),
-			    SCM_FRACTION_DENOMINATOR (x));
-      if (scm_is_false (scm_positive_p (x)))
-	{
-	  /* For negative x, rounding towards zero is correct. */
-	  return q;
-	}
-      else
-	{
-	  /* For positive x, we need to return q+1 unless x is an
-	     integer.  But fractions are never integer, per our
-	     assumptions. */
-	  return scm_sum (q, SCM_INUM1);
-	}
-    }
+    return scm_ceiling_quotient (SCM_FRACTION_NUMERATOR (x),
+				 SCM_FRACTION_DENOMINATOR (x));
   else
     SCM_WTA_DISPATCH_1 (g_scm_ceiling, x, 1, s_scm_ceiling);
 }
diff --git a/libguile/numbers.h b/libguile/numbers.h
index 10a4f17..66cbe63 100644
--- a/libguile/numbers.h
+++ b/libguile/numbers.h
@@ -181,9 +181,21 @@ SCM_API SCM scm_modulo (SCM x, SCM y);
 SCM_API SCM scm_euclidean_divide (SCM x, SCM y);
 SCM_API SCM scm_euclidean_quotient (SCM x, SCM y);
 SCM_API SCM scm_euclidean_remainder (SCM x, SCM y);
+SCM_API SCM scm_floor_divide (SCM x, SCM y);
+SCM_API SCM scm_floor_quotient (SCM x, SCM y);
+SCM_API SCM scm_floor_remainder (SCM x, SCM y);
+SCM_API SCM scm_ceiling_divide (SCM x, SCM y);
+SCM_API SCM scm_ceiling_quotient (SCM x, SCM y);
+SCM_API SCM scm_ceiling_remainder (SCM x, SCM y);
+SCM_API SCM scm_truncate_divide (SCM x, SCM y);
+SCM_API SCM scm_truncate_quotient (SCM x, SCM y);
+SCM_API SCM scm_truncate_remainder (SCM x, SCM y);
 SCM_API SCM scm_centered_divide (SCM x, SCM y);
 SCM_API SCM scm_centered_quotient (SCM x, SCM y);
 SCM_API SCM scm_centered_remainder (SCM x, SCM y);
+SCM_API SCM scm_round_divide (SCM x, SCM y);
+SCM_API SCM scm_round_quotient (SCM x, SCM y);
+SCM_API SCM scm_round_remainder (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/test-suite/tests/numbers.test b/test-suite/tests/numbers.test
index f738189..ef59a02 100644
--- a/test-suite/tests/numbers.test
+++ b/test-suite/tests/numbers.test
@@ -4148,6 +4148,42 @@
                   (test-within-range? 0 <= r < (abs y)))
              (test-eqv? q (/ x y)))))
 
+  (define (valid-floor-answer? x y q r)
+    (and (eq? (exact? q)
+              (exact? r)
+              (and (exact? x) (exact? y)))
+         (test-eqv? r (- x (* q y)))
+         (if (and (finite? x) (finite? y))
+             (and (integer? q)
+                  (if (> y 0)
+                      (test-within-range? 0 <= r < y)
+                      (test-within-range? y < r <= 0)))
+             (test-eqv? q (/ x y)))))
+
+  (define (valid-ceiling-answer? x y q r)
+    (and (eq? (exact? q)
+              (exact? r)
+              (and (exact? x) (exact? y)))
+         (test-eqv? r (- x (* q y)))
+         (if (and (finite? x) (finite? y))
+             (and (integer? q)
+                  (if (> y 0)
+                      (test-within-range? (- y) < r <= 0)
+                      (test-within-range?    0 <= r < (- y))))
+             (test-eqv? q (/ x y)))))
+
+  (define (valid-truncate-answer? x y q r)
+    (and (eq? (exact? q)
+              (exact? r)
+              (and (exact? x) (exact? y)))
+         (test-eqv? r (- x (* q y)))
+         (if (and (finite? x) (finite? y))
+             (and (integer? q)
+                  (if (> x 0)
+                      (test-within-range?          0 <= r < (abs y))
+                      (test-within-range? (- (abs y)) < r <= 0)))
+             (test-eqv? q (/ x y)))))
+
   (define (valid-centered-answer? x y q r)
     (and (eq? (exact? q)
               (exact? r)
@@ -4159,6 +4195,19 @@
                    (* -1/2 (abs y)) <= r < (* +1/2 (abs y))))
              (test-eqv? q (/ x y)))))
 
+  (define (valid-round-answer? x y q r)
+    (and (eq? (exact? q)
+              (exact? r)
+              (and (exact? x) (exact? y)))
+         (test-eqv? r (- x (* q y)))
+         (if (and (finite? x) (finite? y))
+             (and (integer? q)
+                  (let ((ay/2 (/ (abs y) 2)))
+                    (if (even? q)
+                        (test-within-range? (- ay/2) <= r <= ay/2)
+                        (test-within-range? (- ay/2) <  r <  ay/2))))
+             (test-eqv? q (/ x y)))))
+
   (define (for lsts f) (apply for-each f lsts))
 
   (define big (expt 10 (1+ (inexact->exact (ceiling (log10 fixnum-max))))))
@@ -4284,8 +4333,32 @@
                         euclidean-remainder
                         valid-euclidean-answer?))
 
+  (with-test-prefix "floor/"
+    (run-division-tests floor/
+                        floor-quotient
+                        floor-remainder
+                        valid-floor-answer?))
+
+  (with-test-prefix "ceiling/"
+    (run-division-tests ceiling/
+                        ceiling-quotient
+                        ceiling-remainder
+                        valid-ceiling-answer?))
+
+  (with-test-prefix "truncate/"
+    (run-division-tests truncate/
+                        truncate-quotient
+                        truncate-remainder
+                        valid-truncate-answer?))
+
   (with-test-prefix "centered/"
     (run-division-tests centered/
                         centered-quotient
                         centered-remainder
-                        valid-centered-answer?)))
+                        valid-centered-answer?))
+
+  (with-test-prefix "round/"
+    (run-division-tests round/
+                        round-quotient
+                        round-remainder
+                        valid-round-answer?)))
-- 
1.5.6.5


[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #3: Added internal C function to extract from values object --]
[-- Type: text/x-diff, Size: 1799 bytes --]

From a729ef1b8b1905cb6c2ec73cc9b54e2e9453a2b3 Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Thu, 10 Feb 2011 18:03:14 -0500
Subject: [PATCH] Added internal C function to extract from values object

* libguile/values.c (scm_i_extract_values_2): New internal function
  that extracts two values from a values object.

* libguile/values.h: Added prototype.
---
 libguile/values.c |   18 ++++++++++++++++++
 libguile/values.h |    2 ++
 2 files changed, 20 insertions(+), 0 deletions(-)

diff --git a/libguile/values.c b/libguile/values.c
index 967fcd6..b9fe4cc 100644
--- a/libguile/values.c
+++ b/libguile/values.c
@@ -35,6 +35,24 @@
 
 SCM scm_values_vtable;
 
+/* OBJ must be a values object containing exactly two values.
+   scm_i_extract_values_2 puts those two values into *p1 and *p2.  */
+void
+scm_i_extract_values_2 (SCM obj, SCM *p1, SCM *p2)
+{
+  SCM values;
+
+  SCM_ASSERT_TYPE (SCM_VALUESP (obj), obj, SCM_ARG1,
+		   "scm_i_extract_values_2", "values");
+  values = scm_struct_ref (obj, SCM_INUM0);
+  if (!scm_is_null_or_nil (SCM_CDDR (values)))
+    scm_wrong_type_arg_msg
+      ("scm_i_extract_values_2", SCM_ARG1, obj,
+       "a values object containing exactly two values");
+  *p1 = SCM_CAR (values);
+  *p2 = SCM_CADR (values);
+}
+
 static SCM
 print_values (SCM obj, SCM pwps)
 {
diff --git a/libguile/values.h b/libguile/values.h
index 0750aec..65ad8a1 100644
--- a/libguile/values.h
+++ b/libguile/values.h
@@ -30,6 +30,8 @@ SCM_API SCM scm_values_vtable;
 #define SCM_VALUESP(x) (SCM_STRUCTP (x)\
                         && scm_is_eq (scm_struct_vtable (x), scm_values_vtable))
 
+SCM_INTERNAL void scm_i_extract_values_2 (SCM obj, SCM *p1, SCM *p2);
+
 SCM_API SCM scm_values (SCM args);
 SCM_INTERNAL void scm_init_values (void);
 
-- 
1.5.6.5


[-- Warning: decoded text below may be mangled, UTF-8 assumed --]
[-- Attachment #4: Optimize division operators handling of fractions --]
[-- Type: text/x-diff, Size: 54385 bytes --]

From 5a0677e5daafa6dd10dd7612468d1ea1e38692b7 Mon Sep 17 00:00:00 2001
From: Mark H Weaver <mhw@netris.org>
Date: Thu, 10 Feb 2011 18:18:34 -0500
Subject: [PATCH] Optimize division operators handling of fractions

* libguile/numbers.c: (scm_euclidean_quotient, scm_euclidean_remainder,
  scm_euclidean_divide, scm_floor_quotient, scm_floor_remainder,
  scm_floor_divide, scm_ceiling_quotient, scm_ceiling_remainder,
  scm_ceiling_divide, scm_truncate_quotient, scm_truncate_remainder,
  scm_truncate_divide, scm_centered_quotient, scm_centered_remainder,
  scm_centered_divide, scm_round_quotient, scm_round_remainder,
  scm_round_divide): Optimize case where both arguments are exact and at
  least one is a fraction, by reducing to a subproblem involving only
  integers, and then adjusting the resulting remainder as needed.
---
 libguile/numbers.c |  622 +++++++++++++++++++++-------------------------------
 1 files changed, 249 insertions(+), 373 deletions(-)

diff --git a/libguile/numbers.c b/libguile/numbers.c
index 90e449f..3ef8fe5 100644
--- a/libguile/numbers.c
+++ b/libguile/numbers.c
@@ -1070,7 +1070,7 @@ SCM_PRIMITIVE_GENERIC (scm_modulo, "modulo", 2, 0, 0,
 #undef FUNC_NAME
 
 static SCM scm_i_inexact_euclidean_quotient (double x, double y);
-static SCM scm_i_slow_exact_euclidean_quotient (SCM x, SCM y);
+static SCM scm_i_exact_rational_euclidean_quotient (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_euclidean_quotient, "euclidean-quotient", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -1125,7 +1125,7 @@ SCM_PRIMITIVE_GENERIC (scm_euclidean_quotient, "euclidean-quotient", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_euclidean_quotient (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_euclidean_quotient (x, y);
+	return scm_i_exact_rational_euclidean_quotient (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_euclidean_quotient, x, y, SCM_ARG2,
 			    s_scm_euclidean_quotient);
@@ -1171,7 +1171,7 @@ SCM_PRIMITIVE_GENERIC (scm_euclidean_quotient, "euclidean-quotient", 2, 0, 0,
 	return scm_i_inexact_euclidean_quotient
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_euclidean_quotient (x, y);
+	return scm_i_exact_rational_euclidean_quotient (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_euclidean_quotient, x, y, SCM_ARG2,
 			    s_scm_euclidean_quotient);
@@ -1191,8 +1191,11 @@ SCM_PRIMITIVE_GENERIC (scm_euclidean_quotient, "euclidean-quotient", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_euclidean_quotient
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_euclidean_quotient (x, y);
       else
-	return scm_i_slow_exact_euclidean_quotient (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_euclidean_quotient, x, y, SCM_ARG2,
+			    s_scm_euclidean_quotient);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_euclidean_quotient, x, y, SCM_ARG1,
@@ -1213,28 +1216,16 @@ scm_i_inexact_euclidean_quotient (double x, double y)
     return scm_nan ();
 }
 
-/* Compute exact euclidean_quotient 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_euclidean_quotient (SCM x, SCM y)
+scm_i_exact_rational_euclidean_quotient (SCM x, SCM y)
 {
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_euclidean_quotient, x, y, SCM_ARG1,
-			s_scm_euclidean_quotient);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_euclidean_quotient, x, y, SCM_ARG2,
-			s_scm_euclidean_quotient);
-  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_scm_euclidean_quotient);
+  return scm_euclidean_quotient
+    (scm_product (scm_numerator (x), scm_denominator (y)),
+     scm_product (scm_numerator (y), scm_denominator (x)));
 }
 
 static SCM scm_i_inexact_euclidean_remainder (double x, double y);
-static SCM scm_i_slow_exact_euclidean_remainder (SCM x, SCM y);
+static SCM scm_i_exact_rational_euclidean_remainder (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_euclidean_remainder, "euclidean-remainder", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -1294,7 +1285,7 @@ SCM_PRIMITIVE_GENERIC (scm_euclidean_remainder, "euclidean-remainder", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_euclidean_remainder (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_euclidean_remainder (x, y);
+	return scm_i_exact_rational_euclidean_remainder (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_euclidean_remainder, x, y, SCM_ARG2,
 			    s_scm_euclidean_remainder);
@@ -1329,7 +1320,7 @@ SCM_PRIMITIVE_GENERIC (scm_euclidean_remainder, "euclidean-remainder", 2, 0, 0,
 	return scm_i_inexact_euclidean_remainder
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_euclidean_remainder (x, y);
+	return scm_i_exact_rational_euclidean_remainder (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_euclidean_remainder, x, y, SCM_ARG2,
 			    s_scm_euclidean_remainder);
@@ -1349,8 +1340,11 @@ SCM_PRIMITIVE_GENERIC (scm_euclidean_remainder, "euclidean-remainder", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_euclidean_remainder
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_euclidean_remainder (x, y);
       else
-	return scm_i_slow_exact_euclidean_remainder (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_euclidean_remainder, x, y, SCM_ARG2,
+			    s_scm_euclidean_remainder);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_euclidean_remainder, x, y, SCM_ARG1,
@@ -1383,32 +1377,19 @@ scm_i_inexact_euclidean_remainder (double x, double y)
   return scm_from_double (x - q * y);
 }
 
-/* Compute exact euclidean_remainder 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_euclidean_remainder (SCM x, SCM y)
+scm_i_exact_rational_euclidean_remainder (SCM x, SCM y)
 {
-  SCM q;
-
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_euclidean_remainder, x, y, SCM_ARG1,
-			s_scm_euclidean_remainder);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_euclidean_remainder, x, y, SCM_ARG2,
-			s_scm_euclidean_remainder);
-  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_scm_euclidean_remainder);
-  return scm_difference (x, scm_product (y, q));
+  SCM xd = scm_denominator (x);
+  SCM yd = scm_denominator (y);
+  SCM r1 = scm_euclidean_remainder (scm_product (scm_numerator (x), yd),
+				    scm_product (scm_numerator (y), xd));
+  return scm_divide (r1, scm_product (xd, yd));
 }
 
 
 static SCM scm_i_inexact_euclidean_divide (double x, double y);
-static SCM scm_i_slow_exact_euclidean_divide (SCM x, SCM y);
+static SCM scm_i_exact_rational_euclidean_divide (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_euclidean_divide, "euclidean/", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -1477,7 +1458,7 @@ SCM_PRIMITIVE_GENERIC (scm_euclidean_divide, "euclidean/", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_euclidean_divide (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_euclidean_divide (x, y);
+	return scm_i_exact_rational_euclidean_divide (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_euclidean_divide, x, y, SCM_ARG2,
 			    s_scm_euclidean_divide);
@@ -1525,7 +1506,7 @@ SCM_PRIMITIVE_GENERIC (scm_euclidean_divide, "euclidean/", 2, 0, 0,
 	return scm_i_inexact_euclidean_divide
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_euclidean_divide (x, y);
+	return scm_i_exact_rational_euclidean_divide (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_euclidean_divide, x, y, SCM_ARG2,
 			    s_scm_euclidean_divide);
@@ -1545,8 +1526,11 @@ SCM_PRIMITIVE_GENERIC (scm_euclidean_divide, "euclidean/", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_euclidean_divide
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_euclidean_divide (x, y);
       else
-	return scm_i_slow_exact_euclidean_divide (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_euclidean_divide, x, y, SCM_ARG2,
+			    s_scm_euclidean_divide);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_euclidean_divide, x, y, SCM_ARG1,
@@ -1572,32 +1556,23 @@ scm_i_inexact_euclidean_divide (double x, double y)
 				 scm_from_double (r)));
 }
 
-/* Compute exact euclidean quotient and remainder 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_euclidean_divide (SCM x, SCM y)
+scm_i_exact_rational_euclidean_divide (SCM x, SCM y)
 {
-  SCM q, r;
+  SCM q, r, r1;
+  SCM xd = scm_denominator (x);
+  SCM yd = scm_denominator (y);
 
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_euclidean_divide, x, y, SCM_ARG1,
-			s_scm_euclidean_divide);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_euclidean_divide, x, y, SCM_ARG2,
-			s_scm_euclidean_divide);
-  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_scm_euclidean_divide);
-  r = scm_difference (x, scm_product (q, y));
+  scm_i_extract_values_2
+    (scm_euclidean_divide (scm_product (scm_numerator (x), yd),
+			   scm_product (scm_numerator (y), xd)),
+     &q, &r1);
+  r = scm_divide (r1, scm_product (xd, yd));
   return scm_values (scm_list_2 (q, r));
 }
 
 static SCM scm_i_inexact_floor_quotient (double x, double y);
-static SCM scm_i_slow_exact_floor_quotient (SCM x, SCM y);
+static SCM scm_i_exact_rational_floor_quotient (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_floor_quotient, "floor-quotient", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -1647,7 +1622,7 @@ SCM_PRIMITIVE_GENERIC (scm_floor_quotient, "floor-quotient", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_floor_quotient (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_floor_quotient (x, y);
+	return scm_i_exact_rational_floor_quotient (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_floor_quotient, x, y, SCM_ARG2,
 			    s_scm_floor_quotient);
@@ -1688,7 +1663,7 @@ SCM_PRIMITIVE_GENERIC (scm_floor_quotient, "floor-quotient", 2, 0, 0,
 	return scm_i_inexact_floor_quotient
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_floor_quotient (x, y);
+	return scm_i_exact_rational_floor_quotient (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_floor_quotient, x, y, SCM_ARG2,
 			    s_scm_floor_quotient);
@@ -1708,8 +1683,11 @@ SCM_PRIMITIVE_GENERIC (scm_floor_quotient, "floor-quotient", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_floor_quotient
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_floor_quotient (x, y);
       else
-	return scm_i_slow_exact_floor_quotient (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_floor_quotient, x, y, SCM_ARG2,
+			    s_scm_floor_quotient);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_floor_quotient, x, y, SCM_ARG1,
@@ -1726,26 +1704,16 @@ scm_i_inexact_floor_quotient (double x, double y)
     return scm_from_double (floor (x / y));
 }
 
-/* Compute exact floor_quotient 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_floor_quotient (SCM x, SCM y)
+scm_i_exact_rational_floor_quotient (SCM x, SCM y)
 {
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_floor_quotient, x, y, SCM_ARG1,
-			s_scm_floor_quotient);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_floor_quotient, x, y, SCM_ARG2,
-			s_scm_floor_quotient);
-  else if (scm_is_true (scm_zero_p (y)))
-    scm_num_overflow (s_scm_floor_quotient);
-  else
-    return scm_floor (scm_divide (x, y));
+  return scm_floor_quotient
+    (scm_product (scm_numerator (x), scm_denominator (y)),
+     scm_product (scm_numerator (y), scm_denominator (x)));
 }
 
 static SCM scm_i_inexact_floor_remainder (double x, double y);
-static SCM scm_i_slow_exact_floor_remainder (SCM x, SCM y);
+static SCM scm_i_exact_rational_floor_remainder (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_floor_remainder, "floor-remainder", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -1814,7 +1782,7 @@ SCM_PRIMITIVE_GENERIC (scm_floor_remainder, "floor-remainder", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_floor_remainder (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_floor_remainder (x, y);
+	return scm_i_exact_rational_floor_remainder (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_floor_remainder, x, y, SCM_ARG2,
 			    s_scm_floor_remainder);
@@ -1850,7 +1818,7 @@ SCM_PRIMITIVE_GENERIC (scm_floor_remainder, "floor-remainder", 2, 0, 0,
 	return scm_i_inexact_floor_remainder
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_floor_remainder (x, y);
+	return scm_i_exact_rational_floor_remainder (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_floor_remainder, x, y, SCM_ARG2,
 			    s_scm_floor_remainder);
@@ -1870,8 +1838,11 @@ SCM_PRIMITIVE_GENERIC (scm_floor_remainder, "floor-remainder", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_floor_remainder
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_floor_remainder (x, y);
       else
-	return scm_i_slow_exact_floor_remainder (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_floor_remainder, x, y, SCM_ARG2,
+			    s_scm_floor_remainder);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_floor_remainder, x, y, SCM_ARG1,
@@ -1896,28 +1867,19 @@ scm_i_inexact_floor_remainder (double x, double y)
     return scm_from_double (x - y * floor (x / y));
 }
 
-/* Compute exact floor_remainder 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_floor_remainder (SCM x, SCM y)
+scm_i_exact_rational_floor_remainder (SCM x, SCM y)
 {
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_floor_remainder, x, y, SCM_ARG1,
-			s_scm_floor_remainder);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_floor_remainder, x, y, SCM_ARG2,
-			s_scm_floor_remainder);
-  else if (scm_is_true (scm_zero_p (y)))
-    scm_num_overflow (s_scm_floor_remainder);
-  else
-    return scm_difference
-      (x, scm_product (y, scm_floor (scm_divide (x, y))));
+  SCM xd = scm_denominator (x);
+  SCM yd = scm_denominator (y);
+  SCM r1 = scm_floor_remainder (scm_product (scm_numerator (x), yd),
+				scm_product (scm_numerator (y), xd));
+  return scm_divide (r1, scm_product (xd, yd));
 }
 
 
 static SCM scm_i_inexact_floor_divide (double x, double y);
-static SCM scm_i_slow_exact_floor_divide (SCM x, SCM y);
+static SCM scm_i_exact_rational_floor_divide (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_floor_divide, "floor/", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -1998,7 +1960,7 @@ SCM_PRIMITIVE_GENERIC (scm_floor_divide, "floor/", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_floor_divide (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_floor_divide (x, y);
+	return scm_i_exact_rational_floor_divide (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_floor_divide, x, y, SCM_ARG2,
 			    s_scm_floor_divide);
@@ -2042,7 +2004,7 @@ SCM_PRIMITIVE_GENERIC (scm_floor_divide, "floor/", 2, 0, 0,
 	return scm_i_inexact_floor_divide
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_floor_divide (x, y);
+	return scm_i_exact_rational_floor_divide (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_floor_divide, x, y, SCM_ARG2,
 			    s_scm_floor_divide);
@@ -2062,8 +2024,11 @@ SCM_PRIMITIVE_GENERIC (scm_floor_divide, "floor/", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_floor_divide
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_floor_divide (x, y);
       else
-	return scm_i_slow_exact_floor_divide (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_floor_divide, x, y, SCM_ARG2,
+			    s_scm_floor_divide);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_floor_divide, x, y, SCM_ARG1,
@@ -2085,30 +2050,23 @@ scm_i_inexact_floor_divide (double x, double y)
 				 scm_from_double (r)));
 }
 
-/* Compute exact floor quotient and remainder 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_floor_divide (SCM x, SCM y)
+scm_i_exact_rational_floor_divide (SCM x, SCM y)
 {
-  SCM q, r;
+  SCM q, r, r1;
+  SCM xd = scm_denominator (x);
+  SCM yd = scm_denominator (y);
 
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_floor_divide, x, y, SCM_ARG1,
-			s_scm_floor_divide);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_floor_divide, x, y, SCM_ARG2,
-			s_scm_floor_divide);
-  else if (scm_is_true (scm_zero_p (y)))
-    scm_num_overflow (s_scm_floor_divide);
-  else
-    q = scm_floor (scm_divide (x, y));
-  r = scm_difference (x, scm_product (q, y));
+  scm_i_extract_values_2
+    (scm_floor_divide (scm_product (scm_numerator (x), yd),
+		       scm_product (scm_numerator (y), xd)),
+     &q, &r1);
+  r = scm_divide (r1, scm_product (xd, yd));
   return scm_values (scm_list_2 (q, r));
 }
 
 static SCM scm_i_inexact_ceiling_quotient (double x, double y);
-static SCM scm_i_slow_exact_ceiling_quotient (SCM x, SCM y);
+static SCM scm_i_exact_rational_ceiling_quotient (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_ceiling_quotient, "ceiling-quotient", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -2178,7 +2136,7 @@ SCM_PRIMITIVE_GENERIC (scm_ceiling_quotient, "ceiling-quotient", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_ceiling_quotient (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_ceiling_quotient (x, y);
+	return scm_i_exact_rational_ceiling_quotient (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_ceiling_quotient, x, y, SCM_ARG2,
 			    s_scm_ceiling_quotient);
@@ -2219,7 +2177,7 @@ SCM_PRIMITIVE_GENERIC (scm_ceiling_quotient, "ceiling-quotient", 2, 0, 0,
 	return scm_i_inexact_ceiling_quotient
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_ceiling_quotient (x, y);
+	return scm_i_exact_rational_ceiling_quotient (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_ceiling_quotient, x, y, SCM_ARG2,
 			    s_scm_ceiling_quotient);
@@ -2239,8 +2197,11 @@ SCM_PRIMITIVE_GENERIC (scm_ceiling_quotient, "ceiling-quotient", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_ceiling_quotient
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_ceiling_quotient (x, y);
       else
-	return scm_i_slow_exact_ceiling_quotient (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_ceiling_quotient, x, y, SCM_ARG2,
+			    s_scm_ceiling_quotient);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_ceiling_quotient, x, y, SCM_ARG1,
@@ -2257,26 +2218,16 @@ scm_i_inexact_ceiling_quotient (double x, double y)
     return scm_from_double (ceil (x / y));
 }
 
-/* Compute exact ceiling_quotient 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_ceiling_quotient (SCM x, SCM y)
+scm_i_exact_rational_ceiling_quotient (SCM x, SCM y)
 {
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_ceiling_quotient, x, y, SCM_ARG1,
-			s_scm_ceiling_quotient);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_ceiling_quotient, x, y, SCM_ARG2,
-			s_scm_ceiling_quotient);
-  else if (scm_is_true (scm_zero_p (y)))
-    scm_num_overflow (s_scm_ceiling_quotient);
-  else
-    return scm_ceiling (scm_divide (x, y));
+  return scm_ceiling_quotient
+    (scm_product (scm_numerator (x), scm_denominator (y)),
+     scm_product (scm_numerator (y), scm_denominator (x)));
 }
 
 static SCM scm_i_inexact_ceiling_remainder (double x, double y);
-static SCM scm_i_slow_exact_ceiling_remainder (SCM x, SCM y);
+static SCM scm_i_exact_rational_ceiling_remainder (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_ceiling_remainder, "ceiling-remainder", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -2355,7 +2306,7 @@ SCM_PRIMITIVE_GENERIC (scm_ceiling_remainder, "ceiling-remainder", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_ceiling_remainder (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_ceiling_remainder (x, y);
+	return scm_i_exact_rational_ceiling_remainder (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_ceiling_remainder, x, y, SCM_ARG2,
 			    s_scm_ceiling_remainder);
@@ -2391,7 +2342,7 @@ SCM_PRIMITIVE_GENERIC (scm_ceiling_remainder, "ceiling-remainder", 2, 0, 0,
 	return scm_i_inexact_ceiling_remainder
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_ceiling_remainder (x, y);
+	return scm_i_exact_rational_ceiling_remainder (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_ceiling_remainder, x, y, SCM_ARG2,
 			    s_scm_ceiling_remainder);
@@ -2411,8 +2362,11 @@ SCM_PRIMITIVE_GENERIC (scm_ceiling_remainder, "ceiling-remainder", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_ceiling_remainder
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_ceiling_remainder (x, y);
       else
-	return scm_i_slow_exact_ceiling_remainder (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_ceiling_remainder, x, y, SCM_ARG2,
+			    s_scm_ceiling_remainder);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_ceiling_remainder, x, y, SCM_ARG1,
@@ -2437,27 +2391,18 @@ scm_i_inexact_ceiling_remainder (double x, double y)
     return scm_from_double (x - y * ceil (x / y));
 }
 
-/* Compute exact ceiling_remainder 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_ceiling_remainder (SCM x, SCM y)
+scm_i_exact_rational_ceiling_remainder (SCM x, SCM y)
 {
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_ceiling_remainder, x, y, SCM_ARG1,
-			s_scm_ceiling_remainder);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_ceiling_remainder, x, y, SCM_ARG2,
-			s_scm_ceiling_remainder);
-  else if (scm_is_true (scm_zero_p (y)))
-    scm_num_overflow (s_scm_ceiling_remainder);
-  else
-    return scm_difference
-      (x, scm_product (y, scm_ceiling (scm_divide (x, y))));
+  SCM xd = scm_denominator (x);
+  SCM yd = scm_denominator (y);
+  SCM r1 = scm_ceiling_remainder (scm_product (scm_numerator (x), yd),
+				  scm_product (scm_numerator (y), xd));
+  return scm_divide (r1, scm_product (xd, yd));
 }
 
 static SCM scm_i_inexact_ceiling_divide (double x, double y);
-static SCM scm_i_slow_exact_ceiling_divide (SCM x, SCM y);
+static SCM scm_i_exact_rational_ceiling_divide (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_ceiling_divide, "ceiling/", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -2548,7 +2493,7 @@ SCM_PRIMITIVE_GENERIC (scm_ceiling_divide, "ceiling/", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_ceiling_divide (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_ceiling_divide (x, y);
+	return scm_i_exact_rational_ceiling_divide (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_ceiling_divide, x, y, SCM_ARG2,
 			    s_scm_ceiling_divide);
@@ -2592,7 +2537,7 @@ SCM_PRIMITIVE_GENERIC (scm_ceiling_divide, "ceiling/", 2, 0, 0,
 	return scm_i_inexact_ceiling_divide
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_ceiling_divide (x, y);
+	return scm_i_exact_rational_ceiling_divide (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_ceiling_divide, x, y, SCM_ARG2,
 			    s_scm_ceiling_divide);
@@ -2612,8 +2557,11 @@ SCM_PRIMITIVE_GENERIC (scm_ceiling_divide, "ceiling/", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_ceiling_divide
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_ceiling_divide (x, y);
       else
-	return scm_i_slow_exact_ceiling_divide (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_ceiling_divide, x, y, SCM_ARG2,
+			    s_scm_ceiling_divide);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_ceiling_divide, x, y, SCM_ARG1,
@@ -2635,30 +2583,23 @@ scm_i_inexact_ceiling_divide (double x, double y)
 				 scm_from_double (r)));
 }
 
-/* Compute exact ceiling quotient and remainder 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_ceiling_divide (SCM x, SCM y)
+scm_i_exact_rational_ceiling_divide (SCM x, SCM y)
 {
-  SCM q, r;
+  SCM q, r, r1;
+  SCM xd = scm_denominator (x);
+  SCM yd = scm_denominator (y);
 
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_ceiling_divide, x, y, SCM_ARG1,
-			s_scm_ceiling_divide);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_ceiling_divide, x, y, SCM_ARG2,
-			s_scm_ceiling_divide);
-  else if (scm_is_true (scm_zero_p (y)))
-    scm_num_overflow (s_scm_ceiling_divide);
-  else
-    q = scm_ceiling (scm_divide (x, y));
-  r = scm_difference (x, scm_product (q, y));
+  scm_i_extract_values_2
+    (scm_ceiling_divide (scm_product (scm_numerator (x), yd),
+			 scm_product (scm_numerator (y), xd)),
+     &q, &r1);
+  r = scm_divide (r1, scm_product (xd, yd));
   return scm_values (scm_list_2 (q, r));
 }
 
 static SCM scm_i_inexact_truncate_quotient (double x, double y);
-static SCM scm_i_slow_exact_truncate_quotient (SCM x, SCM y);
+static SCM scm_i_exact_rational_truncate_quotient (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_truncate_quotient, "truncate-quotient", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -2706,7 +2647,7 @@ SCM_PRIMITIVE_GENERIC (scm_truncate_quotient, "truncate-quotient", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_truncate_quotient (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_truncate_quotient (x, y);
+	return scm_i_exact_rational_truncate_quotient (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_truncate_quotient, x, y, SCM_ARG2,
 			    s_scm_truncate_quotient);
@@ -2747,7 +2688,7 @@ SCM_PRIMITIVE_GENERIC (scm_truncate_quotient, "truncate-quotient", 2, 0, 0,
 	return scm_i_inexact_truncate_quotient
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_truncate_quotient (x, y);
+	return scm_i_exact_rational_truncate_quotient (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_truncate_quotient, x, y, SCM_ARG2,
 			    s_scm_truncate_quotient);
@@ -2767,8 +2708,11 @@ SCM_PRIMITIVE_GENERIC (scm_truncate_quotient, "truncate-quotient", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_truncate_quotient
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_truncate_quotient (x, y);
       else
-	return scm_i_slow_exact_truncate_quotient (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_truncate_quotient, x, y, SCM_ARG2,
+			    s_scm_truncate_quotient);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_truncate_quotient, x, y, SCM_ARG1,
@@ -2785,26 +2729,16 @@ scm_i_inexact_truncate_quotient (double x, double y)
     return scm_from_double (trunc (x / y));
 }
 
-/* Compute exact truncate_quotient 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_truncate_quotient (SCM x, SCM y)
+scm_i_exact_rational_truncate_quotient (SCM x, SCM y)
 {
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_truncate_quotient, x, y, SCM_ARG1,
-			s_scm_truncate_quotient);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_truncate_quotient, x, y, SCM_ARG2,
-			s_scm_truncate_quotient);
-  else if (scm_is_true (scm_zero_p (y)))
-    scm_num_overflow (s_scm_truncate_quotient);
-  else
-    return scm_truncate_number (scm_divide (x, y));
+  return scm_truncate_quotient
+    (scm_product (scm_numerator (x), scm_denominator (y)),
+     scm_product (scm_numerator (y), scm_denominator (x)));
 }
 
 static SCM scm_i_inexact_truncate_remainder (double x, double y);
-static SCM scm_i_slow_exact_truncate_remainder (SCM x, SCM y);
+static SCM scm_i_exact_rational_truncate_remainder (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_truncate_remainder, "truncate-remainder", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -2848,7 +2782,7 @@ SCM_PRIMITIVE_GENERIC (scm_truncate_remainder, "truncate-remainder", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_truncate_remainder (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_truncate_remainder (x, y);
+	return scm_i_exact_rational_truncate_remainder (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_truncate_remainder, x, y, SCM_ARG2,
 			    s_scm_truncate_remainder);
@@ -2882,7 +2816,7 @@ SCM_PRIMITIVE_GENERIC (scm_truncate_remainder, "truncate-remainder", 2, 0, 0,
 	return scm_i_inexact_truncate_remainder
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_truncate_remainder (x, y);
+	return scm_i_exact_rational_truncate_remainder (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_truncate_remainder, x, y, SCM_ARG2,
 			    s_scm_truncate_remainder);
@@ -2902,8 +2836,11 @@ SCM_PRIMITIVE_GENERIC (scm_truncate_remainder, "truncate-remainder", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_truncate_remainder
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_truncate_remainder (x, y);
       else
-	return scm_i_slow_exact_truncate_remainder (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_truncate_remainder, x, y, SCM_ARG2,
+			    s_scm_truncate_remainder);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_truncate_remainder, x, y, SCM_ARG1,
@@ -2927,28 +2864,19 @@ scm_i_inexact_truncate_remainder (double x, double y)
     return scm_from_double (x - y * trunc (x / y));
 }
 
-/* Compute exact truncate_remainder 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_truncate_remainder (SCM x, SCM y)
+scm_i_exact_rational_truncate_remainder (SCM x, SCM y)
 {
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_truncate_remainder, x, y, SCM_ARG1,
-			s_scm_truncate_remainder);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_truncate_remainder, x, y, SCM_ARG2,
-			s_scm_truncate_remainder);
-  else if (scm_is_true (scm_zero_p (y)))
-    scm_num_overflow (s_scm_truncate_remainder);
-  else
-    return scm_difference
-      (x, scm_product (y, scm_truncate_number (scm_divide (x, y))));
+  SCM xd = scm_denominator (x);
+  SCM yd = scm_denominator (y);
+  SCM r1 = scm_truncate_remainder (scm_product (scm_numerator (x), yd),
+				   scm_product (scm_numerator (y), xd));
+  return scm_divide (r1, scm_product (xd, yd));
 }
 
 
 static SCM scm_i_inexact_truncate_divide (double x, double y);
-static SCM scm_i_slow_exact_truncate_divide (SCM x, SCM y);
+static SCM scm_i_exact_rational_truncate_divide (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_truncate_divide, "truncate/", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -3002,7 +2930,7 @@ SCM_PRIMITIVE_GENERIC (scm_truncate_divide, "truncate/", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_truncate_divide (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_truncate_divide (x, y);
+	return scm_i_exact_rational_truncate_divide (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_truncate_divide, x, y, SCM_ARG2,
 			    s_scm_truncate_divide);
@@ -3047,7 +2975,7 @@ SCM_PRIMITIVE_GENERIC (scm_truncate_divide, "truncate/", 2, 0, 0,
 	return scm_i_inexact_truncate_divide
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_truncate_divide (x, y);
+	return scm_i_exact_rational_truncate_divide (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_truncate_divide, x, y, SCM_ARG2,
 			    s_scm_truncate_divide);
@@ -3067,8 +2995,11 @@ SCM_PRIMITIVE_GENERIC (scm_truncate_divide, "truncate/", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_truncate_divide
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_truncate_divide (x, y);
       else
-	return scm_i_slow_exact_truncate_divide (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_truncate_divide, x, y, SCM_ARG2,
+			    s_scm_truncate_divide);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_truncate_divide, x, y, SCM_ARG1,
@@ -3090,31 +3021,24 @@ scm_i_inexact_truncate_divide (double x, double y)
 				 scm_from_double (r)));
 }
 
-/* Compute exact truncate quotient and remainder 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_truncate_divide (SCM x, SCM y)
+scm_i_exact_rational_truncate_divide (SCM x, SCM y)
 {
-  SCM q, r;
+  SCM q, r, r1;
+  SCM xd = scm_denominator (x);
+  SCM yd = scm_denominator (y);
 
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_truncate_divide, x, y, SCM_ARG1,
-			s_scm_truncate_divide);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_truncate_divide, x, y, SCM_ARG2,
-			s_scm_truncate_divide);
-  else if (scm_is_true (scm_zero_p (y)))
-    scm_num_overflow (s_scm_truncate_divide);
-  else
-    q = scm_truncate_number (scm_divide (x, y));
-  r = scm_difference (x, scm_product (q, y));
+  scm_i_extract_values_2
+    (scm_truncate_divide (scm_product (scm_numerator (x), yd),
+			  scm_product (scm_numerator (y), xd)),
+     &q, &r1);
+  r = scm_divide (r1, scm_product (xd, yd));
   return scm_values (scm_list_2 (q, r));
 }
 
 static SCM scm_i_inexact_centered_quotient (double x, double y);
 static SCM scm_i_bigint_centered_quotient (SCM x, SCM y);
-static SCM scm_i_slow_exact_centered_quotient (SCM x, SCM y);
+static SCM scm_i_exact_rational_centered_quotient (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_centered_quotient, "centered-quotient", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -3184,7 +3108,7 @@ SCM_PRIMITIVE_GENERIC (scm_centered_quotient, "centered-quotient", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_centered_quotient (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_centered_quotient (x, y);
+	return scm_i_exact_rational_centered_quotient (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_centered_quotient, x, y, SCM_ARG2,
 			    s_scm_centered_quotient);
@@ -3233,7 +3157,7 @@ SCM_PRIMITIVE_GENERIC (scm_centered_quotient, "centered-quotient", 2, 0, 0,
 	return scm_i_inexact_centered_quotient
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_centered_quotient (x, y);
+	return scm_i_exact_rational_centered_quotient (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_centered_quotient, x, y, SCM_ARG2,
 			    s_scm_centered_quotient);
@@ -3253,8 +3177,11 @@ SCM_PRIMITIVE_GENERIC (scm_centered_quotient, "centered-quotient", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_centered_quotient
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_centered_quotient (x, y);
       else
-	return scm_i_slow_exact_centered_quotient (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_centered_quotient, x, y, SCM_ARG2,
+			    s_scm_centered_quotient);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_centered_quotient, x, y, SCM_ARG1,
@@ -3318,31 +3245,17 @@ scm_i_bigint_centered_quotient (SCM x, SCM y)
   return scm_i_normbig (q);
 }
 
-/* Compute exact centered quotient 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_centered_quotient (SCM x, SCM y)
+scm_i_exact_rational_centered_quotient (SCM x, SCM y)
 {
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_centered_quotient, x, y, SCM_ARG1,
-			s_scm_centered_quotient);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_centered_quotient, x, y, SCM_ARG2,
-			s_scm_centered_quotient);
-  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_scm_centered_quotient);
+  return scm_centered_quotient
+    (scm_product (scm_numerator (x), scm_denominator (y)),
+     scm_product (scm_numerator (y), scm_denominator (x)));
 }
 
 static SCM scm_i_inexact_centered_remainder (double x, double y);
 static SCM scm_i_bigint_centered_remainder (SCM x, SCM y);
-static SCM scm_i_slow_exact_centered_remainder (SCM x, SCM y);
+static SCM scm_i_exact_rational_centered_remainder (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_centered_remainder, "centered-remainder", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -3409,7 +3322,7 @@ SCM_PRIMITIVE_GENERIC (scm_centered_remainder, "centered-remainder", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_centered_remainder (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_centered_remainder (x, y);
+	return scm_i_exact_rational_centered_remainder (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_centered_remainder, x, y, SCM_ARG2,
 			    s_scm_centered_remainder);
@@ -3450,7 +3363,7 @@ SCM_PRIMITIVE_GENERIC (scm_centered_remainder, "centered-remainder", 2, 0, 0,
 	return scm_i_inexact_centered_remainder
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_centered_remainder (x, y);
+	return scm_i_exact_rational_centered_remainder (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_centered_remainder, x, y, SCM_ARG2,
 			    s_scm_centered_remainder);
@@ -3470,8 +3383,11 @@ SCM_PRIMITIVE_GENERIC (scm_centered_remainder, "centered-remainder", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_centered_remainder
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_centered_remainder (x, y);
       else
-	return scm_i_slow_exact_centered_remainder (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_centered_remainder, x, y, SCM_ARG2,
+			    s_scm_centered_remainder);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_centered_remainder, x, y, SCM_ARG1,
@@ -3544,33 +3460,20 @@ scm_i_bigint_centered_remainder (SCM x, SCM y)
   return scm_i_normbig (r);
 }
 
-/* Compute exact centered_remainder 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_centered_remainder (SCM x, SCM y)
+scm_i_exact_rational_centered_remainder (SCM x, SCM y)
 {
-  SCM q;
-
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_centered_remainder, x, y, SCM_ARG1,
-			s_scm_centered_remainder);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_centered_remainder, x, y, SCM_ARG2,
-			s_scm_centered_remainder);
-  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_scm_centered_remainder);
-  return scm_difference (x, scm_product (y, q));
+  SCM xd = scm_denominator (x);
+  SCM yd = scm_denominator (y);
+  SCM r1 = scm_centered_remainder (scm_product (scm_numerator (x), yd),
+				   scm_product (scm_numerator (y), xd));
+  return scm_divide (r1, scm_product (xd, yd));
 }
 
 
 static SCM scm_i_inexact_centered_divide (double x, double y);
 static SCM scm_i_bigint_centered_divide (SCM x, SCM y);
-static SCM scm_i_slow_exact_centered_divide (SCM x, SCM y);
+static SCM scm_i_exact_rational_centered_divide (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_centered_divide, "centered/", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -3643,7 +3546,7 @@ SCM_PRIMITIVE_GENERIC (scm_centered_divide, "centered/", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_centered_divide (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_centered_divide (x, y);
+	return scm_i_exact_rational_centered_divide (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_centered_divide, x, y, SCM_ARG2,
 			    s_scm_centered_divide);
@@ -3697,7 +3600,7 @@ SCM_PRIMITIVE_GENERIC (scm_centered_divide, "centered/", 2, 0, 0,
 	return scm_i_inexact_centered_divide
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_centered_divide (x, y);
+	return scm_i_exact_rational_centered_divide (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_centered_divide, x, y, SCM_ARG2,
 			    s_scm_centered_divide);
@@ -3708,7 +3611,7 @@ SCM_PRIMITIVE_GENERIC (scm_centered_divide, "centered/", 2, 0, 0,
 	  SCM_BIGP (y) || SCM_FRACTIONP (y))
 	return scm_i_inexact_centered_divide
 	  (SCM_REAL_VALUE (x), scm_to_double (y));
-     else
+      else
 	SCM_WTA_DISPATCH_2 (g_scm_centered_divide, x, y, SCM_ARG2,
 			    s_scm_centered_divide);
     }
@@ -3717,8 +3620,11 @@ SCM_PRIMITIVE_GENERIC (scm_centered_divide, "centered/", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_centered_divide
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_centered_divide (x, y);
       else
-	return scm_i_slow_exact_centered_divide (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_centered_divide, x, y, SCM_ARG2,
+			    s_scm_centered_divide);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_centered_divide, x, y, SCM_ARG1,
@@ -3796,35 +3702,24 @@ scm_i_bigint_centered_divide (SCM x, SCM y)
 				 scm_i_normbig (r)));
 }
 
-/* Compute exact centered quotient and remainder 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_centered_divide (SCM x, SCM y)
+scm_i_exact_rational_centered_divide (SCM x, SCM y)
 {
-  SCM q, r;
+  SCM q, r, r1;
+  SCM xd = scm_denominator (x);
+  SCM yd = scm_denominator (y);
 
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_centered_divide, x, y, SCM_ARG1,
-			s_scm_centered_divide);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_centered_divide, x, y, SCM_ARG2,
-			s_scm_centered_divide);
-  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_scm_centered_divide);
-  r = scm_difference (x, scm_product (q, y));
+  scm_i_extract_values_2
+    (scm_centered_divide (scm_product (scm_numerator (x), yd),
+			  scm_product (scm_numerator (y), xd)),
+     &q, &r1);
+  r = scm_divide (r1, scm_product (xd, yd));
   return scm_values (scm_list_2 (q, r));
 }
 
 static SCM scm_i_inexact_round_quotient (double x, double y);
 static SCM scm_i_bigint_round_quotient (SCM x, SCM y);
-static SCM scm_i_slow_exact_round_quotient (SCM x, SCM y);
+static SCM scm_i_exact_rational_round_quotient (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_round_quotient, "round-quotient", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -3893,7 +3788,7 @@ SCM_PRIMITIVE_GENERIC (scm_round_quotient, "round-quotient", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_round_quotient (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_round_quotient (x, y);
+	return scm_i_exact_rational_round_quotient (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_round_quotient, x, y, SCM_ARG2,
 			    s_scm_round_quotient);
@@ -3944,7 +3839,7 @@ SCM_PRIMITIVE_GENERIC (scm_round_quotient, "round-quotient", 2, 0, 0,
 	return scm_i_inexact_round_quotient
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_round_quotient (x, y);
+	return scm_i_exact_rational_round_quotient (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_round_quotient, x, y, SCM_ARG2,
 			    s_scm_round_quotient);
@@ -3964,8 +3859,11 @@ SCM_PRIMITIVE_GENERIC (scm_round_quotient, "round-quotient", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_round_quotient
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_round_quotient (x, y);
       else
-	return scm_i_slow_exact_round_quotient (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_round_quotient, x, y, SCM_ARG2,
+			    s_scm_round_quotient);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_round_quotient, x, y, SCM_ARG1,
@@ -4014,27 +3912,17 @@ scm_i_bigint_round_quotient (SCM x, SCM y)
   return scm_i_normbig (q);
 }
 
-/* Compute exact round quotient 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_round_quotient (SCM x, SCM y)
+scm_i_exact_rational_round_quotient (SCM x, SCM y)
 {
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_round_quotient, x, y, SCM_ARG1,
-			s_scm_round_quotient);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_round_quotient, x, y, SCM_ARG2,
-			s_scm_round_quotient);
-  else if (scm_is_true (scm_zero_p (y)))
-    scm_num_overflow (s_scm_round_quotient);
-  else
-    return scm_round_number (scm_divide (x, y));
+  return scm_round_quotient
+    (scm_product (scm_numerator (x), scm_denominator (y)),
+     scm_product (scm_numerator (y), scm_denominator (x)));
 }
 
 static SCM scm_i_inexact_round_remainder (double x, double y);
 static SCM scm_i_bigint_round_remainder (SCM x, SCM y);
-static SCM scm_i_slow_exact_round_remainder (SCM x, SCM y);
+static SCM scm_i_exact_rational_round_remainder (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_round_remainder, "round-remainder", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -4104,7 +3992,7 @@ SCM_PRIMITIVE_GENERIC (scm_round_remainder, "round-remainder", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_round_remainder (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_round_remainder (x, y);
+	return scm_i_exact_rational_round_remainder (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_round_remainder, x, y, SCM_ARG2,
 			    s_scm_round_remainder);
@@ -4152,7 +4040,7 @@ SCM_PRIMITIVE_GENERIC (scm_round_remainder, "round-remainder", 2, 0, 0,
 	return scm_i_inexact_round_remainder
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_round_remainder (x, y);
+	return scm_i_exact_rational_round_remainder (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_round_remainder, x, y, SCM_ARG2,
 			    s_scm_round_remainder);
@@ -4172,8 +4060,11 @@ SCM_PRIMITIVE_GENERIC (scm_round_remainder, "round-remainder", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_round_remainder
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_round_remainder (x, y);
       else
-	return scm_i_slow_exact_round_remainder (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_round_remainder, x, y, SCM_ARG2,
+			    s_scm_round_remainder);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_round_remainder, x, y, SCM_ARG1,
@@ -4234,31 +4125,20 @@ scm_i_bigint_round_remainder (SCM x, SCM y)
   return scm_i_normbig (r);
 }
 
-/* Compute exact round_remainder 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_round_remainder (SCM x, SCM y)
+scm_i_exact_rational_round_remainder (SCM x, SCM y)
 {
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_round_remainder, x, y, SCM_ARG1,
-			s_scm_round_remainder);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_round_remainder, x, y, SCM_ARG2,
-			s_scm_round_remainder);
-  else if (scm_is_true (scm_zero_p (y)))
-    scm_num_overflow (s_scm_round_remainder);
-  else
-    {
-      SCM q = scm_round_number (scm_divide (x, y));
-      return scm_difference (x, scm_product (q, y));
-    }
+  SCM xd = scm_denominator (x);
+  SCM yd = scm_denominator (y);
+  SCM r1 = scm_round_remainder (scm_product (scm_numerator (x), yd),
+				scm_product (scm_numerator (y), xd));
+  return scm_divide (r1, scm_product (xd, yd));
 }
 
 
 static SCM scm_i_inexact_round_divide (double x, double y);
 static SCM scm_i_bigint_round_divide (SCM x, SCM y);
-static SCM scm_i_slow_exact_round_divide (SCM x, SCM y);
+static SCM scm_i_exact_rational_round_divide (SCM x, SCM y);
 
 SCM_PRIMITIVE_GENERIC (scm_round_divide, "round/", 2, 0, 0,
 		       (SCM x, SCM y),
@@ -4332,7 +4212,7 @@ SCM_PRIMITIVE_GENERIC (scm_round_divide, "round/", 2, 0, 0,
       else if (SCM_REALP (y))
 	return scm_i_inexact_round_divide (xx, SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_round_divide (x, y);
+	return scm_i_exact_rational_round_divide (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_round_divide, x, y, SCM_ARG2,
 			    s_scm_round_divide);
@@ -4385,7 +4265,7 @@ SCM_PRIMITIVE_GENERIC (scm_round_divide, "round/", 2, 0, 0,
 	return scm_i_inexact_round_divide
 	  (scm_i_big2dbl (x), SCM_REAL_VALUE (y));
       else if (SCM_FRACTIONP (y))
-	return scm_i_slow_exact_round_divide (x, y);
+	return scm_i_exact_rational_round_divide (x, y);
       else
 	SCM_WTA_DISPATCH_2 (g_scm_round_divide, x, y, SCM_ARG2,
 			    s_scm_round_divide);
@@ -4396,7 +4276,7 @@ SCM_PRIMITIVE_GENERIC (scm_round_divide, "round/", 2, 0, 0,
 	  SCM_BIGP (y) || SCM_FRACTIONP (y))
 	return scm_i_inexact_round_divide
 	  (SCM_REAL_VALUE (x), scm_to_double (y));
-     else
+      else
 	SCM_WTA_DISPATCH_2 (g_scm_round_divide, x, y, SCM_ARG2,
 			    s_scm_round_divide);
     }
@@ -4405,8 +4285,11 @@ SCM_PRIMITIVE_GENERIC (scm_round_divide, "round/", 2, 0, 0,
       if (SCM_REALP (y))
 	return scm_i_inexact_round_divide
 	  (scm_i_fraction2double (x), SCM_REAL_VALUE (y));
+      else if (SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y))
+	return scm_i_exact_rational_round_divide (x, y);
       else
-	return scm_i_slow_exact_round_divide (x, y);
+	SCM_WTA_DISPATCH_2 (g_scm_round_divide, x, y, SCM_ARG2,
+			    s_scm_round_divide);
     }
   else
     SCM_WTA_DISPATCH_2 (g_scm_round_divide, x, y, SCM_ARG1,
@@ -4464,26 +4347,19 @@ scm_i_bigint_round_divide (SCM x, SCM y)
 				 scm_i_normbig (r)));
 }
 
-/* Compute exact round quotient and remainder 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_round_divide (SCM x, SCM y)
+scm_i_exact_rational_round_divide (SCM x, SCM y)
 {
-  if (!(SCM_I_INUMP (x) || SCM_BIGP (x) || SCM_FRACTIONP (x)))
-    SCM_WTA_DISPATCH_2 (g_scm_round_divide, x, y, SCM_ARG1,
-			s_scm_round_divide);
-  else if (!(SCM_I_INUMP (y) || SCM_BIGP (y) || SCM_FRACTIONP (y)))
-    SCM_WTA_DISPATCH_2 (g_scm_round_divide, x, y, SCM_ARG2,
-			s_scm_round_divide);
-  else if (scm_is_true (scm_zero_p (y)))
-    scm_num_overflow (s_scm_round_divide);
-  else
-    {
-      SCM q = scm_round_number (scm_divide (x, y));
-      SCM r = scm_difference (x, scm_product (q, y));
-      return scm_values (scm_list_2 (q, r));
-    }
+  SCM q, r, r1;
+  SCM xd = scm_denominator (x);
+  SCM yd = scm_denominator (y);
+
+  scm_i_extract_values_2
+    (scm_round_divide (scm_product (scm_numerator (x), yd),
+		       scm_product (scm_numerator (y), xd)),
+     &q, &r1);
+  r = scm_divide (r1, scm_product (xd, yd));
+  return scm_values (scm_list_2 (q, r));
 }
 
 
-- 
1.5.6.5


^ permalink raw reply related	[flat|nested] 12+ messages in thread
end of thread, other threads:[~2011-02-15 14:49 UTC | newest]

Thread overview: 12+ messages (download: mbox.gz follow: Atom feed
-- links below jump to the message on this page --
     [not found] <860989.4731.qm@web114104.mail.gq1.yahoo.com>
2011-02-12 20:19 ` [PATCH] New division operators, and optimization for fractions Andy Wingo
2011-02-13 14:55   ` Mark H Weaver
2011-02-14 19:44     ` Andy Wingo
     [not found] <595618.30923.qm@web114107.mail.gq1.yahoo.com>
2011-02-12 21:54 ` Andy Wingo
2011-02-12 22:38   ` Mark Weaver
2011-02-10 23:42 Mark H Weaver
2011-02-12 11:55 ` Andy Wingo
2011-02-14 13:56   ` Ken Raeburn
2011-02-14 18:17     ` Mark H Weaver
2011-02-15 10:38 ` Ludovic Courtès
2011-02-15 11:43   ` Mark H Weaver
2011-02-15 14:49     ` Ludovic Courtès

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).