summaryrefslogtreecommitdiff
path: root/big.c
diff options
context:
space:
mode:
Diffstat (limited to 'big.c')
-rw-r--r--big.c914
1 files changed, 379 insertions, 535 deletions
diff --git a/big.c b/big.c
index 705610a..7b58a2f 100644
--- a/big.c
+++ b/big.c
@@ -14,13 +14,6 @@
#include "big.h"
#include <errno.h>
-/* contains remainder from last call to longdiv or shortdiv */
-static word b_rem;
-
-extern char *dicp;
-
-#define maxval (1ll<<60)
-
static double logIBASE, log10IBASE;
word big_one;
@@ -34,94 +27,89 @@ static word shift(word, word);
static word shortdiv(word, word);
static word stimes(word, word);
-void bigsetup(void)
+void bigsetup()
{
logIBASE = log((double)IBASE);
log10IBASE = log10((double)IBASE);
big_one = make(INT, 1, 0);
}
-int isnat(word x)
+int isnat(x)
+word x;
{
return (tag[x] == INT && big_is_positive(x));
}
-/* store C long long as mira bigint */
-word sto_int(long long i)
+word sto_int(i) /* store C long long as mira bigint */
+long long i;
{
- word s = 0;
-
- if (i < 0) {
- s = SIGNBIT;
- i = -i;
- }
-
- word x = make(INT, (s | i) & MAXDIGIT, 0);
-
- for (word *p = &rest(x); i >>= DIGITWIDTH; p = &rest(*p)) {
- *p = make(INT, i & MAXDIGIT, 0);
+ word s, x;
+ if (i < 0)
+ s = SIGNBIT, i = -i;
+ else
+ s = 0;
+ x = make(INT, s | i & MAXDIGIT, 0);
+ if (i >>= DIGITWIDTH) {
+ word *p = &rest(x);
+ *p = make(INT, i & MAXDIGIT, 0), p = &rest(*p);
+ while (i >>= DIGITWIDTH)
+ *p = make(INT, i & MAXDIGIT, 0), p = &rest(*p);
}
+ return (x);
+} /* change to long long, DT Oct 2019 */
- return x;
-}
+#define maxval (1ll<<60)
-/* mira bigint to C long long */
-long long get_int(word x)
+long long get_int(x) /* mira bigint to C long long */
+word x;
{
long long n = digit0(x);
word sign = big_is_negative(x);
-
- if (!(x = rest(x))) {
+ if (!(x = rest(x)))
+ return (sign ? -n : n);
+ {
+ word w = DIGITWIDTH;
+ while (x && w < 60)
+ n += (long long)digit(x) << w, w += DIGITWIDTH, x = rest(x);
+ if (x)
+ n = maxval; /* overflow, return large value */
return (sign ? -n : n);
}
+} /* change to long long, DT Oct 2019 */
- for (word w = DIGITWIDTH; x && w < 60; w += DIGITWIDTH) {
- n += (long long)digit(x) << w;
- x = rest(x);
- }
-
- if (x) {
- n = maxval; /* overflow, return large value */
- }
-
- return (sign ? -n : n);
-}
-
-word bignegate(word x)
+word bignegate(x)
+word x;
{
- if (bigzero(x)) {
- return x;
- }
- return make(INT, (hd[x] & SIGNBIT) ? (hd[x] & MAXDIGIT) : (SIGNBIT | hd[x]),
- tl[x]);
+ if (bigzero(x))
+ return (x);
+ return (make
+ (INT, hd[x] & SIGNBIT ? hd[x] & MAXDIGIT : SIGNBIT | hd[x], tl[x]));
}
-word bigplus(word x, word y)
+word bigplus(x, y)
+word x, y;
{
- if (big_is_positive(x)) {
- if (big_is_positive(y)) {
+ if (big_is_positive(x))
+ if (big_is_positive(y))
return (big_plus(x, y, 0));
- } else {
+ else
return (big_sub(x, y));
- }
-
- } else if (big_is_positive(y)) {
+ else if (big_is_positive(y))
return (big_sub(y, x));
- }
-
- return (big_plus(x, y, SIGNBIT)); /* both negative */
+ else
+ return (big_plus(x, y, SIGNBIT)); /* both negative */
}
-/* ignore input signs, treat x,y as positive */
-word big_plus(word x, word y, int signbit)
+word big_plus(x, y, signbit) /* ignore input signs, treat x,y as positive */
+word x, y;
+int signbit;
{
word d = digit0(x) + digit0(y);
word carry = ((d & IBASE) != 0);
- word r = make(INT, (signbit | d) & MAXDIGIT, 0); /* result */
+ word r = make(INT, signbit | d & MAXDIGIT, 0); /* result */
word *z = &rest(r); /* pointer to rest of result */
x = rest(x);
y = rest(y);
-
while (x && y) { /* this loop has been unwrapped once, see above */
d = carry + digit(x) + digit(y);
carry = ((d & IBASE) != 0);
@@ -130,11 +118,8 @@ word big_plus(word x, word y, int signbit)
y = rest(y);
z = &rest(*z);
}
-
- if (y) {
+ if (y)
x = y; /* by convention x is the longer one */
- }
-
while (x) {
d = carry + digit(x);
carry = ((d & IBASE) != 0);
@@ -142,32 +127,27 @@ word big_plus(word x, word y, int signbit)
x = rest(x);
z = &rest(*z);
}
-
- if (carry) {
+ if (carry)
*z = make(INT, 1, 0);
- }
-
- return r;
+ return (r);
}
-word bigsub(word x, word y)
+word bigsub(x, y)
+word x, y;
{
- if (big_is_positive(x)) {
- if (big_is_positive(y)) {
+ if (big_is_positive(x))
+ if (big_is_positive(y))
return (big_sub(x, y));
- } else {
- return (big_plus(x, y, 0)); /* poz x, negative y */
- }
-
- } else if (big_is_positive(y)) {
- return (big_plus(x, y, SIGNBIT)); /* negative x, poz y */
- }
-
- return (big_sub(y, x)); /* both negative */
+ else
+ return (big_plus(x, y, 0)); /* big_is_positive x, negative y */
+ else if (big_is_positive(y))
+ return (big_plus(x, y, SIGNBIT)); /* negative x, big_is_positive y */
+ else
+ return (big_sub(y, x)); /* both negative */
}
-/* ignore input signs, treat x,y as positive */
-word big_sub(word x, word y)
+word big_sub(x, y) /* ignore input signs, treat x,y as positive */
+word x, y;
{
word d = digit0(x) - digit0(y);
word borrow = (d & IBASE) != 0;
@@ -176,7 +156,6 @@ word big_sub(word x, word y)
word *p = NULL; /* pointer to trailing zeros, if any */
x = rest(x);
y = rest(y);
-
while (x && y) { /* this loop has been unwrapped once, see above */
d = digit(x) - digit(y) - borrow;
borrow = (d & IBASE) != 0;
@@ -190,7 +169,6 @@ word big_sub(word x, word y)
y = rest(y);
z = &rest(*z);
}
-
while (y) { /* at most one of these two loops will be invoked */
d = -digit(y) - borrow;
borrow = ((d & IBASE) != 0);
@@ -203,21 +181,18 @@ word big_sub(word x, word y)
y = rest(y);
z = &rest(*z);
}
-
while (x) { /* alternative loop */
d = digit(x) - borrow;
borrow = ((d & IBASE) != 0);
d = d & MAXDIGIT;
*z = make(INT, d, 0);
- if (d) {
+ if (d)
p = NULL;
- } else if (!p) {
+ else if (!p)
p = z;
- }
x = rest(x);
z = &rest(*z);
}
-
if (borrow) { /* result is negative - take complement and add 1 */
p = NULL;
d = (digit(r) ^ MAXDIGIT) + 1;
@@ -228,141 +203,107 @@ word big_sub(word x, word y)
d = (digit(*z) ^ MAXDIGIT) + borrow;
borrow = ((d & IBASE) != 0);
digit(*z) = d = d & MAXDIGIT;
- if (d) {
+ if (d)
p = NULL;
- } else if (!p) {
+ else if (!p)
p = z;
- }
z = &rest(*z);
}
}
-
- if (p) {
+ if (p)
*p = 0; /* remove redundant (ie trailing) zeros */
- }
-
- return r;
+ return (r);
}
-/* returns +ve,0,-ve as x greater than, equal, less than y */
-int bigcmp(word x, word y)
+int bigcmp(x, y) /* returns +ve,0,-ve as x greater than, equal, less than y */
+word x, y;
{
- word s = big_is_negative(x);
-
- if (big_is_negative(y) != s) {
+ word d, r, s = big_is_negative(x);
+ if (big_is_negative(y) != s)
return (s ? -1 : 1);
- }
-
- word r = digit0(x) - digit0(y);
+ r = digit0(x) - digit0(y);
for (;;) {
x = rest(x);
y = rest(y);
-
- if (!x) {
- if (y) {
+ if (!x)
+ if (y)
return (s ? 1 : -1);
- } else {
+ else
return (s ? -r : r);
- }
- }
-
- if (!y) {
+ if (!y)
return (s ? -1 : 1);
- }
-
- word d = digit(x) - digit(y);
- if (d) {
+ d = digit(x) - digit(y);
+ if (d)
r = d;
- }
}
}
-/* naive multiply - quadratic */
-word bigtimes(word x, word y)
+word bigtimes(x, y) /* naive multiply - quadratic */
+word x, y;
{
- word r = make(INT, 0, 0);
-
if (len(x) < len(y)) {
word hold = x;
x = y;
y = hold;
} /* important optimisation */
-
- if (bigzero(x)) {
- return r; /* short cut */
- }
-
+ word r = make(INT, 0, 0);
word d = digit0(y);
word s = big_is_negative(y);
word n = 0;
-
+ if (bigzero(x))
+ return (r); /* short cut */
for (;;) {
- if (d) {
+ if (d)
r = bigplus(r, shift(n, stimes(x, d)));
- }
-
n++;
y = rest(y);
-
- if (!y) {
+ if (!y)
return (s != big_is_negative(x) ? bignegate(r) : r);
- }
-
d = digit(y);
}
}
-/* multiply big x by n'th power of IBASE */
-/* NB - we assume x non-zero, else unnormalised result */
-word shift(word n, word x)
+word shift(n, x) /* multiply big x by n'th power of IBASE */
+word n, x;
{
- while (n--) {
+ while (n--)
x = make(INT, 0, x);
- }
- return x;
-}
+ return (x);
+} /* NB - we assume x non-zero, else unnormalised result */
-/* multiply big x (>=0) by digit n (>0) */
-word stimes(word x, word n)
+word stimes(x, n) /* multiply big x (>=0) by digit n (>0) */
+word x, n;
{
unsigned d = n * digit0(x); /* ignore sign of x */
word carry = d >> DIGITWIDTH;
word r = make(INT, d & MAXDIGIT, 0);
word *y = &rest(r);
-
- while ((x = rest(x))) {
- d = n * digit(x) + carry;
- *y = make(INT, d & MAXDIGIT, 0);
- y = &rest(*y);
- carry = d >> DIGITWIDTH;
- }
-
- if (carry) {
+ while (x = rest(x))
+ d = n * digit(x) + carry,
+ *y = make(INT, d & MAXDIGIT, 0), y = &rest(*y), carry = d >> DIGITWIDTH;
+ if (carry)
*y = make(INT, carry, 0);
- }
-
- return r;
+ return (r);
}
-/* may assume y~=0 */
-word bigdiv(word x, word y)
+word b_rem; /* contains remainder from last call to longdiv or shortdiv */
+
+word bigdiv(x, y) /* may assume y~=0 */
+word x, y;
{
+ word s1, s2, q;
/* make x,y positive and remember signs */
- word s1 = big_is_negative(y);
- if (s1) {
+ if (s1 = big_is_negative(y))
y = make(INT, digit0(y), rest(y));
- }
-
- word s2 = s1;
- if (big_is_negative(x)) {
- x = make(INT, digit0(x), rest(x));
- s2 = !s1;
- }
-
+ if (big_is_negative(x))
+ x = make(INT, digit0(x), rest(x)), s2 = !s1;
+ else
+ s2 = s1;
/* effect: s1 set iff y negative, s2 set iff signs mixed */
-
- word q = rest(y) ? longdiv(x, y) : shortdiv(x, digit(y));
-
+ if (rest(y))
+ q = longdiv(x, y);
+ else
+ q = shortdiv(x, digit(y));
if (s2) {
if (!bigzero(b_rem)) {
x = q;
@@ -371,50 +312,36 @@ word bigdiv(word x, word y)
if (!rest(x)) {
rest(x) = make(INT, 1, 0);
break;
- } else {
+ } else
x = rest(x);
- }
}
}
-
- if (!bigzero(q)) {
+ if (!bigzero(q))
digit(q) = SIGNBIT | digit(q);
- }
}
-
- return q;
+ return (q);
}
-/* may assume y~=0 */
-word bigmod(word x, word y)
+word bigmod(x, y) /* may assume y~=0 */
+word x, y;
{
- word s1 = big_is_negative(y);
-
+ word s1, s2;
/* make x,y positive and remember signs */
- if (s1) {
+ if (s1 = big_is_negative(y))
y = make(INT, digit0(y), rest(y));
- }
-
- word s2 = s1;
- if (big_is_negative(x)) {
- x = make(INT, digit0(x), rest(x));
- s2 = !s1;
- }
-
+ if (big_is_negative(x))
+ x = make(INT, digit0(x), rest(x)), s2 = !s1;
+ else
+ s2 = s1;
/* effect: s1 set iff y negative, s2 set iff signs mixed */
-
- if (rest(y)) {
+ if (rest(y))
longdiv(x, y);
- } else {
+ else
shortdiv(x, digit(y));
- }
-
if (s2) {
- if (!bigzero(b_rem)) {
+ if (!bigzero(b_rem))
b_rem = bigsub(y, b_rem);
- }
}
-
return (s1 ? bignegate(b_rem) : b_rem);
}
@@ -427,233 +354,184 @@ word bigmod(word x, word y)
magnitudes invariant under change of sign, remainder has sign of
dividend, quotient negative if signs of divi(sor/dend) mixed */
-/* divide big x by single digit n returning big quotient
- and setting external b_rem as side effect */
-
-/* may assume - x>=0,n>0 */
-word shortdiv(word x, word n)
-{
- word d = digit(x);
- word q = 0;
-
- while ((x = rest(x))) { /* reverse rest(x) into q */
- q = make(INT, d, q);
- d = digit(x); /* leaving most sig. digit in d */
- }
-
- x = q;
- d = d / n;
-
- if (d || !q) {
- q = make(INT, d, 0); /* put back first digit (if not leading 0) */
- } else {
- q = 0;
- }
-
- word s_rem = d % n;
- while (x) { /* in situ division of q by n AND destructive reversal */
- d = s_rem * IBASE + digit(x);
- digit(x) = d / n;
+word shortdiv(x, n) /* divide big x by single digit n returning big quotient
+ and setting external b_rem as side effect */
+ /* may assume - x>=0,n>0 */
+word x, n;
+{
+ word d = digit(x), s_rem, q = 0;
+ while (x = rest(x)) /* reverse rest(x) into q */
+ q = make(INT, d, q), d = digit(x); /* leaving most sig. digit in d */
+ {
+ word tmp;
+ x = q;
s_rem = d % n;
- word tmp = x;
- x = rest(x);
- rest(tmp) = q;
- q = tmp;
+ d = d / n;
+ if (d || !q)
+ q = make(INT, d, 0); /* put back first digit (if not leading 0) */
+ else
+ q = 0;
+ while (x) /* in situ division of q by n AND destructive reversal */
+ d = s_rem * IBASE + digit(x), digit(x) = d / n, s_rem = d % n,
+ tmp = x, x = rest(x), rest(tmp) = q, q = tmp;
}
-
b_rem = make(INT, s_rem, 0);
-
- return q;
+ return (q);
}
-/* divide big x by big y returning quotient, leaving
- remainder in extern variable b_rem */
-/* may assume - x>=0,y>0 */
-word longdiv(word x, word y)
+word longdiv(x, y) /* divide big x by big y returning quotient, leaving
+ remainder in extern variable b_rem */
+ /* may assume - x>=0,y>0 */
+word x, y;
{
+ word n, q, ly, y1, scale;
if (bigcmp(x, y) < 0) {
b_rem = x;
return (make(INT, 0, 0));
}
-
- word y1 = msd(y);
- word scale = IBASE / (y1 + 1);
- if (scale > 1) { /* rescale if necessary */
- x = stimes(x, scale);
- y = stimes(y, scale);
- y1 = msd(y);
- }
-
- word n = 0;
- word ly = len(y);
- while (bigcmp(x, y = make(INT, 0, y)) >= 0) {
+ y1 = msd(y);
+ if ((scale = IBASE / (y1 + 1)) > 1) /* rescale if necessary */
+ x = stimes(x, scale), y = stimes(y, scale), y1 = msd(y);
+ n = q = 0;
+ ly = len(y);
+ while (bigcmp(x, y = make(INT, 0, y)) >= 0)
n++;
- }
-
y = rest(y); /* want largest y not exceeding x */
ly += n;
- word q = 0;
for (;;) {
- word d = 0;
- word lx = len(x);
-
- if (lx == ly) {
- if (bigcmp(x, y) >= 0) {
- x = bigsub(x, y);
- d = 1;
- }
-
- } else if (lx > ly) {
+ word d, lx = len(x);
+ if (lx < ly)
+ d = 0;
+ else if (lx == ly)
+ if (bigcmp(x, y) >= 0)
+ x = bigsub(x, y), d = 1;
+ else
+ d = 0;
+ else {
d = ms2d(x) / y1;
- if (d > MAXDIGIT) {
+ if (d > MAXDIGIT)
d = MAXDIGIT;
- }
-
- if ((d -= 2) > 0) {
+ if ((d -= 2) > 0)
x = bigsub(x, stimes(y, d));
- } else {
+ else
d = 0;
- }
-
if (bigcmp(x, y) >= 0) {
- x = bigsub(x, y);
- d++;
- if (bigcmp(x, y) >= 0) {
- x = bigsub(x, y);
- d++;
- }
+ x = bigsub(x, y), d++;
+ if (bigcmp(x, y) >= 0)
+ x = bigsub(x, y), d++;
}
}
-
q = make(INT, d, q);
if (n-- == 0) {
- b_rem = (scale == 1) ? x : shortdiv(x, scale);
- return q;
+ b_rem = scale == 1 ? x : shortdiv(x, scale);
+ return (q);
}
-
ly--;
y = rest(y);
}
} /* see Bird & Wadler p82 for explanation */
-/* no of digits in big x */
-word len(word x)
+word len(x) /* no of digits in big x */
+word x;
{
word n = 1;
- while ((x = rest(x))) {
+ while (x = rest(x))
n++;
- }
- return n;
+ return (n);
}
-/* most significant digit of big x */
-word msd(word x)
+word msd(x) /* most significant digit of big x */
+word x;
{
- while (rest(x)) {
+ while (rest(x))
x = rest(x);
- }
- return digit(x); /* sign? */
+ return (digit(x)); /* sign? */
}
-/* most significant 2 digits of big x (len>=2) */
-word ms2d(word x)
+word ms2d(x) /* most significant 2 digits of big x (len>=2) */
+word x;
{
word d = digit(x);
x = rest(x);
- while (rest(x)) {
- d = digit(x);
- x = rest(x);
- }
+ while (rest(x))
+ d = digit(x), x = rest(x);
return (digit(x) * IBASE + d);
}
-/* assumes y poz */
-word bigpow(word x, word y)
+word bigpow(x, y) /* assumes y big_is_positive */
+word x, y;
{
- word r = make(INT, 1, 0);
- while (rest(y)) {
+ word d, r = make(INT, 1, 0);
+ while (rest(y)) { /* this loop has been unwrapped once, see below */
word i = DIGITWIDTH;
- word d = digit(y);
+ d = digit(y);
while (i--) {
- if (d & 1) {
+ if (d & 1)
r = bigtimes(r, x);
- }
x = bigtimes(x, x);
d >>= 1;
}
y = rest(y);
}
-
- word d = digit(y);
- if (d & 1) {
+ d = digit(y);
+ if (d & 1)
r = bigtimes(r, x);
- }
-
while (d >>= 1) {
x = bigtimes(x, x);
- if (d & 1) {
+ if (d & 1)
r = bigtimes(r, x);
- }
}
-
- return r;
+ return (r);
}
-double bigtodbl(word x)
+double bigtodbl(x)
+word x;
{
- double b = 1.0;
- double r = (double)digit0(x);
-
+ word s = big_is_negative(x);
+ double b = 1.0, r = (double)digit0(x);
x = rest(x);
- while (x) {
- b = b * IBASE;
- r = r + b * digit(x);
- x = rest(x);
- }
-
- return big_is_negative(x) ? -r : r;
+ while (x)
+ b = b * IBASE, r = r + b * digit(x), x = rest(x);
+ if (s)
+ return (-r);
+ return (r);
} /* small end first */
/* note: can return oo, -oo
but is used without surrounding sto_/set)dbl() only in compare() */
-/* not currently used */
-long double bigtoldbl(word x)
-{
- long double b = 1.0L;
- long double r = (long double)digit0(x);
-
+/* not currently used
+long double bigtoldbl(x)
+word x;
+{ int s=big_is_negative(x);
+ long double b=1.0L, r=digit0(x);
x = rest(x);
- while (x) {
- b = b * IBASE;
- r = r + b * digit(x);
- x = rest(x);
- }
-
- return big_is_negative(x) ? -r : r;
+ while(x)b=b*IBASE,r=r+b*digit(x),x=rest(x);
+/*printf("bigtoldbl returns %Le\n",s?-r:r); /* DEBUG
+ if(s)return(-r);
+ return(r);
} /* not compatible with std=c90, lib fns eg sqrtl broken */
-/* entier */
-word dbltobig(double x)
+word dbltobig(x) /* entier */
+double x;
{
+ word s = (x < 0);
word r = make(INT, 0, 0);
word *p = &r;
- double y = fabs(floor(x));
-
- while (y > 0.0) {
+ double y = floor(x);
+/*if(fabs(y-x+1.0)<1e-9)y += 1.0; /* trick due to Peter Bartke, see note */
+ for (y = fabs(y);;) {
double n = fmod(y, (double)IBASE);
digit(*p) = (word) n;
y = (y - n) / (double)IBASE;
- if (y > 0.0) {
+ if (y > 0.0)
rest(*p) = make(INT, 0, 0), p = &rest(*p);
- }
+ else
+ break;
}
-
- if (x < 0) {
+ if (s)
digit(r) = SIGNBIT | digit(r);
- }
-
- return r;
+ return (r);
}
/* produces junk in low order digits if x exceeds range in which integer
@@ -668,125 +546,91 @@ word dbltobig(double x)
the exact integers. There are inherent deficiences in 64 bit fp,
no point in trying to mask this */
-/* logarithm of big x */
-double biglog(word x)
+double biglog(x) /* logarithm of big x */
+word x;
{
- if (big_is_negative(x) || bigzero(x)) {
- errno = EDOM;
- math_error("log");
- }
-
word n = 0;
double r = digit(x);
- while ((x = rest(x))) {
- n++;
- r = digit(x) + r / IBASE;
- }
-
+ if (big_is_negative(x) || bigzero(x))
+ errno = EDOM, math_error("log");
+ while (x = rest(x))
+ n++, r = digit(x) + r / IBASE;
return (log(r) + n * logIBASE);
}
-/* logarithm of big x */
-double biglog10(word x)
+double biglog10(x) /* logarithm of big x */
+word x;
{
- if (big_is_negative(x) || bigzero(x)) {
- errno = EDOM;
- math_error("log10");
- }
-
word n = 0;
double r = digit(x);
- while ((x = rest(x))) {
+ if (big_is_negative(x) || bigzero(x))
+ errno = EDOM, math_error("log10");
+ while (x = rest(x))
n++, r = digit(x) + r / IBASE;
- }
-
return (log10(r) + n * log10IBASE);
}
-/* NB does NOT check for malformed number, assumes already done */
-/* p is a pointer to a null terminated string of digits */
-/* read a big number (in decimal) */
-word bigscan(char *p)
+word bigscan(p) /* read a big number (in decimal) */
+ /* NB does NOT check for malformed number, assumes already done */
+char *p; /* p is a pointer to a null terminated string of digits */
{
- word s = 0;
- word r = make(INT, 0, 0);
- if (*p == '-') {
- s = 1;
- p++; /* optional leading `-' (for NUMVAL) */
- }
-
+ word s = 0, r = make(INT, 0, 0);
+ if (*p == '-')
+ s = 1, p++; /* optional leading `-' (for NUMVAL) */
while (*p) {
- word d = *p - '0';
- word f = 10;
+ word d = *p - '0', f = 10;
p++;
-
- while (*p && f < PTEN) {
- d = 10 * d + *p - '0';
- f = 10 * f;
- p++;
- }
-
+ while (*p && f < PTEN)
+ d = 10 * d + *p - '0', f = 10 * f, p++;
/* rest of loop does r=f*r+d; (in situ) */
d = f * digit(r) + d;
- word carry = d >> DIGITWIDTH;
- word *x = &rest(r);
- digit(r) = d & MAXDIGIT;
-
- while (*x) {
- d = f * digit(*x) + carry;
- digit(*x) = d & MAXDIGIT;
- carry = d >> DIGITWIDTH;
- x = &rest(*x);
- }
-
- if (carry) {
- *x = make(INT, carry, 0);
+ {
+ word carry = d >> DIGITWIDTH;
+ word *x = &rest(r);
+ digit(r) = d & MAXDIGIT;
+ while (*x)
+ d = f * digit(*x) + carry,
+ digit(*x) = d & MAXDIGIT, carry = d >> DIGITWIDTH, x = &rest(*x);
+ if (carry)
+ *x = make(INT, carry, 0);
}
}
-
- if (s && !bigzero(r)) {
+/*if(*p=='e')
+ { int s=bigscan(p+1);
+ r = bigtimes(r,bigpow(make(INT,10,0),s); } */
+ if (s && !bigzero(r))
digit(r) = digit(r) | SIGNBIT;
- }
-
- return r;
+ return (r);
}
/* code to handle (unsigned) exponent commented out */
-/* assumes redundant leading zeros removed */
-/* read unsigned hex number in '\0'-terminated string p to q */
-word bigxscan(char *p, char *q)
+word bigxscan(p, q) /* read unsigned hex number in '\0'-terminated string p to q */
+ /* assumes redundant leading zeros removed */
+char *p, *q;
{
- if (*p == '0' && !p[1]) {
- return make(INT, 0, 0);
- }
-
word r; /* will hold result */
word *x = &r;
-
+ if (*p == '0' && !p[1])
+ return make(INT, 0, 0);
while (q > p) {
unsigned long long hold;
- q = (q - p < 15) ? p : (q - 15); /* read upto 15 hex digits from small end */
+ q = q - p < 15 ? p : q - 15; /* read upto 15 hex digits from small end */
sscanf(q, "%llx", &hold);
*q = '\0';
word count = 4; /* 15 hex digits => 4 bignum digits */
- while (count-- && !(hold == 0 && q == p)) {
- *x = make(INT, hold & MAXDIGIT, 0);
- hold >>= DIGITWIDTH;
- x = &rest(*x);
- }
+ while (count-- && !(hold == 0 && q == p))
+ *x = make(INT, hold & MAXDIGIT, 0), hold >>= DIGITWIDTH, x = &rest(*x);
}
-
return r;
}
-/* read unsigned octal number in '\0'-terminated string p to q */
-word bigoscan(char *p, char *q)
+word bigoscan(p, q) /* read unsigned octal number in '\0'-terminated string p to q */
/* assumes redundant leading zeros removed */
+char *p, *q;
{
word r; /* will hold result */
word *x = &r;
-
while (q > p) {
unsigned hold;
q = q - p < 5 ? p : q - 5; /* read (upto) 5 octal digits from small end */
@@ -794,201 +638,201 @@ word bigoscan(char *p, char *q)
*q = '\0';
*x = make(INT, hold, 0), x = &rest(*x);
}
-
return r;
}
-word digitval(char c)
+word digitval(c)
+char c;
{
return isdigit(c) ? c - '0' : isupper(c) ? 10 + c - 'A' : 10 + c - 'a';
}
-/* numeral (as Miranda string) to big number */
-/* does NOT check for malformed numeral, assumes
- done and that z fully evaluated */
-word strtobig(word z, int base)
+word strtobig(z, base) /* numeral (as Miranda string) to big number */
+ /* does NOT check for malformed numeral, assumes
+ done and that z fully evaluated */
+word z;
+int base;
{
- word s = 0;
- word r = make(INT, 0, 0);
-
- word PBASE = PTEN;
- if (base == 16) {
+ word s = 0, r = make(INT, 0, 0), PBASE = PTEN;
+ if (base == 16)
PBASE = PSIXTEEN;
- } else if (base == 8) {
+ else if (base == 8)
PBASE = PEIGHT;
- }
-
- if (z != NIL && hd[z] == '-') {
- s = 1;
- z = tl[z]; /* optional leading `-' (for NUMVAL) */
- }
-
- if (base != 10) {
+ if (z != NIL && hd[z] == '-')
+ s = 1, z = tl[z]; /* optional leading `-' (for NUMVAL) */
+ if (base != 10)
z = tl[tl[z]]; /* remove "0x" or "0o" */
- }
-
while (z != NIL) {
word d = digitval(hd[z]), f = base;
z = tl[z];
-
- while (z != NIL && f < PBASE) {
- d = base * d + digitval(hd[z]);
- f = base * f;
- z = tl[z];
- }
-
+ while (z != NIL && f < PBASE)
+ d = base * d + digitval(hd[z]), f = base * f, z = tl[z];
/* rest of loop does r=f*r+d; (in situ) */
d = f * digit(r) + d;
- word carry = d >> DIGITWIDTH;
- word *x = &rest(r);
- digit(r) = d & MAXDIGIT;
-
- while (*x) {
- d = f * digit(*x) + carry,
- digit(*x) = d & MAXDIGIT, carry = d >> DIGITWIDTH, x = &rest(*x);
- }
-
- if (carry) {
- *x = make(INT, carry, 0);
+ {
+ word carry = d >> DIGITWIDTH;
+ word *x = &rest(r);
+ digit(r) = d & MAXDIGIT;
+ while (*x)
+ d = f * digit(*x) + carry,
+ digit(*x) = d & MAXDIGIT, carry = d >> DIGITWIDTH, x = &rest(*x);
+ if (carry)
+ *x = make(INT, carry, 0);
}
}
-
- if (s && !bigzero(r)) {
+ if (s && !bigzero(r))
digit(r) = digit(r) | SIGNBIT;
- }
-
- return r;
+ return (r);
}
+extern char *dicp;
-word bigtostr(word x) /* number to decimal string (as Miranda list) */
+word bigtostr(x) /* number to decimal string (as Miranda list) */
+word x;
{
- word sign = big_is_negative(x);
- word s = NIL;
-
+ word x1, sign, s = NIL;
#ifdef DEBUG
extern int debug;
if (debug & 04) { /* print octally */
word d = digit0(x);
+ sign = big_is_negative(x);
for (;;) {
word i = OCTW;
- while (i-- || d) {
+ while (i-- || d)
s = cons('0' + (d & 07), s), d >>= 3;
- }
x = rest(x);
- if (x) {
+ if (x)
s = cons(' ', s), d = digit(x);
- } else {
+ else
return (sign ? cons('-', s) : s);
- }
}
}
#endif
-
if (rest(x) == 0) {
extern char *dicp;
sprintf(dicp, "%ld", getsmallint(x));
return (str_conv(dicp));
}
-
- word x1 = make(INT, digit0(x), 0); /* reverse x into x1 */
- while ((x = rest(x))) {
+ sign = big_is_negative(x);
+ x1 = make(INT, digit0(x), 0); /* reverse x into x1 */
+ while (x = rest(x))
x1 = make(INT, digit(x), x1);
- }
-
x = x1;
for (;;) { /* in situ division of (reversed order) x by PTEN */
- word d = digit(x);
- word rem = d % PTEN;
+ word d = digit(x), rem = d % PTEN;
d = d / PTEN;
x1 = rest(x);
-
- if (d) {
+ if (d)
digit(x) = d;
- } else {
+ else
x = x1; /* remove leading zero from result */
- }
-
- while (x1) {
- d = rem * IBASE + digit(x1);
- digit(x1) = d / PTEN;
- rem = d % PTEN;
- x1 = rest(x1);
- }
-
+ while (x1)
+ d = rem * IBASE + digit(x1),
+ digit(x1) = d / PTEN, rem = d % PTEN, x1 = rest(x1);
/* end of in situ division (also uses x1 as temporary) */
if (x) {
word i = TENW;
- while (i--) {
+ while (i--)
s = cons('0' + rem % 10, s), rem = rem / 10;
- }
} else {
- while (rem) {
+ while (rem)
s = cons('0' + rem % 10, s), rem = rem / 10;
- }
return (sign ? cons('-', s) : s);
}
}
}
-/* integer to hexadecimal string (as Miranda list) */
-word bigtostrx(word x)
+word bigtostrx(x) /* integer to hexadecimal string (as Miranda list) */
+word x;
{
- word r = NIL;
- word sign = big_is_negative(x);
+ word r = NIL, s = big_is_negative(x);
while (x) {
word count = 4; /* 60 bits => 20 octal digits => 4 bignum digits */
unsigned long long factor = 1;
unsigned long long hold = 0;
- while (count-- && x) { /* calculate value of (upto) 4 bignum digits */
+ while (count-- && x) /* calculate value of (upto) 4 bignum digits */
hold = hold + factor * digit0(x),
+ /* printf("::%llx\n",hold), /* DEBUG */
factor <<= 15, x = rest(x);
- }
sprintf(dicp, "%.15llx", hold); /* 15 hex digits = 60 bits */
+ /* printf(":::%s\n",dicp); /* DEBUG */
char *q = dicp + 15;
- while (--q >= dicp) {
+ while (--q >= dicp)
r = cons(*q, r);
- }
}
-
- while (digit(r) == '0' && rest(r) != NIL) {
+ while (digit(r) == '0' && rest(r) != NIL)
r = rest(r); /* remove redundant leading 0's */
- }
-
r = cons('0', cons('x', r));
- if (sign) {
+ if (s)
r = cons('-', r);
- }
-
- return r;
+ return (r);
}
-/* integer to octal string (as Miranda list) */
-word bigtostr8(word x)
+word bigtostr8(x) /* integer to octal string (as Miranda list) */
+word x;
{
- word r = NIL;
- word sign = big_is_negative(x);
-
+ word r = NIL, s = big_is_negative(x);
while (x) {
char *q = dicp + 5;
sprintf(dicp, "%.5lo", digit0(x));
- while (--q >= dicp) {
+ while (--q >= dicp)
r = cons(*q, r);
- }
x = rest(x);
}
-
- while (digit(r) == '0' && rest(r) != NIL) {
+ while (digit(r) == '0' && rest(r) != NIL)
r = rest(r); /* remove redundant leading 0's */
- }
-
r = cons('0', cons('o', r));
-
- if (sign) {
+ if (s)
r = cons('-', r);
- }
+ return (r);
+}
- return r;
+#ifdef DEBUG
+wff(x) /* check for well-formation of integer */
+word x;
+{
+ word y = x;
+ if (tag[x] != INT)
+ printf("BAD TAG %d\n", tag[x]);
+ if (big_is_negative(x) && !digit0(x) && !rest(x))
+ printf("NEGATIVE ZERO!\n");
+ if (digit0(x) & (~MAXDIGIT))
+ printf("OVERSIZED DIGIT!\n");
+ while (x = rest(x))
+ if (tag[x] != INT)
+ printf("BAD INTERNAL TAG %d\n", tag[x]);
+ else if (digit(x) & (~MAXDIGIT))
+ printf("OVERSIZED DIGIT!\n");
+ else if (!digit(x) && !rest(x))
+ printf("TRAILING ZERO!\n");
+ return (y);
+}
+
+normalise(x) /* remove trailing zeros */
+word x;
+{
+ if (rest(x))
+ rest(x) = norm1(rest(x));
+ return (wff(x));
}
+norm1(x)
+word x;
+{
+ if (rest(x))
+ rest(x) = norm1(rest(x));
+ return (!digit(x) && !rest(x) ? 0 : x);
+}
+
+#endif
+
+/* stall(s)
+char *s;
+{ fprintf(stderr,"big integer %s not yet implemented\n",s);
+ exit(0);
+}
+
+#define destrev(x,y,z) while(x)z=x,x=rest(x),rest(z)=y,y=z;
+/* destructively reverse x into y using z as temp */
+
/* END OF MIRANDA INTEGER PACKAGE */