/* -*-c-*-
*
- * $Id: mptext.c,v 1.1 1999/11/17 18:02:16 mdw Exp $
+ * $Id: mptext.c,v 1.4 1999/12/22 15:56:56 mdw Exp $
*
* Textual representation of multiprecision numbers
*
/*----- Revision history --------------------------------------------------*
*
* $Log: mptext.c,v $
+ * Revision 1.4 1999/12/22 15:56:56 mdw
+ * Use clever recursive algorithm for writing numbers out.
+ *
+ * Revision 1.3 1999/12/10 23:23:26 mdw
+ * Allocate slightly less memory.
+ *
+ * Revision 1.2 1999/11/20 22:24:15 mdw
+ * Use function versions of MPX_UMULN and MPX_UADDN.
+ *
* Revision 1.1 1999/11/17 18:02:16 mdw
* New multiprecision integer arithmetic suite.
*
#include <ctype.h>
#include <stdio.h>
-#include <mLib/darray.h>
-
#include "mp.h"
#include "mptext.h"
-
-/*----- Data structures ---------------------------------------------------*/
-
-#ifndef CHARV
- DA_DECL(charv, char);
-# define CHARV
-#endif
+#include "paranoia.h"
/*----- Main code ---------------------------------------------------------*/
/* --- Initialize the destination number --- */
- MP_MODIFY(m, 16);
+ MP_MODIFY(m, 4);
m->vl = m->v;
m->f &= ~MP_UNDEF;
/* --- Stick the character on the end of my integer --- */
- MP_ENSURE(m, MP_LEN(m) + 1);
- MPX_UMULN(m->v, m->vl, m->v, m->vl - 1, radix);
- MPX_UADDN(m->v, m->vl, x);
- MP_SHRINK(m);
+ mp_ensure(m, MP_LEN(m) + 1);
+ mpx_umuln(m->v, m->vl, m->v, m->vl - 1, radix);
+ mpx_uaddn(m->v, m->vl, x);
+ mp_shrink(m);
f |= f_ok;
}
* Use: Writes a large integer in textual form.
*/
+/* --- 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.
+ */
+
+static int simple(mp *m, int radix, unsigned z,
+ const mptext_ops *ops, void *p)
+{
+ int rc = 0;
+ char buf[64];
+ unsigned i = sizeof(buf);
+
+ do {
+ int ch;
+ mpw x;
+
+ x = mpx_udivn(m->v, m->vl, m->v, m->vl, radix);
+ MP_SHRINK(m);
+ if (x < 10)
+ ch = '0' + x;
+ else
+ ch = 'a' + x - 10;
+ buf[--i] = ch;
+ if (z)
+ z--;
+ } while (i && MP_LEN(m));
+
+ if (MP_LEN(m))
+ rc = simple(m, radix, z, ops, p);
+ else {
+ static const char zero[32] = "00000000000000000000000000000000";
+ while (!rc && z >= sizeof(zero)) {
+ rc = ops->put(zero, sizeof(zero), p);
+ z -= sizeof(zero);
+ }
+ if (!rc && z)
+ rc = ops->put(zero, z, p);
+ }
+ if (!rc)
+ ops->put(buf + i, sizeof(buf) - i, p);
+ if (m->f & MP_BURN)
+ BURN(buf);
+ return (rc);
+}
+
+/* --- Complicated case --- *
+ *
+ * If the number is small, fall back to the simple case above. Otherwise
+ * divide and take remainder by current large power of the radix, and emit
+ * each separately. Don't emit a zero quotient. Be very careful about
+ * leading zeroes on the remainder part, because they're deeply significant.
+ */
+
+static int complicated(mp *m, int radix, mp **pr, unsigned i, unsigned z,
+ const mptext_ops *ops, void *p)
+{
+ int rc = 0;
+ mp *q = MP_NEW;
+ unsigned d = 1 << i;
+
+ if (MP_LEN(m) < 8)
+ return (simple(m, radix, z, ops, p));
+
+ mp_div(&q, &m, m, pr[i]);
+ if (!MP_LEN(q))
+ d = z;
+ else {
+ if (z > d)
+ z -= d;
+ else
+ z = 0;
+ rc = complicated(q, radix, pr, i - 1, z, ops, p);
+ }
+ if (!rc)
+ rc = complicated(m, radix, pr, i - 1, d, ops, p);
+ mp_drop(q);
+ return (rc);
+}
+
+/* --- Main driver code --- */
+
int mp_write(mp *m, int radix, const mptext_ops *ops, void *p)
{
- charv v = DA_INIT;
- mpw b = radix;
- mp bb;
+ int rc;
/* --- Set various things up --- */
m = MP_COPY(m);
- mp_build(&bb, &b, &b + 1);
+ MP_SPLIT(m);
/* --- If the number is negative, sort that out --- */
if (m->f & MP_NEG) {
if (ops->put("-", 1, p))
return (EOF);
- MP_SPLIT(m);
- m->f &= ~MP_NEG;
}
- /* --- Write digits to a temporary array --- */
-
- do {
- mp *q = MP_NEW, *r = MP_NEW;
- int ch;
- mpw x;
+ /* --- If the number is small, do it the easy way --- */
+
+ if (MP_LEN(m) < 8)
+ rc = simple(m, radix, 0, ops, p);
+
+ /* --- Use a clever algorithm --- *
+ *
+ * Square the radix repeatedly, remembering old results, until I get
+ * something more than half the size of the number @m@. Use this to divide
+ * the number: the quotient and remainder will be approximately the same
+ * size, and I'll have split them on a digit boundary, so I can just emit
+ * the quotient and remainder recursively, in order.
+ *
+ * The array size copes with the largest number possibly representable on
+ * the host machine. Such a large number shouldn't ever arise in real use.
+ */
+
+ else {
+ mp *pr[CHAR_BIT * sizeof(size_t)];
+ size_t target = MP_LEN(m) / 2;
+ unsigned i = 0;
+ mp *z = mp_create(1);
+
+ /* --- Set up the exponent table --- */
+
+ z->v[0] = radix;
+ z->f = 0;
+ for (;;) {
+ assert(((void)"Number is too unimaginably huge",
+ i < sizeof(pr) / sizeof(pr[0])));
+ pr[i++] = z;
+ if (MP_LEN(z) > target)
+ break;
+ z = mp_sqr(MP_NEW, z);
+ }
- mp_div(&q, &r, m, &bb);
- x = r->v[0];
- if (x < 10)
- ch = '0' + x;
- else
- ch = 'a' + x - 10;
- DA_UNSHIFT(&v, ch);
- MP_DROP(m);
- MP_DROP(r);
- m = q;
- } while (MP_CMP(m, !=, MP_ZERO));
+ /* --- Write out the answer --- */
- /* --- Finished that --- */
+ rc = complicated(m, radix, pr, i - 1, 0, ops, p);
- MP_DROP(m);
+ /* --- Tidy away the array --- */
- {
- int rc = ops->put(DA(&v), DA_LEN(&v), p);
- DA_DESTROY(&v);
- return (rc ? EOF : 0);
+ while (i > 0)
+ mp_drop(pr[--i]);
}
+
+ /* --- Tidying up code --- */
+
+ MP_DROP(m);
+ return (rc);
}
/*----- Test rig ----------------------------------------------------------*/
}
dstr_destroy(&d);
+ assert(mparena_count(MPARENA_GLOBAL) == 0);
return (ok);
}