Mercurial > hg > xemacs-beta
diff src/number-mp.c @ 1983:9c872f33ecbe
[xemacs-hg @ 2004-04-05 22:49:31 by james]
Add bignum, ratio, and bigfloat support.
author | james |
---|---|
date | Mon, 05 Apr 2004 22:50:11 +0000 |
parents | |
children | 4f5e0297b73f |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/src/number-mp.c Mon Apr 05 22:50:11 2004 +0000 @@ -0,0 +1,509 @@ +/* Numeric types for XEmacs using the MP library. + Copyright (C) 2004 Jerry James. + +This file is part of XEmacs. + +XEmacs is free software; you can redistribute it and/or modify it +under the terms of the GNU General Public License as published by the +Free Software Foundation; either version 2, or (at your option) any +later version. + +XEmacs 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 General Public License +for more details. + +You should have received a copy of the GNU General Public License +along with XEmacs; see the file COPYING. If not, write to +the Free Software Foundation, Inc., 59 Temple Place - Suite 330, +Boston, MA 02111-1307, USA. */ + +/* Synched up with: Not in FSF. */ + +#include <config.h> +#include <limits.h> +#include <math.h> +#include "lisp.h" + +static MINT *bignum_bytesize, *bignum_long_sign_bit, *bignum_one, *bignum_two; +MINT *bignum_zero, *intern_bignum; +MINT *bignum_min_int, *bignum_max_int, *bignum_max_uint; +MINT *bignum_min_long, *bignum_max_long, *bignum_max_ulong; +short div_rem; + +char * +bignum_to_string (bignum b, int base) +{ + REGISTER unsigned int i; + unsigned int bufsize = 128U, index = 0U; + int sign; + char *buffer = xnew_array (char, 128), *retval; + MINT *quo = MP_ITOM (0); + short rem; + + /* FIXME: signal something if base is < 2 or doesn't fit into a short. */ + + /* Save the sign for later */ + sign = MP_MCMP (b, bignum_zero); + + if (sign == 0) + { + XREALLOC_ARRAY (buffer, char, 2); + buffer[0] = '0'; + buffer[1] = '\0'; + return buffer; + } + /* Copy abs(b) into quo for destructive modification */ + else if (sign < 0) + MP_MSUB (bignum_zero, b, quo); + else + MP_MOVE (b, quo); + + quo = MP_ITOM (0); + + /* Loop over the digits of b (in BASE) and place each one into buffer */ + for (i = 0U; MP_MCMP(quo, bignum_zero) > 0; i++) + { + MP_SDIV (quo, base, quo, &rem); + if (index == bufsize) + { + bufsize <<= 1; + XREALLOC_ARRAY (buffer, char, bufsize); + } + buffer[index++] = rem < 10 ? rem + '0' : rem - 10 + 'a'; + } + MP_MFREE (quo); + + /* Reverse the digits, maybe add a minus sign, and add a null terminator */ + bufsize = index + (sign < 0 ? 1 : 0) + 1; + retval = xnew_array (char, bufsize); + if (sign < 0) + { + retval[0] = '-'; + i = 1; + } + else + i = 0; + for (; i < bufsize - 1; i++) + retval[i] = buffer[--index]; + retval[bufsize - 1] = '\0'; + xfree (buffer, char *); + return retval; +} + +#define BIGNUM_TO_TYPE(type,accumtype) do { \ + MP_MULT (b, quo, quo); \ + for (i = 0U; i < sizeof(type); i++) \ + { \ + MP_SDIV (quo, 256, quo, &rem); \ + retval |= ((accumtype) rem) << (8 * i); \ + } \ + MP_MFREE (quo); \ +} while (0) + +int +bignum_to_int (bignum b) +{ + short rem, sign; + unsigned int retval = 0; + REGISTER unsigned int i; + MINT *quo; + + sign = MP_MCMP (b, bignum_zero) < 0 ? -1 : 1; + quo = MP_ITOM (sign); + BIGNUM_TO_TYPE (int, unsigned int); + return ((int) retval) * sign; +} + +unsigned int +bignum_to_uint (bignum b) +{ + short rem; + unsigned int retval = 0U; + REGISTER unsigned int i; + MINT *quo; + + quo = MP_ITOM (MP_MCMP (b, bignum_zero) < 0 ? -1 : 1); + BIGNUM_TO_TYPE (unsigned int, unsigned int); + return retval; +} + +long +bignum_to_long (bignum b) +{ + short rem, sign; + unsigned long retval = 0L; + REGISTER unsigned int i; + MINT *quo; + + sign = MP_MCMP (b, bignum_zero) < 0 ? -1 : 1; + quo = MP_ITOM (sign); + BIGNUM_TO_TYPE (long, unsigned long); + return ((long) retval) * sign; +} + +unsigned long +bignum_to_ulong (bignum b) +{ + short rem; + unsigned long retval = 0UL; + REGISTER unsigned int i; + MINT *quo; + + quo = MP_ITOM (MP_MCMP (b, bignum_zero) < 0 ? -1 : 1); + BIGNUM_TO_TYPE (unsigned long, unsigned long); + return retval; +} + +double +bignum_to_double (bignum b) +{ + short rem, sign; + double retval = 0.0; + REGISTER unsigned int i; + MINT *quo; + + sign = MP_MCMP (b, bignum_zero) < 0 ? -1 : 1; + quo = MP_ITOM (sign); + MP_MULT (b, quo, quo); + for (i = 0U; MP_MCMP(quo, bignum_zero) > 0; i++) + { + MP_SDIV (quo, 256, quo, &rem); + retval += rem * i * 256; + } + MP_MFREE (quo); + return retval * sign; +} + +static short +char_to_number (char c) +{ + if (c >= '0' && c <= '9') + return c - '0'; + if (c >= 'a' && c <= 'z') + return c - 'a' + 10; + if (c >= 'A' && c <= 'Z') + return c - 'A' + 10; + return -1; +} + +int +bignum_set_string (bignum b, const char *s, int base) +{ + MINT *mbase; + short digit; + + if (base == 0) + { + if (s[0] == '0' && (s[1] == 'x' || s[1] == 'X')) + { + base = 16; + s += 2; + } + else if (*s == '0') + { + base = 8; + s++; + } + else + base = 10; + } + + /* FIXME: signal something if base is < 2 or doesn't fit into a short. */ + + mbase = MP_ITOM ((short) base); + MP_MOVE (bignum_zero, b); + + for (digit = char_to_number (*s); digit >= 0 && digit < base; + digit = char_to_number (*++s)) + { + MINT *temp; + + MP_MULT (b, mbase, b); + temp = MP_ITOM (digit); + MP_MADD (b, temp, b); + MP_MFREE (temp); + } + + return (digit >= 0) ? -1 : 0; +} + +void +bignum_set_long (MINT *b, long l) +{ + /* Negative l is hard, not least because -LONG_MIN == LONG_MIN. We pretend + that l is unsigned, then subtract off the amount equal to the sign bit. */ + bignum_set_ulong (b, (unsigned long) l); + if (l < 0L) + MP_MSUB (b, bignum_long_sign_bit, b); +} + +void +bignum_set_ulong (bignum b, unsigned long l) +{ + REGISTER unsigned int i; + MINT *multiplier = MP_ITOM (1); + + MP_MOVE (bignum_zero, b); + for (i = 0UL; l > 0UL; l >>= 8, i++) + { + MINT *temp = MP_ITOM ((short) (l & 255)); + MP_MULT (multiplier, temp, temp); + MP_MADD (b, temp, b); + MP_MULT (multiplier, bignum_bytesize, multiplier); + MP_MFREE (temp); + } + MP_MFREE (multiplier); +} + +void +bignum_set_double (bignum b, double d) +{ + REGISTER unsigned int i; + int negative = (d < 0) ? 1 : 0; + MINT *multiplier = MP_ITOM (1); + + MP_MOVE (bignum_zero, b); + if (negative) + d = -d; + for (i = 0UL; d > 0.0; d /= 256, i++) + { + MINT *temp = MP_ITOM ((short) fmod (d, 256.0)); + MP_MULT (multiplier, temp, temp); + MP_MADD (b, temp, b); + MP_MULT (multiplier, bignum_bytesize, multiplier); + MP_MFREE (temp); + } + MP_MFREE (multiplier); + if (negative) + MP_MSUB (bignum_zero, b, b); +} + +/* Return nonzero if b1 is exactly divisible by b2 */ +int +bignum_divisible_p (bignum b1, bignum b2) +{ + int retval; + MINT *rem = MP_ITOM (0); + MP_MDIV (b1, b2, intern_bignum, rem); + retval = (MP_MCMP (rem, bignum_zero) == 0); + MP_MFREE (rem); + return retval; +} + +void bignum_ceil (bignum quotient, bignum N, bignum D) +{ + MP_MDIV (N, D, quotient, intern_bignum); + if (MP_MCMP (intern_bignum, bignum_zero) > 0 && + MP_MCMP (quotient, bignum_zero) > 0) + MP_MADD (quotient, bignum_one, quotient); +} + +void bignum_floor (bignum quotient, bignum N, bignum D) +{ + MP_MDIV (N, D, quotient, intern_bignum); + if (MP_MCMP (intern_bignum, bignum_zero) > 0 && + MP_MCMP (quotient, bignum_zero) < 0) + MP_MSUB (quotient, bignum_one, quotient); +} + +/* RESULT = N to the POWth power */ +void +bignum_pow (bignum result, bignum n, unsigned long pow) +{ + MP_MOVE (bignum_one, result); + for ( ; pow > 0UL; pow--) + MP_MULT (result, n, result); +} + +/* lcm(b1,b2) = b1 * b2 / gcd(b1, b2) */ +void +bignum_lcm (bignum result, bignum b1, bignum b2) +{ + MP_MULT (b1, b2, result); + MP_GCD (b1, b2, intern_bignum); + MP_MDIV (result, intern_bignum, result, intern_bignum); +} + +/* FIXME: We can't handle negative args, so right now we just make them + positive before doing anything else. How should we really handle negative + args? */ +#define bignum_bit_op(result, b1, b2, op) \ + REGISTER unsigned int i; \ + MINT *multiplier = MP_ITOM (1), *n1 = MP_ITOM (0), *n2 = MP_ITOM (0); \ + \ + if (MP_MCMP (bignum_zero, b1) > 0) \ + MP_MSUB (bignum_zero, b1, n1); \ + else \ + MP_MOVE (b1, n1); \ + if (MP_MCMP (bignum_zero, b2) > 0) \ + MP_MSUB (bignum_zero, b2, n2); \ + else \ + MP_MOVE (b2, n2); \ + \ + MP_MOVE (bignum_zero, result); \ + \ + for (i = 0UL; MP_MCMP (bignum_zero, n1) < 0 && \ + MP_MCMP (bignum_zero, n2) < 0; i++) \ + { \ + short byte1, byte2; \ + MINT *temp; \ + \ + MP_SDIV (n1, 256, n1, &byte1); \ + MP_SDIV (n2, 256, n2, &byte2); \ + temp = MP_ITOM (byte1 op byte2); \ + MP_MULT (multiplier, temp, temp); \ + MP_MADD (result, temp, result); \ + MP_MULT (multiplier, bignum_bytesize, multiplier); \ + MP_MFREE (temp); \ + } \ + MP_MFREE (n2); \ + MP_MFREE (n1); \ + MP_MFREE (multiplier) + +void +bignum_and (bignum result, bignum b1, bignum b2) +{ + bignum_bit_op (result, b1, b2, &); +} + +void +bignum_ior (bignum result, bignum b1, bignum b2) +{ + bignum_bit_op (result, b1, b2, |); +} + +void +bignum_xor (bignum result, bignum b1, bignum b2) +{ + bignum_bit_op (result, b1, b2, ^); +} + +/* NOT is not well-defined for bignums ... where do you stop flipping bits? + We just flip until we see the last one. This is probably a bad idea. */ +void +bignum_not (bignum result, bignum b) +{ + REGISTER unsigned int i; + MINT *multiplier = MP_ITOM (1), *n = MP_ITOM (0); + + if (MP_MCMP (bignum_zero, b) > 0) + MP_MSUB (bignum_zero, b, n); + else + MP_MOVE (b, n); + + MP_MOVE (bignum_zero, result); + + for (i = 0UL; MP_MCMP (bignum_zero, n) < 0; i++) + { + short byte; + MINT *temp; + + MP_SDIV (n, 256, n, &byte); + temp = MP_ITOM (~byte); + MP_MULT (multiplier, temp, temp); + MP_MADD (result, temp, result); + MP_MULT (multiplier, bignum_bytesize, multiplier); + MP_MFREE (temp); + } + MP_MFREE (n); + MP_MFREE (multiplier); +} + +void +bignum_setbit (bignum b, unsigned long bit) +{ + bignum_pow (intern_bignum, bignum_two, bit); + bignum_ior (b, b, intern_bignum); +} + +/* This is so evil, even I feel queasy. */ +void +bignum_clrbit (bignum b, unsigned long bit) +{ + MINT *num = MP_ITOM (0); + + /* See if the bit is already set, and subtract it off if not */ + MP_MOVE (b, intern_bignum); + bignum_pow (num, bignum_two, bit); + bignum_ior (intern_bignum, intern_bignum, num); + if (MP_MCMP (b, intern_bignum) == 0) + MP_MSUB (b, num, b); + MP_MFREE (num); +} + +int +bignum_testbit (bignum b, unsigned long bit) +{ + bignum_pow (intern_bignum, bignum_two, bit); + bignum_and (intern_bignum, b, intern_bignum); + return MP_MCMP (intern_bignum, bignum_zero); +} + +void +bignum_lshift (bignum result, bignum b, unsigned long bits) +{ + bignum_pow (intern_bignum, bignum_two, bits); + MP_MULT (b, intern_bignum, result); +} + +void +bignum_rshift (bignum result, bignum b, unsigned long bits) +{ + bignum_pow (intern_bignum, bignum_two, bits); + MP_MDIV (b, intern_bignum, result, intern_bignum); +} + +void bignum_random_seed(unsigned long seed) +{ + /* FIXME: Implement me */ +} + +void bignum_random(bignum result, bignum limit) +{ + /* FIXME: Implement me */ + MP_MOVE (bignum_zero, result); +} + +void +init_number_mp () +{ + REGISTER unsigned int i; + + bignum_zero = MP_ITOM (0); + bignum_one = MP_ITOM (1); + bignum_two = MP_ITOM (2); + + /* intern_bignum holds throwaway values from macro expansions in + number-mp.h. Its value is immaterial. */ + intern_bignum = MP_ITOM (0); + + /* bignum_bytesize holds the number of bits in a byte. */ + bignum_bytesize = MP_ITOM (256); + + /* bignum_long_sign_bit holds an adjustment for negative longs. */ + bignum_long_sign_bit = MP_ITOM (256); + for (i = 1UL; i < sizeof (long); i++) + MP_MULT (bignum_bytesize, bignum_long_sign_bit, bignum_long_sign_bit); + + /* The MP interface only supports turning short ints into MINTs, so we have + to set these the hard way. */ + + bignum_min_int = MP_ITOM (0); + bignum_set_long (bignum_min_int, INT_MIN); + + bignum_max_int = MP_ITOM (0); + bignum_set_long (bignum_max_int, INT_MAX); + + bignum_max_uint = MP_ITOM (0); + bignum_set_ulong (bignum_max_uint, UINT_MAX); + + bignum_min_long = MP_ITOM (0); + bignum_set_long (bignum_min_long, LONG_MIN); + + bignum_max_long = MP_ITOM (0); + bignum_set_long (bignum_max_long, LONG_MAX); + + bignum_max_ulong = MP_ITOM (0); + bignum_set_ulong (bignum_max_ulong, ULONG_MAX); +}