X-Git-Url: https://git.distorted.org.uk/u/mdw/catacomb/blobdiff_plain/2b26f2d79b3c5342a0daaba359cc00e13de2a9ca..2685767a6125c1620719c7de6234aedf41857b7e:/mptext.c

diff --git a/mptext.c b/mptext.c
index f408b90..55ef579 100644
--- a/mptext.c
+++ b/mptext.c
@@ -1,6 +1,6 @@
 /* -*-c-*-
  *
- * $Id: mptext.c,v 1.5 2000/06/17 11:46:19 mdw Exp $
+ * $Id: mptext.c,v 1.9 2001/02/03 16:05:17 mdw Exp $
  *
  * Textual representation of multiprecision numbers
  *
@@ -30,6 +30,21 @@
 /*----- Revision history --------------------------------------------------* 
  *
  * $Log: mptext.c,v $
+ * Revision 1.9  2001/02/03 16:05:17  mdw
+ * Make flags be unsigned.  Improve the write algorithm: recurse until the
+ * parts are one word long and use single-precision arithmetic from there.
+ * Fix off-by-one bug when breaking the number apart.
+ *
+ * Revision 1.8  2000/12/06 20:32:42  mdw
+ * Reduce binary bytes (to allow marker bits to be ignored).  Fix error
+ * message string a bit.  Allow leading `+' signs.
+ *
+ * Revision 1.7  2000/07/15 10:01:08  mdw
+ * Bug fix in binary input.
+ *
+ * Revision 1.6  2000/06/25 12:58:23  mdw
+ * Fix the derivation of `depth' commentary.
+ *
  * Revision 1.5  2000/06/17 11:46:19  mdw
  * New and much faster stack-based algorithm for reading integers.  Support
  * reading and writing binary integers in bases between 2 and 256.
@@ -64,12 +79,13 @@
  *
  * This is the number of bits in a @size_t@ object.  Why? 
  *
- * Just to convince yourself that this is correct: let @b = MPW_MAX + 1@.
- * Then the largest possible @mp@ is %$M - 1$% where %$M = b^Z$%.  Let %$r$%
- * be a radix to read or write.  Since the recursion squares the radix at
- * each step, the highest number reached by the recursion is %$d$%, where:
+ * To see this, let %$b = \mathit{MPW\_MAX} + 1$% and let %$Z$% be the
+ * largest @size_t@ value.  Then the largest possible @mp@ is %$M - 1$% where
+ * %$M = b^Z$%.  Let %$r$% be a radix to read or write.  Since the recursion
+ * squares the radix at each step, the highest number reached by the
+ * recursion is %$d$%, where:
  *
- *   %$r^(2^d) = b^Z$%.
+ *   %$r^{2^d} = b^Z$%.
  *
  * Solving gives that %$d = \lg \log_r b^Z$%.  If %$r = 2$%, this is maximum,
  * so choosing %$d = \lg \lg b^Z = \lg (Z \lg b) = \lg Z + \lg \lg b$%.
@@ -136,10 +152,8 @@ mp *mp_read(mp *m, int radix, const mptext_ops *ops, void *p)
 
   /* --- Flags --- */
 
-  enum {
-    f_neg = 1u,
-    f_ok = 2u
-  };
+#define f_neg 1u
+#define f_ok 2u
 
   /* --- Initialize the stacks --- */
 
@@ -162,11 +176,10 @@ mp *mp_read(mp *m, int radix, const mptext_ops *ops, void *p)
 
   /* --- Handle an initial sign --- */
 
