+++ /dev/null
-/* -*-c-*-
- *
- * $Id$
- *
- * Textual representation of multiprecision numbers
- *
- * (c) 1999 Straylight/Edgeware
- */
-
-/*----- Licensing notice --------------------------------------------------*
- *
- * This file is part of Catacomb.
- *
- * Catacomb is free software; you can redistribute it and/or modify
- * it under the terms of the GNU Library General Public License as
- * published by the Free Software Foundation; either version 2 of the
- * License, or (at your option) any later version.
- *
- * Catacomb is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- * GNU Library General Public License for more details.
- *
- * You should have received a copy of the GNU Library General Public
- * License along with Catacomb; if not, write to the Free
- * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
- * MA 02111-1307, USA.
- */
-
-/*----- Header files ------------------------------------------------------*/
-
-#include <ctype.h>
-#include <limits.h>
-#include <stdio.h>
-
-#include "mp.h"
-#include "mptext.h"
-#include "paranoia.h"
-
-/*----- Magical numbers ---------------------------------------------------*/
-
-/* --- Maximum recursion depth --- *
- *
- * This is the number of bits in a @size_t@ object. Why?
- *
- * To see this, let %$b = \textit{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$%.
- *
- * 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$%.
- *
- * Expressing %$\lg Z$% as @CHAR_BIT * sizeof(size_t)@ yields an
- * overestimate, since a @size_t@ representation may contain `holes'.
- * Choosing to represent %$\lg \lg b$% by 10 is almost certainly sufficient
- * for `some time to come'.
- */
-
-#define DEPTH (CHAR_BIT * sizeof(size_t) + 10)
-
-/*----- Main code ---------------------------------------------------------*/
-
-/* --- @mp_read@ --- *
- *
- * Arguments: @mp *m@ = destination multiprecision number
- * @int radix@ = base to assume for data (or zero to guess)
- * @const mptext_ops *ops@ = pointer to operations block
- * @void *p@ = data for the operations block
- *
- * Returns: The integer read, or zero if it didn't work.
- *
- * Use: Reads an integer from some source. If the @radix@ is
- * specified, the number is assumed to be given in that radix,
- * with the letters `a' (either upper- or lower-case) upwards
- * standing for digits greater than 9. Otherwise, base 10 is
- * assumed unless the number starts with `0' (octal), `0x' (hex)
- * or `nnn_' (base `nnn'). An arbitrary amount of whitespace
- * before the number is ignored.
- */
-
-/* --- About the algorithm --- *
- *
- * The algorithm here is rather aggressive. I maintain an array of
- * successive squarings of the radix, and a stack of partial results, each
- * with a counter attached indicating which radix square to multiply by.
- * Once the item at the top of the stack reaches the same counter level as
- * the next item down, they are combined together and the result is given a
- * counter level one higher than either of the results.
- *
- * Gluing the results together at the end is slightly tricky. Pay attention
- * to the code.
- *
- * This is more complicated because of the need to handle the slightly
- * bizarre syntax.
- */
-
-mp *mp_read(mp *m, int radix, const mptext_ops *ops, void *p)
-{
- int ch; /* Current char being considered */
- unsigned f = 0; /* Flags about the current number */
- int r; /* Radix to switch over to */
- mpw rd; /* Radix as an @mp@ digit */
- mp rr; /* The @mp@ for the radix */
- unsigned nf = m ? m->f & MP_BURN : 0; /* New @mp@ flags */
-
- /* --- Stacks --- */
-
- mp *pow[DEPTH]; /* List of powers */
- unsigned pows; /* Next index to fill */
- struct { unsigned i; mp *m; } s[DEPTH]; /* Main stack */
- unsigned sp; /* Current stack pointer */
-
- /* --- Flags --- */
-
-#define f_neg 1u
-#define f_ok 2u
-#define f_start 4u
-
- /* --- Initialize the stacks --- */
-
- mp_build(&rr, &rd, &rd + 1);
- pow[0] = &rr;
- pows = 1;
-
- sp = 0;
-
- /* --- Initialize the destination number --- */
-
- if (m)
- MP_DROP(m);
-
- /* --- Read an initial character --- */
-
- ch = ops->get(p);
- if (radix >= 0) {
- while (isspace(ch))
- ch = ops->get(p);
- }
-
- /* --- Handle an initial sign --- */
-
- 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 --- */
-
- if (radix > 0) {
- assert(((void)"ascii radix must be <= 62", radix <= 62));
- rd = radix;
- r = -1;
- } else if (radix < 0) {
- rd = -radix;
- assert(((void)"binary radix must fit in a byte", rd <= UCHAR_MAX));
- r = -1;
- } else if (ch != '0') {
- rd = 10;
- r = 0;
- } else {
- ch = ops->get(p);
- switch (ch) {
- case 'x':
- rd = 16;
- goto prefix;
- case 'o':
- rd = 8;
- goto prefix;
- case 'b':
- rd = 2;
- goto prefix;
- prefix:
- ch = ops->get(p);
- break;
- default:
- rd = 8;
- f |= f_ok;
- }
- r = -1;
- }
-
- /* --- Use fast algorithm for binary radix --- *
- *
- * This is the restart point after having parsed a radix number from the
- * input. We check whether the radix is binary, and if so use a fast
- * algorithm which just stacks the bits up in the right order.
- */
-
-restart:
- switch (rd) {
- unsigned bit;
-
- case 2: bit = 1; goto bin;
- case 4: bit = 2; goto bin;
- case 8: bit = 3; goto bin;
- case 16: bit = 4; goto bin;
- case 32: bit = 5; goto bin;
- case 64: bit = 6; goto bin;
- case 128: bit = 7; goto bin;
- default:
- break;
-
- /* --- The fast binary algorithm --- *
- *
- * We stack bits up starting at the top end of a word. When one word is
- * full, we write it to the integer, and start another with the left-over
- * bits. When the array in the integer is full, we resize using low-level
- * calls and copy the current data to the top end. Finally, we do a single
- * bit-shift when we know where the end of the number is.
- */
-
- bin: {
- mpw a = 0;
- unsigned b = MPW_BITS;
- size_t len, n;
- mpw *v;
-
- m = mp_dest(MP_NEW, 1, nf);
- len = n = m->sz;
- n = len;
- v = m->v + n;
- for (;; ch = ops->get(p)) {
- unsigned x;
-
- if (ch < 0)
- break;
-
- /* --- Check that the character is a digit and in range --- */
-
- if (radix < 0)
- x = ch % rd;
- else {
- if (!isalnum(ch))
- break;
- if (ch >= '0' && ch <= '9')
- x = ch - '0';
- else {
- if (rd <= 36)
- ch = tolower(ch);
- if (ch >= 'a' && ch <= 'z') /* ASCII dependent! */
- x = ch - 'a' + 10;
- else if (ch >= 'A' && ch <= 'Z')
- x = ch - 'A' + 36;
- else
- break;
- }
- }
- if (x >= rd)
- break;
-
- /* --- Feed the digit into the accumulator --- */
-
- f |= f_ok;
- if (!x && !(f & f_start))
- continue;
- f |= f_start;
- if (b > bit) {
- b -= bit;
- a |= MPW(x) << b;
- } else {
- a |= MPW(x) >> (bit - b);
- b += MPW_BITS - bit;
- *--v = MPW(a);
- n--;
- if (!n) {
- n = len;
- len <<= 1;
- v = mpalloc(m->a, len);
- memcpy(v + n, m->v, MPWS(n));
- mpfree(m->a, m->v);
- m->v = v;
- v = m->v + n;
- }
- a = (b < MPW_BITS) ? MPW(x) << b : 0;
- }
- }
-
- /* --- Finish up --- */
-
- if (!(f & f_ok)) {
- mp_drop(m);
- m = 0;
- } else {
- *--v = MPW(a);
- n--;
- m->sz = len;
- m->vl = m->v + len;
- m->f &= ~MP_UNDEF;
- m = mp_lsr(m, m, (unsigned long)n * MPW_BITS + b);
- }
- ops->unget(ch, p);
- goto done;
- }}
-
- /* --- Time to start --- */
-
- for (;; ch = ops->get(p)) {
- unsigned x;
-
- if (ch < 0)
- break;
-
- /* --- An underscore indicates a numbered base --- */
-
- if (ch == '_' && r > 0 && r <= 62) {
- unsigned i;
-
- /* --- Clear out the stacks --- */
-
- for (i = 1; i < pows; i++)
- MP_DROP(pow[i]);
- pows = 1;
- for (i = 0; i < sp; i++)
- MP_DROP(s[i].m);
- sp = 0;
-
- /* --- Restart the search --- */
-
- rd = r;
- r = -1;
- f &= ~f_ok;
- ch = ops->get(p);
- goto restart;
- }
-
- /* --- Check that the character is a digit and in range --- */
-
- if (radix < 0)
- x = ch % rd;
- else {
- if (!isalnum(ch))
- break;
- if (ch >= '0' && ch <= '9')
- x = ch - '0';
- else {
- if (rd <= 36)
- ch = tolower(ch);
- if (ch >= 'a' && ch <= 'z') /* ASCII dependent! */
- x = ch - 'a' + 10;
- else if (ch >= 'A' && ch <= 'Z')
- x = ch - 'A' + 36;
- else
- break;
- }
- }
-
- /* --- Sort out what to do with the character --- */
-
- if (x >= 10 && r >= 0)
- r = -1;
- if (x >= rd)
- break;
-
- if (r >= 0)
- r = r * 10 + x;
-
- /* --- Stick the character on the end of my integer --- */
-
- assert(((void)"Number is too unimaginably huge", sp < DEPTH));
- s[sp].m = m = mp_new(1, nf);
- m->v[0] = x;
- s[sp].i = 0;
-
- /* --- Now grind through the stack --- */
-
- while (sp > 0 && s[sp - 1].i == s[sp].i) {
-
- /* --- Combine the top two items --- */
-
- sp--;
- m = s[sp].m;
- m = mp_mul(m, m, pow[s[sp].i]);
- m = mp_add(m, m, s[sp + 1].m);
- s[sp].m = m;
- MP_DROP(s[sp + 1].m);
- s[sp].i++;
-
- /* --- Make a new radix power if necessary --- */
-
- if (s[sp].i >= pows) {
- assert(((void)"Number is too unimaginably huge", pows < DEPTH));
- pow[pows] = mp_sqr(MP_NEW, pow[pows - 1]);
- pows++;
- }
- }
- f |= f_ok;
- sp++;
- }
-
- ops->unget(ch, p);
-
- /* --- If we're done, compute the rest of the number --- */
-
- if (f & f_ok) {
- if (!sp)
- return (MP_ZERO);
- else {
- mp *z = MP_ONE;
- sp--;
-
- while (sp > 0) {
-
- /* --- Combine the top two items --- */
-
- sp--;
- m = s[sp].m;
- z = mp_mul(z, z, pow[s[sp + 1].i]);
- m = mp_mul(m, m, z);
- m = mp_add(m, m, s[sp + 1].m);
- s[sp].m = m;
- MP_DROP(s[sp + 1].m);
-
- /* --- Make a new radix power if necessary --- */
-
- if (s[sp].i >= pows) {
- assert(((void)"Number is too unimaginably huge", pows < DEPTH));
- pow[pows] = mp_sqr(MP_NEW, pow[pows - 1]);
- pows++;
- }
- }
- MP_DROP(z);
- m = s[0].m;
- }
- } else {
- unsigned i;
- for (i = 0; i < sp; i++)
- MP_DROP(s[i].m);
- }
-
- /* --- Clear the radix power list --- */
-
- {
- unsigned i;
- for (i = 1; i < pows; i++)
- MP_DROP(pow[i]);
- }
-
- /* --- Bail out if the number was bad --- */
-
-done:
- if (!(f & f_ok))
- return (0);
-
- /* --- Set the sign and return --- */
-
- if (f & f_neg)
- m->f |= MP_NEG;
- MP_SHRINK(m);
- return (m);
-
-#undef f_start
-#undef f_neg
-#undef f_ok
-}
-
-/* --- @mp_write@ --- *
- *
- * Arguments: @mp *m@ = pointer to a multi-precision integer
- * @int radix@ = radix to use when writing the number out
- * @const mptext_ops *ops@ = pointer to an operations block
- * @void *p@ = data for the operations block
- *
- * Returns: Zero if it worked, nonzero otherwise.
- *
- * Use: Writes a large integer in textual form.
- */
-
-/* --- Simple case --- *
- *
- * 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(mpw n, int radix, unsigned z,
- const mptext_ops *ops, void *p)
-{
- int rc = 0;
- char buf[64];
- unsigned i = sizeof(buf);
- int rd = radix > 0 ? radix : -radix;
-
- do {
- int ch;
- mpw x;
-
- x = n % rd;
- n /= rd;
- if (radix < 0)
- ch = x;
- else if (x < 10)
- ch = '0' + x;
- else if (x < 36) /* Ascii specific */
- ch = 'a' + x - 10;
- else
- ch = 'A' + x - 36;
- buf[--i] = ch;
- if (z)
- z--;
- } while (i && n);
-
- if (n)
- rc = simple(n, radix, z, ops, p);
- else {
- char zbuf[32];
- memset(zbuf, (radix < 0) ? 0 : '0', sizeof(zbuf));
- while (!rc && z >= sizeof(zbuf)) {
- rc = ops->put(zbuf, sizeof(zbuf), p);
- z -= sizeof(zbuf);
- }
- if (!rc && z)
- rc = ops->put(zbuf, z, p);
- }
- if (!rc)
- rc = ops->put(buf + i, sizeof(buf) - i, p);
- 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) < 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_ZEROP(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);
-}
-
-/* --- Binary case --- *
- *
- * Special case for binary output. Goes much faster.
- */
-
-static int binary(mp *m, int bit, int radix, const mptext_ops *ops, void *p)
-{
- mpw *v;
- mpw a;
- int rc = 0;
- unsigned b;
- unsigned mask;
- unsigned long n;
- unsigned f = 0;
- char buf[8], *q;
- unsigned x;
- int ch;
-
-#define f_out 1u
-
- /* --- Work out where to start --- */
-
- n = mp_bits(m);
- if (n % bit)
- n += bit - (n % bit);
- b = n % MPW_BITS;
- n /= MPW_BITS;
-
- if (n >= MP_LEN(m)) {
- n--;
- b += MPW_BITS;
- }
-
- v = m->v + n;
- a = *v;
- mask = (1 << bit) - 1;
- q = buf;
-
- /* --- Main code --- */
-
- for (;;) {
- if (b > bit) {
- b -= bit;
- x = a >> b;
- } else {
- x = a << (bit - b);
- b += MPW_BITS - bit;
- if (v == m->v)
- break;
- a = *--v;
- if (b < MPW_BITS)
- x |= a >> b;
- }
- x &= mask;
- if (!x && !(f & f_out))
- continue;
-
- if (radix < 0)
- ch = x;
- else if (x < 10)
- ch = '0' + x;
- else if (x < 36)
- ch = 'a' + x - 10; /* Ascii specific */
- else
- ch = 'A' + x - 36;
- *q++ = ch;
- if (q >= buf + sizeof(buf)) {
- if ((rc = ops->put(buf, sizeof(buf), p)) != 0)
- goto done;
- q = buf;
- }
- f |= f_out;
- }
-
- x &= mask;
- if (radix < 0)
- ch = x;
- else if (x < 10)
- ch = '0' + x;
- else if (x < 36)
- ch = 'a' + x - 10; /* Ascii specific */
- else
- ch = 'A' + x - 36;
- *q++ = ch;
- rc = ops->put(buf, q - buf, p);
-
-done:
- mp_drop(m);
- return (rc);
-
-#undef f_out
-}
-
-/* --- Main driver code --- */
-
-int mp_write(mp *m, int radix, const mptext_ops *ops, void *p)
-{
- int rc;
-
- if (MP_EQ(m, MP_ZERO))
- return (ops->put(radix > 0 ? "0" : "\0", 1, p));
-
- /* --- Set various things up --- */
-
- m = MP_COPY(m);
- MP_SPLIT(m);
-
- /* --- Check the radix for sensibleness --- */
-
- if (radix > 0)
- assert(((void)"ascii radix must be <= 62", radix <= 62));
- else if (radix < 0)
- assert(((void)"binary radix must fit in a byte", -radix <= UCHAR_MAX));
- else
- assert(((void)"radix can't be zero in mp_write", 0));
-
- /* --- If the number is negative, sort that out --- */
-
- if (MP_NEGP(m)) {
- assert(radix > 0);
- if (ops->put("-", 1, p))
- return (EOF);
- m->f &= ~MP_NEG;
- }
-
- /* --- Handle binary radix --- */
-
- switch (radix) {
- case 2: case -2: return (binary(m, 1, radix, ops, p));
- case 4: case -4: return (binary(m, 2, radix, ops, p));
- case 8: case -8: return (binary(m, 3, radix, ops, p));
- case 16: case -16: return (binary(m, 4, radix, ops, p));
- case 32: case -32: return (binary(m, 5, radix, ops, p));
- case -64: return (binary(m, 6, radix, ops, p));
- case -128: return (binary(m, 7, radix, ops, p));
- }
-
- /* --- If the number is small, do it the easy way --- */
-
- if (MP_LEN(m) < 2)
- rc = simple(MP_LEN(m) ? m->v[0] : 0, 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.
- */
-
- else {
- mp *pr[DEPTH];
- size_t target = (MP_LEN(m) + 1) / 2;
- unsigned i = 0;
- mp *z = mp_new(1, 0);
-
- /* --- Set up the exponent table --- */
-
- z->v[0] = (radix > 0 ? radix : -radix);
- z->f = 0;
- for (;;) {
- assert(((void)"Number is too unimaginably huge", i < DEPTH));
- pr[i++] = z;
- if (MP_LEN(z) > target)
- break;
- z = mp_sqr(MP_NEW, z);
- }
-
- /* --- Write out the answer --- */
-
- rc = complicated(m, radix, pr, i - 1, 0, ops, p);
-
- /* --- Tidy away the array --- */
-
- while (i > 0)
- mp_drop(pr[--i]);
- }
-
- /* --- Tidying up code --- */
-
- MP_DROP(m);
- return (rc);
-}
-
-/*----- Test rig ----------------------------------------------------------*/
-
-#ifdef TEST_RIG
-
-#include <mLib/testrig.h>
-
-static int verify(dstr *v)
-{
- int ok = 1;
- int ib = *(int *)v[0].buf, ob = *(int *)v[2].buf;
- dstr d = DSTR_INIT;
- size_t off = 0;
- mp *m = mp_readdstr(MP_NEW, &v[1], &off, ib);
- if (m) {
- if (!ob) {
- fprintf(stderr, "*** unexpected successful parse\n"
- "*** input [%2i] = ", ib);
- if (ib < 0)
- type_hex.dump(&v[1], stderr);
- else
- fputs(v[1].buf, stderr);
- mp_writedstr(m, &d, 10);
- fprintf(stderr, "\n*** (value = %s)\n", d.buf);
- ok = 0;
- } else {
- mp_writedstr(m, &d, ob);
- if (d.len != v[3].len || memcmp(d.buf, v[3].buf, d.len) != 0) {
- fprintf(stderr, "*** failed read or write\n"
- "*** input [%2i] = ", ib);
- if (ib < 0)
- type_hex.dump(&v[1], stderr);
- else
- fputs(v[1].buf, stderr);
- fprintf(stderr, "\n*** output [%2i] = ", ob);
- if (ob < 0)
- type_hex.dump(&d, stderr);
- else
- fputs(d.buf, stderr);
- fprintf(stderr, "\n*** expected [%2i] = ", ob);
- if (ob < 0)
- type_hex.dump(&v[3], stderr);
- else
- fputs(v[3].buf, stderr);
- fputc('\n', stderr);
- ok = 0;
- }
- }
- mp_drop(m);
- } else {
- if (ob) {
- fprintf(stderr, "*** unexpected parse failure\n"
- "*** input [%2i] = ", ib);
- if (ib < 0)
- type_hex.dump(&v[1], stderr);
- else
- fputs(v[1].buf, stderr);
- fprintf(stderr, "\n*** expected [%2i] = ", ob);
- if (ob < 0)
- type_hex.dump(&v[3], stderr);
- else
- fputs(v[3].buf, stderr);
- fputc('\n', stderr);
- ok = 0;
- }
- }
-
- if (v[1].len - off != v[4].len ||
- memcmp(v[1].buf + off, v[4].buf, v[4].len) != 0) {
- fprintf(stderr, "*** leftovers incorrect\n"
- "*** input [%2i] = ", ib);
- if (ib < 0)
- type_hex.dump(&v[1], stderr);
- else
- fputs(v[1].buf, stderr);
- fprintf(stderr, "\n*** expected `%s'\n"
- "*** found `%s'\n",
- v[4].buf, v[1].buf + off);
- ok = 0;
- }
-
- dstr_destroy(&d);
- assert(mparena_count(MPARENA_GLOBAL) == 0);
- return (ok);
-}
-
-static test_chunk tests[] = {
- { "mptext-ascii", verify,
- { &type_int, &type_string, &type_int, &type_string, &type_string, 0 } },
- { "mptext-bin-in", verify,
- { &type_int, &type_hex, &type_int, &type_string, &type_string, 0 } },
- { "mptext-bin-out", verify,
- { &type_int, &type_string, &type_int, &type_hex, &type_string, 0 } },
- { 0, 0, { 0 } }
-};
-
-int main(int argc, char *argv[])
-{
- sub_init();
- test_run(argc, argv, tests, SRCDIR "/tests/mptext");
- return (0);
-}
-
-#endif
-
-/*----- That's all, folks -------------------------------------------------*/
projects / u / mdw / catacomb / blobdiff