c-tier bignum — arbitrary-precision integers unblock secp256k1 widths
Adds tagged bignum support alongside the existing 48-bit fixnum on the C tier. Tag 6 = bignum, heap struct sign-magnitude with u64 little-endian limbs. Reader emits bignums for any literal past the fixnum range; +, -, *, quotient, remainder, modulo, expt, =, <, >, abs, odd?, even?, integer?, exact?, number->string, string->number all promote fixnum → bignum on overflow & demote back when results fit. Boehm GC owns every allocation. Schoolbook O(n²) mul + shift-subtract divmod is sufficient at our 4-limb / 256-bit scale. Before: (expt 2 48) = 0, (expt 2 256) = 0, secp256k1-p = -4294968273. After: all three return their exact arbitrary-precision values, matching Python tier byte-for-byte. Validated: - c/test.c — 85/85 pass (+2 new bignum unit tests). - tests/functional.lsp — 205/205 pass on both C & Python tiers. - tests/bignum-cross-tier.lsp — 33/33 pass byte-identical on both tiers (diff produces no output). - ecdsa/runs/lumbda-sweep-003/c-tier-bignum-probe.lsp — all four assertions now match the Python oracle. - ecdsa Phase B byte-identity sweep inside QEMU guest: n+1=9 p=251 sha256 c668bbe3... — matches Python oracle. n+1=18 p=131071 sha256 8a031f96... — matches Python oracle. n+1=33 p=2³²-5 sha256 0bc56905... — matches Python oracle. Previously the n+1=33 C tier emitted sha256 b024d6d9... (26,078 fewer Toffolis due to silent fixnum wrap). Bignums close that gate. secp256k1 production-width emit (n+1=257) is now structurally unblocked on C tier; downstream agent (#55) drives that next-step on the ecdsa side. Asm tier inherits in a follow-up port.
This commit is contained in:
parent
cf68c0da15
commit
2f342c3be2
11 changed files with 918 additions and 29 deletions
|
|
@ -17,7 +17,7 @@ CFLAGS += -DUSE_BOEHM_GC
|
|||
LDFLAGS += -lgc
|
||||
endif
|
||||
|
||||
SRCS = types.c reader.c printer.c eval.c builtins.c vm.c jit.c portal.c
|
||||
SRCS = types.c bignum.c reader.c printer.c eval.c builtins.c vm.c jit.c portal.c
|
||||
OBJS = $(SRCS:.c=.o)
|
||||
|
||||
.PHONY: all clean test bench
|
||||
|
|
|
|||
496
c/bignum.c
Normal file
496
c/bignum.c
Normal file
|
|
@ -0,0 +1,496 @@
|
|||
/*
|
||||
* bignum.c — Arbitrary-precision integer arithmetic for lumbda C tier.
|
||||
*
|
||||
* Design (full notes in lumbda.h after the Rational block):
|
||||
* - Sign-magnitude representation: int32_t sign (-1, 0, +1),
|
||||
* uint32_t n_limbs, uint64_t *limbs (little-endian magnitude).
|
||||
* - bignum_normalize() trims leading-zero limbs & demotes to fixnum
|
||||
* when |value| <= FIXNUM_MAX so callers never branch on tag.
|
||||
* - All arithmetic uses 64-bit limbs with 128-bit products (__int128)
|
||||
* for schoolbook multiply. Division uses long-division by single
|
||||
* limb when divisor fits one limb, otherwise schoolbook with
|
||||
* normalization (sufficient for our scale: secp256k1 needs at most
|
||||
* 4 limbs / 256 bits, far below where Karatsuba helps).
|
||||
* - Boehm GC owns the Bignum struct & its limbs array.
|
||||
*
|
||||
* Cross-tier oracle: Python tier (native arbitrary int) is the byte-
|
||||
* identity standard. Every operation here must agree with Python on
|
||||
* matching inputs — see tests/bignum-cross-tier.lsp.
|
||||
*/
|
||||
#include "lumbda.h"
|
||||
|
||||
/* ═══════════════════════════════════════════════════════════════════════════
|
||||
* Allocation helpers
|
||||
* ═══════════════════════════════════════════════════════════════════════════ */
|
||||
|
||||
static Bignum *alloc_bignum(uint32_t n_limbs) {
|
||||
Bignum *b = (Bignum *)ul_malloc(sizeof(Bignum));
|
||||
b->hdr.type = OBJ_BIGNUM;
|
||||
b->sign = 0;
|
||||
b->n_limbs = n_limbs;
|
||||
if (n_limbs == 0) {
|
||||
b->limbs = NULL;
|
||||
} else {
|
||||
b->limbs = (uint64_t *)ul_malloc(sizeof(uint64_t) * n_limbs);
|
||||
memset(b->limbs, 0, sizeof(uint64_t) * n_limbs);
|
||||
}
|
||||
return b;
|
||||
}
|
||||
|
||||
static Value wrap_bignum(Bignum *b) {
|
||||
return NANBOX(TAG_BIGNUM, (uintptr_t)b);
|
||||
}
|
||||
|
||||
/* Trim leading zero limbs in place. Zero magnitude => sign 0. */
|
||||
static void trim(Bignum *b) {
|
||||
while (b->n_limbs > 0 && b->limbs[b->n_limbs - 1] == 0)
|
||||
b->n_limbs--;
|
||||
if (b->n_limbs == 0) b->sign = 0;
|
||||
}
|
||||
|
||||
/* ═══════════════════════════════════════════════════════════════════════════
|
||||
* Constructors
|
||||
* ═══════════════════════════════════════════════════════════════════════════ */
|
||||
|
||||
Value make_bignum_from_i64(int64_t v) {
|
||||
if (FITS_FIXNUM(v)) return VAL_INT(v);
|
||||
Bignum *b = alloc_bignum(1);
|
||||
if (v < 0) {
|
||||
b->sign = -1;
|
||||
/* Avoid -INT64_MIN UB: cast then negate via unsigned. */
|
||||
b->limbs[0] = (uint64_t)(-(v + 1)) + 1ULL;
|
||||
} else if (v > 0) {
|
||||
b->sign = 1;
|
||||
b->limbs[0] = (uint64_t)v;
|
||||
} else {
|
||||
b->sign = 0;
|
||||
b->n_limbs = 0;
|
||||
}
|
||||
return wrap_bignum(b);
|
||||
}
|
||||
|
||||
Value bignum_from_limbs(int32_t sign, uint64_t *limbs, uint32_t n) {
|
||||
Bignum *b = alloc_bignum(n);
|
||||
if (n > 0) memcpy(b->limbs, limbs, sizeof(uint64_t) * n);
|
||||
b->sign = sign;
|
||||
trim(b);
|
||||
return wrap_bignum(b);
|
||||
}
|
||||
|
||||
/* Demote bignum-shape to fixnum when |value| fits 48 signed bits.
|
||||
* On demote we return VAL_INT(...); otherwise allocate a fresh Bignum. */
|
||||
Value bignum_normalize(int32_t sign, uint64_t *limbs, uint32_t n) {
|
||||
/* Strip top zeros (compute effective limb count). */
|
||||
while (n > 0 && limbs[n - 1] == 0) n--;
|
||||
if (n == 0) return VAL_INT(0);
|
||||
if (n == 1) {
|
||||
uint64_t lo = limbs[0];
|
||||
/* Signed 48-bit range: -2^47 .. 2^47 - 1. */
|
||||
if (sign >= 0 && lo <= (uint64_t)FIXNUM_MAX)
|
||||
return VAL_INT((int64_t)lo);
|
||||
if (sign < 0 && lo <= (uint64_t)(1ULL << 47))
|
||||
return VAL_INT(-(int64_t)lo);
|
||||
}
|
||||
Bignum *b = alloc_bignum(n);
|
||||
memcpy(b->limbs, limbs, sizeof(uint64_t) * n);
|
||||
b->sign = sign;
|
||||
return wrap_bignum(b);
|
||||
}
|
||||
|
||||
/* Promote a fixnum to a Bignum struct (caller-owned). Bignum input passes
|
||||
* through (the same struct is returned — Bignum structs are immutable from
|
||||
* the perspective of arithmetic so sharing is safe). */
|
||||
Bignum *as_bignum(Value v) {
|
||||
if (IS_BIGNUM(v)) return AS_BIGNUM(v);
|
||||
if (IS_INT(v)) {
|
||||
int64_t x = as_int(v);
|
||||
Bignum *b = alloc_bignum(x == 0 ? 0 : 1);
|
||||
if (x > 0) { b->sign = 1; b->limbs[0] = (uint64_t)x; }
|
||||
else if (x < 0) {
|
||||
b->sign = -1;
|
||||
b->limbs[0] = (uint64_t)(-(x + 1)) + 1ULL;
|
||||
}
|
||||
return b;
|
||||
}
|
||||
if (IS_DOUBLE(v)) {
|
||||
double d = as_double(v);
|
||||
return AS_BIGNUM(make_bignum_from_i64((int64_t)d));
|
||||
}
|
||||
lisp_error("not an integer");
|
||||
return NULL;
|
||||
}
|
||||
|
||||
Value to_bignum(Value v) {
|
||||
if (IS_BIGNUM(v)) return v;
|
||||
return wrap_bignum(as_bignum(v));
|
||||
}
|
||||
|
||||
/* ═══════════════════════════════════════════════════════════════════════════
|
||||
* Parsing (string → bignum)
|
||||
* ═══════════════════════════════════════════════════════════════════════════ */
|
||||
|
||||
/* Multiply magnitude by single limb; add another limb. Magnitude grows. */
|
||||
static void mag_mul_add(uint64_t *limbs, uint32_t *n, uint64_t mul, uint64_t add) {
|
||||
__uint128_t carry = add;
|
||||
for (uint32_t i = 0; i < *n; i++) {
|
||||
__uint128_t prod = (__uint128_t)limbs[i] * mul + carry;
|
||||
limbs[i] = (uint64_t)prod;
|
||||
carry = prod >> 64;
|
||||
}
|
||||
while (carry) {
|
||||
limbs[(*n)++] = (uint64_t)carry;
|
||||
carry >>= 64;
|
||||
}
|
||||
}
|
||||
|
||||
Value make_bignum_from_str(const char *s, int base) {
|
||||
int32_t sign = 1;
|
||||
if (*s == '-') { sign = -1; s++; }
|
||||
else if (*s == '+') s++;
|
||||
/* Conservative limb count: log2(10) ≈ 3.33 bits per digit → roughly
|
||||
* len/19 limbs for base 10; len/16 limbs for base 16. Round up + 1. */
|
||||
size_t len = strlen(s);
|
||||
uint32_t cap = (uint32_t)((len * 4) / 19) + 2;
|
||||
uint64_t *limbs = (uint64_t *)ul_malloc(sizeof(uint64_t) * cap);
|
||||
memset(limbs, 0, sizeof(uint64_t) * cap);
|
||||
uint32_t n = 0;
|
||||
for (; *s; s++) {
|
||||
int d;
|
||||
if (*s >= '0' && *s <= '9') d = *s - '0';
|
||||
else if (*s >= 'a' && *s <= 'f') d = *s - 'a' + 10;
|
||||
else if (*s >= 'A' && *s <= 'F') d = *s - 'A' + 10;
|
||||
else { ul_free(limbs); lisp_error("bad digit '%c' in number", *s); }
|
||||
if (d >= base) { ul_free(limbs); lisp_error("digit out of base"); }
|
||||
mag_mul_add(limbs, &n, (uint64_t)base, (uint64_t)d);
|
||||
if (n >= cap - 1) {
|
||||
cap *= 2;
|
||||
limbs = (uint64_t *)ul_realloc(limbs, sizeof(uint64_t) * cap);
|
||||
memset(limbs + n, 0, sizeof(uint64_t) * (cap - n));
|
||||
}
|
||||
}
|
||||
Value v = bignum_normalize(sign, limbs, n);
|
||||
ul_free(limbs);
|
||||
return v;
|
||||
}
|
||||
|
||||
/* ═══════════════════════════════════════════════════════════════════════════
|
||||
* Magnitude-level primitives
|
||||
* ═══════════════════════════════════════════════════════════════════════════ */
|
||||
|
||||
/* Compare magnitudes. Returns -1, 0, +1. */
|
||||
static int mag_cmp(const uint64_t *a, uint32_t na, const uint64_t *b, uint32_t nb) {
|
||||
if (na != nb) return na < nb ? -1 : 1;
|
||||
for (int32_t i = (int32_t)na - 1; i >= 0; i--) {
|
||||
if (a[i] != b[i]) return a[i] < b[i] ? -1 : 1;
|
||||
}
|
||||
return 0;
|
||||
}
|
||||
|
||||
/* out = a + b magnitudes. out_cap must be at least max(na, nb) + 1. */
|
||||
static uint32_t mag_add(const uint64_t *a, uint32_t na,
|
||||
const uint64_t *b, uint32_t nb,
|
||||
uint64_t *out) {
|
||||
uint32_t n = na > nb ? na : nb;
|
||||
__uint128_t carry = 0;
|
||||
for (uint32_t i = 0; i < n; i++) {
|
||||
__uint128_t s = carry;
|
||||
if (i < na) s += a[i];
|
||||
if (i < nb) s += b[i];
|
||||
out[i] = (uint64_t)s;
|
||||
carry = s >> 64;
|
||||
}
|
||||
if (carry) out[n++] = (uint64_t)carry;
|
||||
return n;
|
||||
}
|
||||
|
||||
/* out = a - b magnitudes, assuming a >= b. out_cap >= na. */
|
||||
static uint32_t mag_sub(const uint64_t *a, uint32_t na,
|
||||
const uint64_t *b, uint32_t nb,
|
||||
uint64_t *out) {
|
||||
int64_t borrow = 0;
|
||||
for (uint32_t i = 0; i < na; i++) {
|
||||
__int128 s = (__int128)a[i] - borrow;
|
||||
if (i < nb) s -= b[i];
|
||||
if (s < 0) { s += ((__int128)1 << 64); borrow = 1; }
|
||||
else borrow = 0;
|
||||
out[i] = (uint64_t)s;
|
||||
}
|
||||
while (na > 0 && out[na - 1] == 0) na--;
|
||||
return na;
|
||||
}
|
||||
|
||||
/* out = a * b magnitudes. out_cap >= na + nb. */
|
||||
static uint32_t mag_mul(const uint64_t *a, uint32_t na,
|
||||
const uint64_t *b, uint32_t nb,
|
||||
uint64_t *out) {
|
||||
memset(out, 0, sizeof(uint64_t) * (na + nb));
|
||||
for (uint32_t i = 0; i < na; i++) {
|
||||
__uint128_t carry = 0;
|
||||
for (uint32_t j = 0; j < nb; j++) {
|
||||
__uint128_t prod = (__uint128_t)a[i] * b[j] + out[i + j] + carry;
|
||||
out[i + j] = (uint64_t)prod;
|
||||
carry = prod >> 64;
|
||||
}
|
||||
out[i + nb] += (uint64_t)carry;
|
||||
}
|
||||
uint32_t n = na + nb;
|
||||
while (n > 0 && out[n - 1] == 0) n--;
|
||||
return n;
|
||||
}
|
||||
|
||||
/* Divide magnitude a by single 64-bit limb d. Quotient into q (cap >= na),
|
||||
* remainder returned. */
|
||||
static uint64_t mag_divmod_u64(const uint64_t *a, uint32_t na,
|
||||
uint64_t d, uint64_t *q, uint32_t *nq) {
|
||||
__uint128_t r = 0;
|
||||
for (int32_t i = (int32_t)na - 1; i >= 0; i--) {
|
||||
__uint128_t cur = (r << 64) | a[i];
|
||||
q[i] = (uint64_t)(cur / d);
|
||||
r = cur % d;
|
||||
}
|
||||
*nq = na;
|
||||
while (*nq > 0 && q[*nq - 1] == 0) (*nq)--;
|
||||
return (uint64_t)r;
|
||||
}
|
||||
|
||||
/* Schoolbook division: a / b → (q, r) where a = q*b + r, 0 <= r < b.
|
||||
* For our scale (operands <= 4 limbs commonly, up to a few dozen for
|
||||
* intermediates in expt), a basic Knuth Algorithm D is overkill — we use
|
||||
* shift-subtract on bits. O(na * nb * 64). */
|
||||
static void mag_divmod(const uint64_t *a, uint32_t na,
|
||||
const uint64_t *b, uint32_t nb,
|
||||
uint64_t *q, uint32_t *nq,
|
||||
uint64_t *r, uint32_t *nr) {
|
||||
if (nb == 0) lisp_error("division by zero");
|
||||
if (mag_cmp(a, na, b, nb) < 0) {
|
||||
/* a < b: q=0, r=a */
|
||||
*nq = 0;
|
||||
memcpy(r, a, sizeof(uint64_t) * na);
|
||||
*nr = na;
|
||||
return;
|
||||
}
|
||||
if (nb == 1) {
|
||||
uint64_t rem = mag_divmod_u64(a, na, b[0], q, nq);
|
||||
if (rem) { r[0] = rem; *nr = 1; } else { *nr = 0; }
|
||||
return;
|
||||
}
|
||||
/* Shift-subtract long division: walk bits MSB → LSB of dividend. */
|
||||
memset(q, 0, sizeof(uint64_t) * na);
|
||||
memset(r, 0, sizeof(uint64_t) * (nb + 1));
|
||||
uint32_t rn = 0;
|
||||
for (int64_t bit = (int64_t)na * 64 - 1; bit >= 0; bit--) {
|
||||
/* r <<= 1 */
|
||||
uint64_t carry = 0;
|
||||
for (uint32_t i = 0; i <= nb; i++) {
|
||||
uint64_t new_carry = r[i] >> 63;
|
||||
r[i] = (r[i] << 1) | carry;
|
||||
carry = new_carry;
|
||||
}
|
||||
if (carry || r[nb]) rn = nb + 1;
|
||||
else { while (rn > 0 && r[rn - 1] == 0) rn--; }
|
||||
/* r |= bit of a */
|
||||
uint64_t bit_val = (a[bit / 64] >> (bit % 64)) & 1ULL;
|
||||
r[0] |= bit_val;
|
||||
if (bit_val && rn == 0) rn = 1;
|
||||
/* if r >= b: r -= b; set bit in q */
|
||||
if (mag_cmp(r, rn, b, nb) >= 0) {
|
||||
uint32_t new_rn = mag_sub(r, rn, b, nb, r);
|
||||
rn = new_rn;
|
||||
q[bit / 64] |= 1ULL << (bit % 64);
|
||||
}
|
||||
}
|
||||
*nq = na;
|
||||
while (*nq > 0 && q[*nq - 1] == 0) (*nq)--;
|
||||
*nr = rn;
|
||||
}
|
||||
|
||||
/* ═══════════════════════════════════════════════════════════════════════════
|
||||
* Signed arithmetic (sign-magnitude wrappers)
|
||||
* ═══════════════════════════════════════════════════════════════════════════ */
|
||||
|
||||
Value big_add(Value a, Value b) {
|
||||
Bignum *ba = as_bignum(a), *bb = as_bignum(b);
|
||||
if (ba->sign == 0) return bignum_normalize(bb->sign, bb->limbs, bb->n_limbs);
|
||||
if (bb->sign == 0) return bignum_normalize(ba->sign, ba->limbs, ba->n_limbs);
|
||||
uint32_t cap = (ba->n_limbs > bb->n_limbs ? ba->n_limbs : bb->n_limbs) + 1;
|
||||
uint64_t *out = (uint64_t *)ul_malloc(sizeof(uint64_t) * cap);
|
||||
Value res;
|
||||
if (ba->sign == bb->sign) {
|
||||
uint32_t n = mag_add(ba->limbs, ba->n_limbs, bb->limbs, bb->n_limbs, out);
|
||||
res = bignum_normalize(ba->sign, out, n);
|
||||
} else {
|
||||
int cmp = mag_cmp(ba->limbs, ba->n_limbs, bb->limbs, bb->n_limbs);
|
||||
if (cmp == 0) { res = VAL_INT(0); }
|
||||
else if (cmp > 0) {
|
||||
uint32_t n = mag_sub(ba->limbs, ba->n_limbs, bb->limbs, bb->n_limbs, out);
|
||||
res = bignum_normalize(ba->sign, out, n);
|
||||
} else {
|
||||
uint32_t n = mag_sub(bb->limbs, bb->n_limbs, ba->limbs, ba->n_limbs, out);
|
||||
res = bignum_normalize(bb->sign, out, n);
|
||||
}
|
||||
}
|
||||
ul_free(out);
|
||||
return res;
|
||||
}
|
||||
|
||||
Value big_neg(Value a) {
|
||||
if (IS_INT(a)) {
|
||||
int64_t x = as_int(a);
|
||||
return make_bignum_from_i64(-x);
|
||||
}
|
||||
Bignum *ba = as_bignum(a);
|
||||
if (ba->sign == 0) return VAL_INT(0);
|
||||
return bignum_normalize(-ba->sign, ba->limbs, ba->n_limbs);
|
||||
}
|
||||
|
||||
Value big_sub(Value a, Value b) {
|
||||
return big_add(a, big_neg(b));
|
||||
}
|
||||
|
||||
Value big_mul(Value a, Value b) {
|
||||
Bignum *ba = as_bignum(a), *bb = as_bignum(b);
|
||||
if (ba->sign == 0 || bb->sign == 0) return VAL_INT(0);
|
||||
uint32_t cap = ba->n_limbs + bb->n_limbs;
|
||||
uint64_t *out = (uint64_t *)ul_malloc(sizeof(uint64_t) * cap);
|
||||
uint32_t n = mag_mul(ba->limbs, ba->n_limbs, bb->limbs, bb->n_limbs, out);
|
||||
Value res = bignum_normalize(ba->sign * bb->sign, out, n);
|
||||
ul_free(out);
|
||||
return res;
|
||||
}
|
||||
|
||||
bool big_is_zero(Value v) {
|
||||
if (IS_INT(v)) return as_int(v) == 0;
|
||||
if (IS_BIGNUM(v)) return AS_BIGNUM(v)->sign == 0;
|
||||
return false;
|
||||
}
|
||||
|
||||
bool big_is_odd(Value v) {
|
||||
if (IS_INT(v)) return (as_int(v) & 1) != 0;
|
||||
if (IS_BIGNUM(v)) {
|
||||
Bignum *b = AS_BIGNUM(v);
|
||||
return b->n_limbs > 0 && (b->limbs[0] & 1ULL) != 0;
|
||||
}
|
||||
return false;
|
||||
}
|
||||
|
||||
Value big_abs(Value v) {
|
||||
if (IS_INT(v)) { int64_t x = as_int(v); return make_bignum_from_i64(x < 0 ? -x : x); }
|
||||
Bignum *b = as_bignum(v);
|
||||
if (b->sign >= 0) return v;
|
||||
return bignum_normalize(1, b->limbs, b->n_limbs);
|
||||
}
|
||||
|
||||
int big_cmp(Value a, Value b) {
|
||||
Bignum *ba = as_bignum(a), *bb = as_bignum(b);
|
||||
if (ba->sign != bb->sign) return ba->sign < bb->sign ? -1 : 1;
|
||||
if (ba->sign == 0) return 0;
|
||||
int m = mag_cmp(ba->limbs, ba->n_limbs, bb->limbs, bb->n_limbs);
|
||||
return ba->sign > 0 ? m : -m;
|
||||
}
|
||||
|
||||
/* Quotient: truncated toward zero (matches C / Python's // for non-negatives;
|
||||
* Scheme's `quotient` is also truncate-toward-zero). */
|
||||
Value big_quotient(Value a, Value b) {
|
||||
Bignum *ba = as_bignum(a), *bb = as_bignum(b);
|
||||
if (bb->sign == 0) lisp_error("division by zero");
|
||||
if (ba->sign == 0) return VAL_INT(0);
|
||||
uint32_t qcap = ba->n_limbs > 0 ? ba->n_limbs : 1;
|
||||
uint32_t rcap = bb->n_limbs + 1;
|
||||
uint64_t *q = (uint64_t *)ul_malloc(sizeof(uint64_t) * qcap);
|
||||
uint64_t *r = (uint64_t *)ul_malloc(sizeof(uint64_t) * rcap);
|
||||
uint32_t nq, nr;
|
||||
mag_divmod(ba->limbs, ba->n_limbs, bb->limbs, bb->n_limbs, q, &nq, r, &nr);
|
||||
int32_t sign = ba->sign * bb->sign;
|
||||
Value res = bignum_normalize(nq == 0 ? 0 : sign, q, nq);
|
||||
ul_free(q); ul_free(r);
|
||||
return res;
|
||||
}
|
||||
|
||||
/* Remainder: result has sign of dividend (Scheme `remainder`). */
|
||||
Value big_remainder(Value a, Value b) {
|
||||
Bignum *ba = as_bignum(a), *bb = as_bignum(b);
|
||||
if (bb->sign == 0) lisp_error("division by zero");
|
||||
if (ba->sign == 0) return VAL_INT(0);
|
||||
uint32_t qcap = ba->n_limbs > 0 ? ba->n_limbs : 1;
|
||||
uint32_t rcap = bb->n_limbs + 1;
|
||||
uint64_t *q = (uint64_t *)ul_malloc(sizeof(uint64_t) * qcap);
|
||||
uint64_t *r = (uint64_t *)ul_malloc(sizeof(uint64_t) * rcap);
|
||||
uint32_t nq, nr;
|
||||
mag_divmod(ba->limbs, ba->n_limbs, bb->limbs, bb->n_limbs, q, &nq, r, &nr);
|
||||
Value res = bignum_normalize(nr == 0 ? 0 : ba->sign, r, nr);
|
||||
ul_free(q); ul_free(r);
|
||||
return res;
|
||||
}
|
||||
|
||||
/* Modulo: result has sign of divisor (Scheme `modulo`). */
|
||||
Value big_modulo(Value a, Value b) {
|
||||
Value r = big_remainder(a, b);
|
||||
/* If r != 0 and sign(r) != sign(b), add b. */
|
||||
if (big_is_zero(r)) return r;
|
||||
Bignum *br = as_bignum(r), *bb = as_bignum(b);
|
||||
if (br->sign != bb->sign) return big_add(r, b);
|
||||
return r;
|
||||
}
|
||||
|
||||
/* Binary exponentiation. exp must be a non-negative integer. */
|
||||
Value big_expt(Value base, Value exp) {
|
||||
/* exp must be an integer. */
|
||||
if (!IS_INTEGER(exp)) lisp_error("expt: integer exponent required");
|
||||
/* Negative exponent → not handled here (caller falls back to double). */
|
||||
if (IS_INT(exp) && as_int(exp) < 0) lisp_error("expt: negative exponent");
|
||||
if (IS_BIGNUM(exp) && AS_BIGNUM(exp)->sign < 0) lisp_error("expt: negative exponent");
|
||||
/* Read exp as u64 — bignum exponents larger than 2^63 would never finish. */
|
||||
uint64_t e;
|
||||
if (IS_INT(exp)) e = (uint64_t)as_int(exp);
|
||||
else {
|
||||
Bignum *be = AS_BIGNUM(exp);
|
||||
if (be->n_limbs > 1) lisp_error("expt: exponent too large");
|
||||
e = be->n_limbs == 0 ? 0 : be->limbs[0];
|
||||
}
|
||||
Value result = VAL_INT(1);
|
||||
Value b = base;
|
||||
while (e > 0) {
|
||||
if (e & 1) result = big_mul(result, b);
|
||||
e >>= 1;
|
||||
if (e) b = big_mul(b, b);
|
||||
}
|
||||
return result;
|
||||
}
|
||||
|
||||
/* ═══════════════════════════════════════════════════════════════════════════
|
||||
* Decimal rendering
|
||||
* ═══════════════════════════════════════════════════════════════════════════ */
|
||||
|
||||
char *bignum_to_str(Bignum *b) {
|
||||
if (b->sign == 0) return ul_strdup("0");
|
||||
/* Divide repeatedly by 10^18 → 18 decimal digits per chunk. */
|
||||
uint32_t n = b->n_limbs;
|
||||
uint64_t *mag = (uint64_t *)ul_malloc(sizeof(uint64_t) * n);
|
||||
memcpy(mag, b->limbs, sizeof(uint64_t) * n);
|
||||
/* Generous size: log10(2) ≈ 0.30103 → n*64*0.31 + sign + null. */
|
||||
size_t cap = (size_t)n * 21 + 4;
|
||||
char *buf = (char *)ul_malloc(cap);
|
||||
size_t pos = cap;
|
||||
buf[--pos] = '\0';
|
||||
const uint64_t CHUNK = 1000000000000000000ULL; /* 10^18 fits in u64. */
|
||||
while (n > 0) {
|
||||
uint64_t rem = mag_divmod_u64(mag, n, CHUNK, mag, &n);
|
||||
if (n > 0) {
|
||||
/* Full 18-digit chunk (zero-pad). */
|
||||
for (int i = 0; i < 18; i++) {
|
||||
buf[--pos] = '0' + (int)(rem % 10);
|
||||
rem /= 10;
|
||||
}
|
||||
} else {
|
||||
/* Top chunk: no leading zeros. */
|
||||
if (rem == 0) buf[--pos] = '0';
|
||||
else while (rem) { buf[--pos] = '0' + (int)(rem % 10); rem /= 10; }
|
||||
}
|
||||
}
|
||||
if (b->sign < 0) buf[--pos] = '-';
|
||||
/* Move string to start. */
|
||||
size_t len = cap - pos;
|
||||
memmove(buf, buf + pos, len);
|
||||
ul_free(mag);
|
||||
return buf;
|
||||
}
|
||||
65
c/builtins.c
65
c/builtins.c
|
|
@ -113,11 +113,13 @@ static Value bi_atan(Value *a, int n, Env *e) {
|
|||
|
||||
static Value bi_quotient(Value *a, int n, Env *e) {
|
||||
(void)e; CHECK_ARITY("quotient", 2);
|
||||
if (IS_INTEGER(a[0]) && IS_INTEGER(a[1])) return big_quotient(a[0], a[1]);
|
||||
return VAL_INT(as_number_int(a[0]) / as_number_int(a[1]));
|
||||
}
|
||||
|
||||
static Value bi_remainder(Value *a, int n, Env *e) {
|
||||
(void)e; CHECK_ARITY("remainder", 2);
|
||||
if (IS_INTEGER(a[0]) && IS_INTEGER(a[1])) return big_remainder(a[0], a[1]);
|
||||
int64_t x = as_number_int(a[0]), y = as_number_int(a[1]);
|
||||
int64_t r = x % y;
|
||||
/* Remainder has sign of dividend */
|
||||
|
|
@ -126,6 +128,7 @@ static Value bi_remainder(Value *a, int n, Env *e) {
|
|||
|
||||
static Value bi_modulo(Value *a, int n, Env *e) {
|
||||
(void)e; CHECK_ARITY("modulo", 2);
|
||||
if (IS_INTEGER(a[0]) && IS_INTEGER(a[1])) return big_modulo(a[0], a[1]);
|
||||
int64_t x = as_number_int(a[0]), y = as_number_int(a[1]);
|
||||
int64_t r = x % y;
|
||||
if ((r > 0 && y < 0) || (r < 0 && y > 0)) r += y;
|
||||
|
|
@ -134,18 +137,25 @@ static Value bi_modulo(Value *a, int n, Env *e) {
|
|||
|
||||
static Value bi_expt(Value *a, int n, Env *e) {
|
||||
(void)e; CHECK_ARITY("expt", 2);
|
||||
if (IS_INT(a[0]) && IS_INT(a[1]) && as_int(a[1]) >= 0) {
|
||||
int64_t base = as_int(a[0]), exp = as_int(a[1]);
|
||||
int64_t result = 1;
|
||||
for (int64_t i = 0; i < exp; i++) result *= base;
|
||||
return VAL_INT(result);
|
||||
/* Integer-base / non-negative integer-exp → bignum-aware path.
|
||||
* This catches (expt 2 256) and friends — fixnum * fixnum would
|
||||
* silently overflow inside the old int64 loop. */
|
||||
if (IS_INTEGER(a[0]) && IS_INTEGER(a[1])) {
|
||||
bool nonneg = (IS_INT(a[1]) && as_int(a[1]) >= 0) ||
|
||||
(IS_BIGNUM(a[1]) && AS_BIGNUM(a[1])->sign >= 0);
|
||||
if (nonneg) return big_expt(a[0], a[1]);
|
||||
}
|
||||
return make_double(pow(as_number_double(a[0]), as_number_double(a[1])));
|
||||
}
|
||||
|
||||
static Value bi_abs(Value *a, int n, Env *e) {
|
||||
(void)e; CHECK_ARITY("abs", 1);
|
||||
if (IS_INT(a[0])) { int64_t v = as_int(a[0]); return VAL_INT(v < 0 ? -v : v); }
|
||||
if (IS_INT(a[0])) {
|
||||
int64_t v = as_int(a[0]);
|
||||
if (v == FIXNUM_MIN) return big_abs(a[0]); /* avoid -INT_MIN overflow */
|
||||
return VAL_INT(v < 0 ? -v : v);
|
||||
}
|
||||
if (IS_BIGNUM(a[0])) return big_abs(a[0]);
|
||||
return make_double(fabs(as_number_double(a[0])));
|
||||
}
|
||||
|
||||
|
|
@ -231,6 +241,15 @@ static Value bi_number_to_string(Value *a, int n, Env *e) {
|
|||
char buf[128];
|
||||
if (n > 1) {
|
||||
int base = (int)as_number_int(a[1]);
|
||||
/* Bignum decimal path uses our renderer; other bases fall through to
|
||||
* int64 — at our scale (256-bit secp256k1) base-10 covers every site
|
||||
* that actually needs precision. */
|
||||
if (IS_BIGNUM(a[0]) && base == 10) {
|
||||
char *s = bignum_to_str(AS_BIGNUM(a[0]));
|
||||
Value r = make_string_from_cstr(s);
|
||||
ul_free(s);
|
||||
return r;
|
||||
}
|
||||
int64_t val = as_number_int(a[0]);
|
||||
switch (base) {
|
||||
case 2: {
|
||||
|
|
@ -280,8 +299,8 @@ NUM_CMP(num_ge, num_ge)
|
|||
NUM_PRED(zero_p, is_number(a[0]) && num_eq(a[0], VAL_INT(0)))
|
||||
NUM_PRED(positive_p, is_number(a[0]) && num_gt(a[0], VAL_INT(0)))
|
||||
NUM_PRED(negative_p, is_number(a[0]) && num_lt(a[0], VAL_INT(0)))
|
||||
NUM_PRED(odd_p, IS_INT(a[0]) && (as_int(a[0]) % 2 != 0))
|
||||
NUM_PRED(even_p, IS_INT(a[0]) && (as_int(a[0]) % 2 == 0))
|
||||
NUM_PRED(odd_p, IS_INTEGER(a[0]) && big_is_odd(a[0]))
|
||||
NUM_PRED(even_p, IS_INTEGER(a[0]) && !big_is_odd(a[0]))
|
||||
|
||||
static Value bi_nan_p(Value *a, int n, Env *e) {
|
||||
(void)e; CHECK_ARITY("nan?", 1);
|
||||
|
|
@ -357,10 +376,10 @@ static Value bi_equal_p(Value *a, int n, Env *e) {
|
|||
}
|
||||
|
||||
TYPE_PRED(number_p, is_number(a[0]))
|
||||
TYPE_PRED(integer_p, IS_INT(a[0]) || (IS_DOUBLE(a[0]) && as_double(a[0]) == floor(as_double(a[0]))))
|
||||
TYPE_PRED(integer_p, IS_INT(a[0]) || IS_BIGNUM(a[0]) || (IS_DOUBLE(a[0]) && as_double(a[0]) == floor(as_double(a[0]))))
|
||||
TYPE_PRED(real_p, is_number(a[0]))
|
||||
TYPE_PRED(rational_p, IS_INT(a[0]) || IS_RATIONAL(a[0]) || (IS_DOUBLE(a[0]) && isfinite(as_double(a[0]))))
|
||||
TYPE_PRED(exact_p, IS_INT(a[0]) || IS_RATIONAL(a[0]))
|
||||
TYPE_PRED(rational_p, IS_INT(a[0]) || IS_BIGNUM(a[0]) || IS_RATIONAL(a[0]) || (IS_DOUBLE(a[0]) && isfinite(as_double(a[0]))))
|
||||
TYPE_PRED(exact_p, IS_INT(a[0]) || IS_BIGNUM(a[0]) || IS_RATIONAL(a[0]))
|
||||
TYPE_PRED(inexact_p, IS_DOUBLE(a[0]))
|
||||
TYPE_PRED(pair_p, IS_PAIR(a[0]))
|
||||
TYPE_PRED(null_p, IS_NIL(a[0]))
|
||||
|
|
@ -796,6 +815,30 @@ static Value bi_string_to_number(Value *a, int n, Env *e) {
|
|||
(void)e; CHECK_MIN_ARITY("string->number", 1); check_string(a[0]);
|
||||
const char *s = AS_STRING(a[0])->data;
|
||||
int base = n > 1 ? (int)as_number_int(a[1]) : 10;
|
||||
/* Detect all-digit (optionally signed) shape → integer (fixnum or
|
||||
* bignum on overflow). */
|
||||
{
|
||||
const char *p = s;
|
||||
if (*p == '-' || *p == '+') p++;
|
||||
bool all_digits = (*p != '\0');
|
||||
for (const char *q = p; *q; q++) {
|
||||
char c = *q;
|
||||
int dval;
|
||||
if (c >= '0' && c <= '9') dval = c - '0';
|
||||
else if (base > 10 && c >= 'a' && c <= 'z') dval = c - 'a' + 10;
|
||||
else if (base > 10 && c >= 'A' && c <= 'Z') dval = c - 'A' + 10;
|
||||
else { all_digits = false; break; }
|
||||
if (dval >= base) { all_digits = false; break; }
|
||||
}
|
||||
if (all_digits) {
|
||||
char *end;
|
||||
errno = 0;
|
||||
long long val = strtoll(s, &end, base);
|
||||
if (*end == '\0' && errno == 0 && FITS_FIXNUM(val))
|
||||
return VAL_INT(val);
|
||||
return make_bignum_from_str(s, base);
|
||||
}
|
||||
}
|
||||
char *end;
|
||||
errno = 0;
|
||||
long long val = strtoll(s, &end, base);
|
||||
|
|
|
|||
6
c/eval.c
6
c/eval.c
|
|
@ -480,7 +480,7 @@ static bool is_literal(SyntaxTransformer *st, const char *name) {
|
|||
static bool sr_match(SyntaxTransformer *st, Value pat, Value form, Env *bindings) {
|
||||
if (IS_NIL(pat)) return IS_NIL(form);
|
||||
if (pat == VAL_TRUE || pat == VAL_FALSE) return pat == form;
|
||||
if (IS_INT(pat) || IS_DOUBLE(pat)) return values_equal(pat, form);
|
||||
if (IS_INT(pat) || IS_DOUBLE(pat) || IS_BIGNUM(pat)) return values_equal(pat, form);
|
||||
|
||||
if (IS_SYM(pat)) {
|
||||
const char *pname = sym_name(pat);
|
||||
|
|
@ -573,7 +573,7 @@ static bool sr_match(SyntaxTransformer *st, Value pat, Value form, Env *bindings
|
|||
|
||||
static Value sr_expand(SyntaxTransformer *st, Value tmpl, Env *bindings) {
|
||||
if (IS_NIL(tmpl) || tmpl == VAL_TRUE || tmpl == VAL_FALSE) return tmpl;
|
||||
if (IS_INT(tmpl) || IS_DOUBLE(tmpl)) return tmpl;
|
||||
if (IS_INT(tmpl) || IS_DOUBLE(tmpl) || IS_BIGNUM(tmpl)) return tmpl;
|
||||
if (IS_STRING(tmpl)) return tmpl;
|
||||
|
||||
if (IS_SYM(tmpl)) {
|
||||
|
|
@ -669,7 +669,7 @@ Value leval(Value expr, Env *env) {
|
|||
/* Self-evaluating */
|
||||
if (IS_NIL(expr) || IS_VOID(expr) || IS_TRUE(expr) || IS_FALSE(expr) ||
|
||||
IS_EOF(expr) || IS_CHAR(expr)) return expr;
|
||||
if (IS_INT(expr) || IS_DOUBLE(expr) || IS_RATIONAL(expr)) return expr;
|
||||
if (IS_INT(expr) || IS_DOUBLE(expr) || IS_RATIONAL(expr) || IS_BIGNUM(expr)) return expr;
|
||||
if (IS_STRING(expr)) return expr;
|
||||
if (IS_VECTOR(expr)) return expr;
|
||||
if (IS_BUILTIN(expr)) return expr;
|
||||
|
|
|
|||
77
c/lumbda.h
77
c/lumbda.h
|
|
@ -62,6 +62,10 @@ static inline char *ul_strdup(const char *s) {
|
|||
* Tag 3 = special (NIL, VOID, TRUE, FALSE, EOF, char in payload)
|
||||
* Tag 4 = builtin function pointer
|
||||
* Tag 5 = rational (pointer to Rational struct)
|
||||
* Tag 6 = bignum (pointer to Bignum struct) — arbitrary-precision integer.
|
||||
* Used for any integer outside the 48-bit fixnum range. Arithmetic
|
||||
* promotes fixnum to bignum on overflow & demotes bignum back to
|
||||
* fixnum when result fits, so callers never need to distinguish.
|
||||
*
|
||||
* Plain doubles (no NaN payload) are floating-point numbers.
|
||||
* ═══════════════════════════════════════════════════════════════════════════ */
|
||||
|
|
@ -81,6 +85,12 @@ typedef uint64_t Value;
|
|||
#define TAG_SPECIAL 3ULL
|
||||
#define TAG_BUILTIN 4ULL
|
||||
#define TAG_RATIONAL 5ULL
|
||||
#define TAG_BIGNUM 6ULL
|
||||
|
||||
/* Fixnum range: signed 48-bit (low 48 bits of payload). */
|
||||
#define FIXNUM_MIN ((int64_t)(-(1LL << 47)))
|
||||
#define FIXNUM_MAX ((int64_t)((1LL << 47) - 1))
|
||||
#define FITS_FIXNUM(v) ((v) >= FIXNUM_MIN && (v) <= FIXNUM_MAX)
|
||||
|
||||
/* Construct a NaN-boxed value */
|
||||
#define NANBOX(tag, payload) (QNAN | ((uint64_t)(tag) << TAG_SHIFT) | ((uint64_t)(payload) & PAYLOAD_MASK))
|
||||
|
|
@ -137,6 +147,9 @@ static inline int64_t as_int(Value v) {
|
|||
#define IS_SPECIAL(v) (!IS_DOUBLE(v) && GET_TAG(v) == TAG_SPECIAL)
|
||||
#define IS_BUILTIN(v) (!IS_DOUBLE(v) && GET_TAG(v) == TAG_BUILTIN)
|
||||
#define IS_RATIONAL(v) (!IS_DOUBLE(v) && GET_TAG(v) == TAG_RATIONAL)
|
||||
#define IS_BIGNUM(v) (!IS_DOUBLE(v) && GET_TAG(v) == TAG_BIGNUM)
|
||||
/* Integer-shaped: either fixnum or bignum. */
|
||||
#define IS_INTEGER(v) (IS_INT(v) || IS_BIGNUM(v))
|
||||
|
||||
#define IS_NIL(v) ((v) == VAL_NIL)
|
||||
#define IS_VOID(v) ((v) == VAL_VOID)
|
||||
|
|
@ -170,6 +183,7 @@ typedef enum {
|
|||
OBJ_ERROR,
|
||||
OBJ_SYNTAX_TRANSFORMER,
|
||||
OBJ_RATIONAL,
|
||||
OBJ_BIGNUM,
|
||||
} ObjType;
|
||||
|
||||
typedef struct ObjHeader {
|
||||
|
|
@ -288,6 +302,65 @@ Value make_rational(int64_t num, int64_t den);
|
|||
/* Returns integer Value if den==1, otherwise rational */
|
||||
Value rational_normalize(int64_t num, int64_t den);
|
||||
|
||||
/* ── Bignum — arbitrary-precision integer ────────────────────────────────── */
|
||||
|
||||
/*
|
||||
* Bignum representation: sign-magnitude.
|
||||
* sign: -1 (negative), 0 (zero — must have n_limbs == 0), +1 (positive)
|
||||
* n_limbs: number of u64 magnitude limbs (always trimmed of leading zeros)
|
||||
* limbs[i]: little-endian magnitude (limb 0 = least-significant 64 bits)
|
||||
*
|
||||
* Boehm GC owns every Bignum struct & its limbs array. Construction goes
|
||||
* through bignum_normalize() which trims leading-zero limbs & forces zero
|
||||
* to canonical form (sign=0, n_limbs=0).
|
||||
*
|
||||
* Result-producing helpers ALWAYS return a Value: if the magnitude fits in
|
||||
* 48 signed bits, a fixnum Value is returned; otherwise a bignum Value.
|
||||
* Callers therefore never branch on which one they get back.
|
||||
*/
|
||||
typedef struct Bignum {
|
||||
ObjHeader hdr;
|
||||
int32_t sign; /* -1, 0, +1 */
|
||||
uint32_t n_limbs; /* magnitude limb count after trim */
|
||||
uint64_t *limbs; /* little-endian magnitude */
|
||||
} Bignum;
|
||||
|
||||
#define AS_BIGNUM(v) ((Bignum *)(uintptr_t)GET_PAYLOAD(v))
|
||||
|
||||
/* Constructors */
|
||||
Value make_bignum_from_i64(int64_t v);
|
||||
Value make_bignum_from_str(const char *s, int base); /* base 10 or 16 */
|
||||
Value bignum_from_limbs(int32_t sign, uint64_t *limbs, uint32_t n);
|
||||
|
||||
/* Demote a bignum to fixnum when magnitude fits 48 signed bits. Caller passes
|
||||
* a (sign, limbs, n_limbs) triple; result is fixnum Value or bignum Value. */
|
||||
Value bignum_normalize(int32_t sign, uint64_t *limbs, uint32_t n);
|
||||
|
||||
/* Promote a fixnum Value to a bignum Value (no-op for already-bignum). */
|
||||
Value to_bignum(Value v);
|
||||
|
||||
/* Decimal-string render (caller frees via ul_free). */
|
||||
char *bignum_to_str(Bignum *b);
|
||||
|
||||
/* Arithmetic — accept fixnum or bignum, return fixnum-or-bignum
|
||||
* (demoted when possible). Sign/magnitude internal. */
|
||||
Value big_add(Value a, Value b);
|
||||
Value big_sub(Value a, Value b);
|
||||
Value big_mul(Value a, Value b);
|
||||
Value big_neg(Value a);
|
||||
Value big_quotient(Value a, Value b);
|
||||
Value big_remainder(Value a, Value b); /* sign of dividend */
|
||||
Value big_modulo(Value a, Value b); /* sign of divisor */
|
||||
Value big_expt(Value base, Value exp); /* exp >= 0, integer */
|
||||
int big_cmp(Value a, Value b); /* -1, 0, +1 */
|
||||
bool big_is_zero(Value v);
|
||||
bool big_is_odd(Value v);
|
||||
Value big_abs(Value v);
|
||||
|
||||
/* Coerce any number Value to a fresh Bignum struct (caller owns) for ops
|
||||
* that want explicit bignum arguments. */
|
||||
Bignum *as_bignum(Value v);
|
||||
|
||||
/* ── Environment ─────────────────────────────────────────────────────────── */
|
||||
|
||||
typedef struct EnvBinding {
|
||||
|
|
@ -719,9 +792,9 @@ static inline bool is_callable(Value v) {
|
|||
return false;
|
||||
}
|
||||
|
||||
/* Number extraction — works for int, double, rational */
|
||||
/* Number extraction — works for int, double, rational, bignum */
|
||||
static inline bool is_number(Value v) {
|
||||
return IS_INT(v) || IS_DOUBLE(v) || IS_RATIONAL(v);
|
||||
return IS_INT(v) || IS_DOUBLE(v) || IS_RATIONAL(v) || IS_BIGNUM(v);
|
||||
}
|
||||
|
||||
double as_number_double(Value v); /* coerce any number to double */
|
||||
|
|
|
|||
14
c/printer.c
14
c/printer.c
|
|
@ -136,6 +136,13 @@ static void show_value(Value v, bool display, StringBuilder *sb) {
|
|||
return;
|
||||
}
|
||||
|
||||
if (IS_BIGNUM(v)) {
|
||||
char *s = bignum_to_str(AS_BIGNUM(v));
|
||||
sb_appendz(sb, s);
|
||||
ul_free(s);
|
||||
return;
|
||||
}
|
||||
|
||||
if (IS_BUILTIN(v)) {
|
||||
sb_appendz(sb, "#<builtin>");
|
||||
return;
|
||||
|
|
@ -243,6 +250,13 @@ static void show_value(Value v, bool display, StringBuilder *sb) {
|
|||
sb_appendz(sb, buf);
|
||||
break;
|
||||
}
|
||||
|
||||
case OBJ_BIGNUM: {
|
||||
char *s = bignum_to_str(AS_BIGNUM(v));
|
||||
sb_appendz(sb, s);
|
||||
ul_free(s);
|
||||
break;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
|
|
|
|||
23
c/reader.c
23
c/reader.c
|
|
@ -230,17 +230,22 @@ Value parse_atom(const char *tok) {
|
|||
return v;
|
||||
}
|
||||
|
||||
/* Integer */
|
||||
/* Integer (fixnum or bignum). Match optional sign + all-digits. */
|
||||
{
|
||||
char *end;
|
||||
errno = 0;
|
||||
long long val = strtoll(tok, &end, 10);
|
||||
if (*end == '\0' && errno == 0) {
|
||||
/* Check if fits in 48-bit int */
|
||||
if (val >= -(1LL << 47) && val < (1LL << 47))
|
||||
const char *p = tok;
|
||||
if (*p == '-' || *p == '+') p++;
|
||||
bool all_digits = (*p != '\0');
|
||||
for (const char *q = p; *q; q++) {
|
||||
if (!isdigit((unsigned char)*q)) { all_digits = false; break; }
|
||||
}
|
||||
if (all_digits) {
|
||||
char *end;
|
||||
errno = 0;
|
||||
long long val = strtoll(tok, &end, 10);
|
||||
if (*end == '\0' && errno == 0 && FITS_FIXNUM(val))
|
||||
return VAL_INT(val);
|
||||
else
|
||||
return make_double((double)val);
|
||||
/* Out of fixnum range (or overflow) → bignum. */
|
||||
return make_bignum_from_str(tok, 10);
|
||||
}
|
||||
}
|
||||
|
||||
|
|
|
|||
47
c/test.c
47
c/test.c
|
|
@ -135,6 +135,51 @@ TEST(rational_creation) {
|
|||
ASSERT_EQ_INT(as_int(r2), 2);
|
||||
}
|
||||
|
||||
TEST(bignum_basic) {
|
||||
/* (expt 2 48) — first value past fixnum cap → bignum. */
|
||||
char *got = run_show("(expt 2 48)");
|
||||
ASSERT_EQ_STR(got, "281474976710656");
|
||||
ul_free(got);
|
||||
|
||||
/* (expt 2 256) — 78-digit hex world. */
|
||||
got = run_show("(expt 2 256)");
|
||||
ASSERT_EQ_STR(got, "115792089237316195423570985008687907853269984665640564039457584007913129639936");
|
||||
ul_free(got);
|
||||
|
||||
/* secp256k1 prime — full triple-subtract. */
|
||||
got = run_show("(- (expt 2 256) (expt 2 32) 977)");
|
||||
ASSERT_EQ_STR(got, "115792089237316195423570985008687907853269984665640564039457584007908834671663");
|
||||
ul_free(got);
|
||||
|
||||
/* Demotion back to fixnum when result fits. */
|
||||
got = run_show("(- (expt 2 64) (expt 2 64))");
|
||||
ASSERT_EQ_STR(got, "0");
|
||||
ul_free(got);
|
||||
}
|
||||
|
||||
TEST(bignum_arith) {
|
||||
/* Modular arithmetic across 78-digit operands. */
|
||||
char *got = run_show(
|
||||
"(define p (- (expt 2 256) (expt 2 32) 977))"
|
||||
"(modulo (- 0 1) p)");
|
||||
/* p - 1 */
|
||||
ASSERT_EQ_STR(got, "115792089237316195423570985008687907853269984665640564039457584007908834671662");
|
||||
ul_free(got);
|
||||
|
||||
/* quotient + remainder identity for bignum scale. */
|
||||
got = run_show(
|
||||
"(define p (- (expt 2 256) (expt 2 32) 977))"
|
||||
"(let* ((a (* p 7)) (q (quotient a p)) (r (remainder a p)))"
|
||||
" (list q r))");
|
||||
ASSERT_EQ_STR(got, "(7 0)");
|
||||
ul_free(got);
|
||||
|
||||
/* Predicate roundtrip. */
|
||||
got = run_show("(list (integer? (expt 2 100)) (exact? (expt 2 100)) (odd? (- (expt 2 100) 1)))");
|
||||
ASSERT_EQ_STR(got, "(#t #t #t)");
|
||||
ul_free(got);
|
||||
}
|
||||
|
||||
TEST(environment) {
|
||||
Env *g = make_env(NULL);
|
||||
g->global = g;
|
||||
|
|
@ -886,6 +931,8 @@ int main(void) {
|
|||
run_test_string_creation();
|
||||
run_test_vector_creation();
|
||||
run_test_rational_creation();
|
||||
run_test_bignum_basic();
|
||||
run_test_bignum_arith();
|
||||
run_test_environment();
|
||||
|
||||
printf("\n[reader]\n");
|
||||
|
|
|
|||
69
c/types.c
69
c/types.c
|
|
@ -304,6 +304,13 @@ double as_number_double(Value v) {
|
|||
Rational *r = AS_RATIONAL(v);
|
||||
return (double)r->num / (double)r->den;
|
||||
}
|
||||
if (IS_BIGNUM(v)) {
|
||||
Bignum *b = AS_BIGNUM(v);
|
||||
double d = 0.0;
|
||||
for (int32_t i = (int32_t)b->n_limbs - 1; i >= 0; i--)
|
||||
d = d * 18446744073709551616.0 + (double)b->limbs[i];
|
||||
return b->sign < 0 ? -d : d;
|
||||
}
|
||||
lisp_error("not a number");
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -315,6 +322,12 @@ int64_t as_number_int(Value v) {
|
|||
Rational *r = AS_RATIONAL(v);
|
||||
return r->num / r->den;
|
||||
}
|
||||
if (IS_BIGNUM(v)) {
|
||||
Bignum *b = AS_BIGNUM(v);
|
||||
uint64_t lo = b->n_limbs > 0 ? b->limbs[0] : 0;
|
||||
int64_t s = (int64_t)lo;
|
||||
return b->sign < 0 ? -s : s;
|
||||
}
|
||||
lisp_error("not a number");
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -327,11 +340,37 @@ static void to_rational(Value v, int64_t *num, int64_t *den) {
|
|||
}
|
||||
|
||||
static bool is_exact(Value v) {
|
||||
return IS_INT(v) || IS_RATIONAL(v);
|
||||
return IS_INT(v) || IS_RATIONAL(v) || IS_BIGNUM(v);
|
||||
}
|
||||
|
||||
/* Fixnum add/sub/mul that promote on overflow. Result is fixnum or bignum. */
|
||||
static Value fixnum_add(int64_t a, int64_t b) {
|
||||
int64_t r;
|
||||
if (__builtin_add_overflow(a, b, &r) || !FITS_FIXNUM(r))
|
||||
return big_add(make_bignum_from_i64(a), make_bignum_from_i64(b));
|
||||
return VAL_INT(r);
|
||||
}
|
||||
static Value fixnum_sub(int64_t a, int64_t b) {
|
||||
int64_t r;
|
||||
if (__builtin_sub_overflow(a, b, &r) || !FITS_FIXNUM(r))
|
||||
return big_sub(make_bignum_from_i64(a), make_bignum_from_i64(b));
|
||||
return VAL_INT(r);
|
||||
}
|
||||
static Value fixnum_mul(int64_t a, int64_t b) {
|
||||
int64_t r;
|
||||
if (__builtin_mul_overflow(a, b, &r) || !FITS_FIXNUM(r))
|
||||
return big_mul(make_bignum_from_i64(a), make_bignum_from_i64(b));
|
||||
return VAL_INT(r);
|
||||
}
|
||||
|
||||
Value num_add(Value a, Value b) {
|
||||
/* Pure integer path (fixnum or bignum) — never produces a rational. */
|
||||
if (IS_INTEGER(a) && IS_INTEGER(b)) {
|
||||
if (IS_INT(a) && IS_INT(b)) return fixnum_add(as_int(a), as_int(b));
|
||||
return big_add(a, b);
|
||||
}
|
||||
if (is_exact(a) && is_exact(b)) {
|
||||
/* Rational mix — fall back to int64_t numerator/denominator path. */
|
||||
int64_t an, ad, bn, bd;
|
||||
to_rational(a, &an, &ad);
|
||||
to_rational(b, &bn, &bd);
|
||||
|
|
@ -341,6 +380,10 @@ Value num_add(Value a, Value b) {
|
|||
}
|
||||
|
||||
Value num_sub(Value a, Value b) {
|
||||
if (IS_INTEGER(a) && IS_INTEGER(b)) {
|
||||
if (IS_INT(a) && IS_INT(b)) return fixnum_sub(as_int(a), as_int(b));
|
||||
return big_sub(a, b);
|
||||
}
|
||||
if (is_exact(a) && is_exact(b)) {
|
||||
int64_t an, ad, bn, bd;
|
||||
to_rational(a, &an, &ad);
|
||||
|
|
@ -351,6 +394,10 @@ Value num_sub(Value a, Value b) {
|
|||
}
|
||||
|
||||
Value num_mul(Value a, Value b) {
|
||||
if (IS_INTEGER(a) && IS_INTEGER(b)) {
|
||||
if (IS_INT(a) && IS_INT(b)) return fixnum_mul(as_int(a), as_int(b));
|
||||
return big_mul(a, b);
|
||||
}
|
||||
if (is_exact(a) && is_exact(b)) {
|
||||
int64_t an, ad, bn, bd;
|
||||
to_rational(a, &an, &ad);
|
||||
|
|
@ -361,6 +408,19 @@ Value num_mul(Value a, Value b) {
|
|||
}
|
||||
|
||||
Value num_div(Value a, Value b) {
|
||||
/* Integer / integer that divides cleanly stays integer (matches Python /
|
||||
* Scheme semantics for `/` between exacts when the quotient is exact).
|
||||
* For bignum case we treat as exact division (truncating to integer when
|
||||
* quotient is exact, else fall back to double for now — secp256k1 ops
|
||||
* never use fractional bignum). */
|
||||
if (IS_INTEGER(a) && IS_INTEGER(b)) {
|
||||
if (big_is_zero(b)) lisp_error("division by zero");
|
||||
Value q = big_quotient(a, b);
|
||||
Value r = big_remainder(a, b);
|
||||
if (big_is_zero(r)) return q;
|
||||
/* Inexact fallback. */
|
||||
return make_double(as_number_double(a) / as_number_double(b));
|
||||
}
|
||||
if (is_exact(a) && is_exact(b)) {
|
||||
int64_t an, ad, bn, bd;
|
||||
to_rational(a, &an, &ad);
|
||||
|
|
@ -374,7 +434,11 @@ Value num_div(Value a, Value b) {
|
|||
}
|
||||
|
||||
Value num_neg(Value a) {
|
||||
if (IS_INT(a)) return VAL_INT(-as_int(a));
|
||||
if (IS_INT(a)) {
|
||||
int64_t v = as_int(a);
|
||||
return fixnum_sub(0, v); /* handles -FIXNUM_MIN safely */
|
||||
}
|
||||
if (IS_BIGNUM(a)) return big_neg(a);
|
||||
if (IS_RATIONAL(a)) {
|
||||
Rational *r = AS_RATIONAL(a);
|
||||
return make_rational(-r->num, r->den);
|
||||
|
|
@ -383,6 +447,7 @@ Value num_neg(Value a) {
|
|||
}
|
||||
|
||||
static int num_cmp(Value a, Value b) {
|
||||
if (IS_INTEGER(a) && IS_INTEGER(b)) return big_cmp(a, b);
|
||||
if (is_exact(a) && is_exact(b)) {
|
||||
int64_t an, ad, bn, bd;
|
||||
to_rational(a, &an, &ad);
|
||||
|
|
|
|||
2
c/vm.c
2
c/vm.c
|
|
@ -148,7 +148,7 @@ void bc_compile(Value expr, CodeObj *code, Env *env, bool tail) {
|
|||
if (IS_NIL(expr) || IS_TRUE(expr) || IS_FALSE(expr) || IS_EOF(expr) || IS_CHAR(expr)) {
|
||||
code_emit(code, OP_CONST, expr); return;
|
||||
}
|
||||
if (IS_INT(expr) || IS_DOUBLE(expr) || IS_RATIONAL(expr)) {
|
||||
if (IS_INT(expr) || IS_DOUBLE(expr) || IS_RATIONAL(expr) || IS_BIGNUM(expr)) {
|
||||
code_emit(code, OP_CONST, expr); return;
|
||||
}
|
||||
if (IS_STRING(expr) || IS_VECTOR(expr)) {
|
||||
|
|
|
|||
146
tests/bignum-cross-tier.lsp
Normal file
146
tests/bignum-cross-tier.lsp
Normal file
|
|
@ -0,0 +1,146 @@
|
|||
;;; bignum-cross-tier.lsp — Cross-tier byte-identity for arbitrary-precision
|
||||
;;; integers. Runs identically under Python tier & C tier. Python tier is
|
||||
;;; our oracle (native int handles every magnitude); C tier joined the
|
||||
;;; bignum club here & must agree byte-for-byte with Python.
|
||||
;;;
|
||||
;;; Python: python3 lumbda.py --fast tests/bignum-cross-tier.lsp
|
||||
;;; C: ./c/lumbda tests/bignum-cross-tier.lsp
|
||||
;;;
|
||||
;;; Every assertion logs PASS or FAIL; a final summary lines tallies them.
|
||||
;;; Identity discipline: every PASS line on both tiers must match
|
||||
;;; byte-for-byte.
|
||||
|
||||
(define *pass* 0)
|
||||
(define *fail* 0)
|
||||
|
||||
(define (assert-equal name got expected)
|
||||
(if (equal? got expected)
|
||||
(begin (set! *pass* (+ *pass* 1))
|
||||
(display "PASS: ") (display name) (newline))
|
||||
(begin (set! *fail* (+ *fail* 1))
|
||||
(display "FAIL: ") (display name)
|
||||
(display " got=") (write got)
|
||||
(display " expected=") (write expected) (newline))))
|
||||
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
;;; Literals & basic reader/printer round-trip
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
|
||||
(assert-equal "fixnum-literal-1" 42 42)
|
||||
(assert-equal "fixnum-literal-2" -7 -7)
|
||||
(assert-equal "bignum-literal-pos"
|
||||
281474976710656 ; 2^48
|
||||
281474976710656)
|
||||
(assert-equal "bignum-literal-neg"
|
||||
-281474976710656
|
||||
-281474976710656)
|
||||
(assert-equal "bignum-literal-256"
|
||||
115792089237316195423570985008687907853269984665640564039457584007913129639936
|
||||
115792089237316195423570985008687907853269984665640564039457584007913129639936)
|
||||
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
;;; expt — promotion through fixnum boundary
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
|
||||
(assert-equal "expt-2-30" (expt 2 30) 1073741824)
|
||||
(assert-equal "expt-2-47" (expt 2 47) 140737488355328)
|
||||
(assert-equal "expt-2-48" (expt 2 48) 281474976710656)
|
||||
(assert-equal "expt-2-64" (expt 2 64) 18446744073709551616)
|
||||
(assert-equal "expt-2-128" (expt 2 128)
|
||||
340282366920938463463374607431768211456)
|
||||
(assert-equal "expt-2-256" (expt 2 256)
|
||||
115792089237316195423570985008687907853269984665640564039457584007913129639936)
|
||||
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
;;; secp256k1 prime — 78-digit triple-subtraction lands precisely
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
|
||||
(define secp-p (- (expt 2 256) (expt 2 32) 977))
|
||||
(assert-equal "secp256k1-prime"
|
||||
secp-p
|
||||
115792089237316195423570985008687907853269984665640564039457584007908834671663)
|
||||
(assert-equal "secp-p length" (string-length (number->string secp-p)) 78)
|
||||
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
;;; Arithmetic at bignum scale
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
|
||||
(assert-equal "add-bignum"
|
||||
(+ secp-p secp-p)
|
||||
231584178474632390847141970017375815706539969331281128078915168015817669343326)
|
||||
|
||||
(assert-equal "sub-bignum"
|
||||
(- (* secp-p 2) secp-p)
|
||||
secp-p)
|
||||
|
||||
;;; (* secp-p secp-p) — magnitude check via number->string length; exact
|
||||
;;; value differs across tier baselines only if their bignum products
|
||||
;;; diverge — Python tier serves as oracle (literal computed there).
|
||||
(assert-equal "mul-bignum"
|
||||
(* secp-p secp-p)
|
||||
13407807929942597099574024998205846127479365820592393377723561443720769383374469661147847687812952081302854773939601805382211292725060150247698793015185569)
|
||||
|
||||
;;; Multi-step modular: a * b mod p
|
||||
(define a 12345678901234567890123456789012345678901234567890)
|
||||
(define b 98765432109876543210987654321098765432109876543210)
|
||||
(assert-equal "mod-bignum"
|
||||
(modulo (* a b) secp-p)
|
||||
(modulo (* a b) secp-p))
|
||||
|
||||
(assert-equal "quotient-bignum"
|
||||
(quotient (* secp-p 7) secp-p) 7)
|
||||
|
||||
(assert-equal "remainder-bignum"
|
||||
(remainder (+ (* secp-p 3) 42) secp-p) 42)
|
||||
|
||||
;;; modulo wraps negative dividend by divisor sign
|
||||
(assert-equal "modulo-neg-bignum"
|
||||
(modulo (- 0 1) secp-p)
|
||||
(- secp-p 1))
|
||||
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
;;; Comparisons across fixnum/bignum boundary
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
|
||||
(assert-equal "lt-bignum-fixnum" (< 1 (expt 2 100)) #t)
|
||||
(assert-equal "gt-bignum-fixnum" (> (expt 2 100) 999) #t)
|
||||
(assert-equal "eq-bignum-bignum" (= (expt 2 64) 18446744073709551616) #t)
|
||||
(assert-equal "neq-bignum-fixnum" (= (expt 2 50) 0) #f)
|
||||
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
;;; Predicates
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
|
||||
(assert-equal "integer?-bignum" (integer? (expt 2 100)) #t)
|
||||
(assert-equal "exact?-bignum" (exact? (expt 2 100)) #t)
|
||||
(assert-equal "odd?-bignum" (odd? (- (expt 2 100) 1)) #t)
|
||||
(assert-equal "even?-bignum" (even? (expt 2 100)) #t)
|
||||
(assert-equal "zero?-bignum" (zero? (- (expt 2 100) (expt 2 100))) #t)
|
||||
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
;;; Promotion & demotion
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
|
||||
;;; Result of (- bignum bignum) that fits a fixnum demotes back.
|
||||
(assert-equal "demote" (- (expt 2 64) (expt 2 64)) 0)
|
||||
(assert-equal "demote-small" (- (expt 2 50) (- (expt 2 50) 7)) 7)
|
||||
|
||||
;;; number->string round-trip
|
||||
(assert-equal "n->s-bignum"
|
||||
(number->string (expt 2 100))
|
||||
"1267650600228229401496703205376")
|
||||
|
||||
;;; string->number round-trip
|
||||
(assert-equal "s->n-bignum"
|
||||
(string->number "1267650600228229401496703205376")
|
||||
(expt 2 100))
|
||||
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
;;; Summary
|
||||
;;; ═══════════════════════════════════════════════════════════════
|
||||
|
||||
(display "════════════════════════════════════════") (newline)
|
||||
(display "bignum cross-tier: ")
|
||||
(display *pass*) (display " passed, ")
|
||||
(display *fail*) (display " failed")
|
||||
(newline)
|
||||
Loading…
Add table
Add a link
Reference in a new issue