-  if (ch == '-') {
-    f |= f_neg;
-    ch = ops->get(p);
-    while (isspace(ch))
-      ch = ops->get(p);
+  if (radix >= 0 && (ch == '-' || ch == '+')) {
+    if (ch == '-')
+      f |= f_neg;
+    do ch = ops->get(p); while isspace(ch);
   }
 
   /* --- If the radix is zero, look for leading zeros --- */
@@ -177,7 +190,7 @@ mp *mp_read(mp *m, int radix, const mptext_ops *ops, void *p)
     r = -1;
   } else if (radix < 0) {
     rd = -radix;
-    assert(((void)"binary radix must fit in a byte ", rd < UCHAR_MAX));
+    assert(((void)"binary radix must fit in a byte", rd < UCHAR_MAX));
     r = -1;
   } else if (ch != '0') {
     rd = 10;
@@ -199,6 +212,9 @@ mp *mp_read(mp *m, int radix, const mptext_ops *ops, void *p)
   for (;; ch = ops->get(p)) {
     int x;
 
+    if (ch < 0)
+      break;
+
     /* --- An underscore indicates a numbered base --- */
 
     if (ch == '_' && r > 0 && r <= 36) {
@@ -224,7 +240,7 @@ mp *mp_read(mp *m, int radix, const mptext_ops *ops, void *p)
     /* --- Check that the character is a digit and in range --- */
 
     if (radix < 0)
-      x = ch;
+      x = ch % rd;
     else {
       if (!isalnum(ch))
 	break;
@@ -340,6 +356,9 @@ mp *mp_read(mp *m, int radix, const mptext_ops *ops, void *p)
   if (f & f_neg)
     m->f |= MP_NEG;
   return (m);
+
+#undef f_neg
+#undef f_ok
 }
 
 /* --- @mp_write@ --- *
@@ -356,14 +375,14 @@ mp *mp_read(mp *m, int radix, const mptext_ops *ops, void *p)
 
 /* --- Simple case --- *
  *
- * Use a fixed-sized buffer and the simple single-precision division
- * algorithm to pick off low-order digits.  Put each digit in a buffer,
- * working backwards from the end.  If the buffer becomes full, recurse to
- * get another one.  Ensure that there are at least @z@ digits by writing
- * leading zeroes if there aren't enough real digits.
+ * Use a fixed-sized buffer and single-precision arithmetic to pick off
+ * low-order digits.  Put each digit in a buffer, working backwards from the
+ * end.  If the buffer becomes full, recurse to get another one.  Ensure that
+ * there are at least @z@ digits by writing leading zeroes if there aren't
+ * enough real digits.
  */
 
-static int simple(mp *m, int radix, unsigned z,
+static int simple(mpw n, int radix, unsigned z,
 		  const mptext_ops *ops, void *p)
 {
   int rc = 0;
@@ -375,23 +394,21 @@ static int simple(mp *m, int radix, unsigned z,
     int ch;
     mpw x;
 
-    x = mpx_udivn(m->v, m->vl, m->v, m->vl, rd);
-    MP_SHRINK(m);
+    x = n % rd;
+    n /= rd;
     if (radix < 0)
       ch = x;
-    else {
-      if (x < 10)
-	ch = '0' + x;
-      else
-	ch = 'a' + x - 10;
-    }
+    else if (x < 10)
+      ch = '0' + x;
+    else
+      ch = 'a' + x - 10;
     buf[--i] = ch;
     if (z)
       z--;
-  } while (i && MP_LEN(m));
+  } while (i && n);
 
-  if (MP_LEN(m))
-    rc = simple(m, radix, z, ops, p);
+  if (n)
+    rc = simple(n, radix, z, ops, p);
   else {
     static const char zero[32] = "00000000000000000000000000000000";
     while (!rc && z >= sizeof(zero)) {
@@ -402,9 +419,8 @@ static int simple(mp *m, int radix, unsigned z,
       rc = ops->put(zero, z, p);
   }
   if (!rc)
-    ops->put(buf + i, sizeof(buf) - i, p);
-  if (m->f & MP_BURN)
-    BURN(buf);
+    rc = ops->put(buf + i, sizeof(buf) - i, p);
+  BURN(buf);
   return (rc);
 }
 
@@ -423,9 +439,10 @@ static int complicated(mp *m, int radix, mp **pr, unsigned i, unsigned z,
   mp *q = MP_NEW;
   unsigned d = 1 << i;
 
-  if (MP_LEN(m) < 8)
-    return (simple(m, radix, z, ops, p));
+  if (MP_LEN(m) < 2)
+    return (simple(MP_LEN(m) ? m->v[0] : 0, radix, z, ops, p));
 
+  assert(i);
   mp_div(&q, &m, m, pr[i]);
   if (!MP_LEN(q))
     d = z;
@@ -472,8 +489,8 @@ int mp_write(mp *m, int radix, const mptext_ops *ops, void *p)
 
   /* --- If the number is small, do it the easy way --- */
 
-  if (MP_LEN(m) < 8)
-    rc = simple(m, radix, 0, ops, p);
+  if (MP_LEN(m) < 2)
+    rc = simple(MP_LEN(m) ? m->v[0] : 0, radix, 0, ops, p);
 
   /* --- Use a clever algorithm --- *
    *
@@ -486,7 +503,7 @@ int mp_write(mp *m, int radix, const mptext_ops *ops, void *p)
 
   else {
     mp *pr[DEPTH];
-    size_t target = MP_LEN(m) / 2;
+    size_t target = (MP_LEN(m) + 1) / 2;
     unsigned i = 0;
     mp *z = mp_new(1, 0);