;; -*- Mode: Lisp; -*- ;; --------------------------------------------------------------------------- ;; Title: Most Simple Bignum Implementation in Most Simple Lisp ;; Created: 2025-02-08 ;; Author: Gilbert Baumann ;; License: MIT style (see below) ;; --------------------------------------------------------------------------- ;; (c) copyright 2025 by Gilbert Baumann ;; Permission is hereby granted, free of charge, to any person obtaining ;; a copy of this software and associated documentation files (the ;; "Software"), to deal in the Software without restriction, including ;; without limitation the rights to use, copy, modify, merge, publish, ;; distribute, sublicense, and/or sell copies of the Software, and to ;; permit persons to whom the Software is furnished to do so, subject to ;; the following conditions: ;; ;; The above copyright notice and this permission notice shall be ;; included in all copies or substantial portions of the Software. ;; ;; THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, ;; EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF ;; MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. ;; IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY ;; CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, ;; TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE ;; SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. ;;;; -- Purpose ------------------------------------------------------------------------------- ;; In case you find yourself stranded in like the 80s left with an itty ;; machine and an even more itty LISP that doesn't even bother to implement ;; bignums, but at least has fixnums, this little file might be helpful. ;; We implement bignum addition, subtraction, multiplication, and division in ;; the most straight-forward way. Bignums are represented in B's complement ;; representation as lists (sic!). LSB digit first. ;; Bignums must be normalized. The sign of a bignum is negative if the most ;; significant digit is not less than half the the base. Another 0 or B-1 is ;; added as needed. Further a bignum may not be longer than needed, but must ;; feature at least one digit. (That is 0 is (0) not NIL). ;; Suppose B = 10: ;; 0 = (0) ;; 123 = (3 2 1) ;; 567 = (7 6 5 0) extra 0 to make it positive. ;; -1 = (9) ;; -99 = (1 0 9) ;; -433 = (7 6 5) ;;;; ------------------------------------------------------------------------------------------ (defpackage :bignum (:use :cl)) (in-package :bignum) (defparameter +base+ ;; Doesn't matter. For debugging 10 is a good choice. But it better is even! (ash 1 (floor (integer-length most-positive-fixnum) 2))) (define-symbol-macro +halfbase+ (floor +base+ 2)) (defparameter +z-zero+ '(0) "A bignum zero") (defun z-cons (new yet) "Pushes a new digit, and ensures the result is still normalized." (if (and yet (null (cdr yet)) (= (car yet) (digit-sign new))) (list new) (cons new yet))) (defun z-add (a b &optional (c 0) &aux (sum 0) a-digit b-digit res) "Adds two integers in B's complement representation" (loop (setq a-digit (if a (car a) (digit-sign a-digit))) (setq b-digit (if b (car b) (digit-sign b-digit))) (setf (values c sum) (truncate (+ c a-digit b-digit) +base+)) (push sum res) (when (and (null a) (null b)) (return)) (pop a) (pop b)) (norm-and-reverse res)) (defun z-sub (a b &optional (c 1) &aux (sum 0) a-digit b-digit res) "Substracts two integers in B's complement representation" (loop (setq a-digit (if a (car a) (digit-sign a-digit))) (setq b-digit (if b (car b) (digit-sign b-digit))) (setf (values c sum) (truncate (+ a-digit (digit-neg b-digit) c) +base+)) (push sum res) (when (and (null a) (null b)) (return)) (pop a) (pop b)) (norm-and-reverse res)) (defun norm-and-reverse (res) "Return the reverse and normalized integer from /res/." (do () ((not (and (cdr res) (= (car res) (digit-sign (cadr res)))))) (pop res)) (nreverse res)) (defun digit-sign (x) "Sign of a digit. 0 if positive, B/2 if negative." (if (< x +halfbase+) 0 (1- +base+))) (defun digit-neg (d) "Complement of a single digit" (- (1- +base+) d)) ;;; Complement (defun magnitude-sign (x) "-> |x| ; minusp" (if (z-minusp x) (values (z-neg x) t) (values x nil))) (defun z-neg (b &aux res b-digit (c 1) sum) "Complement of an integer" (loop (setq b-digit (if b (car b) (digit-sign b-digit))) (setf (values c sum) (floor (+ (digit-neg b-digit) c) +base+)) (push sum res) (when (null b) (return)) (pop b)) (norm-and-reverse res)) (defun z-minusp (a) (do ((d (car a) (car a)) (a (cdr a) (cdr a))) ((null a) (>= d +halfbase+)))) ;;; Multiplication (defun z-mul (a b) "Product of two signed integers." (multiple-value-bind (am as) (magnitude-sign a) (multiple-value-bind (bm bs) (magnitude-sign b) (let ((p (n-mul am bm))) (if (if as bs (not bs)) p (z-neg p)))))) (defun n-mul (a b) "Product of two unsigned integers." (do ((acc +z-zero+ (z-add acc (z-scale-by-k**base s (z-scale (car k) b)))) (i 0 (+ i 1)) (k a (cdr k)) (s 0 (+ s 1))) ((null k) acc))) (defun z-scale-by-k**base (k z) "The number /z/ multiplied by B**k" (do ((k k (- k 1)) (z z (cons 0 z))) ((= k 0) z))) (defun z-scale (k a &optional (c 0) &aux (res nil) sum) "The number /a/ multiplied by the single digit /k/." (let (a-digit) (loop (setq a-digit (if a (car a) (digit-sign a-digit))) (setf (values c sum) (floor (+ (* k a-digit) c) +base+)) (push sum res) (unless a (return)) (pop a)) (norm-and-reverse res))) ;;; Division (defun z-div (a b) "Signed integer division" (multiple-value-bind (am as) (magnitude-sign a) (multiple-value-bind (bm bs) (magnitude-sign b) (multiple-value-bind (q r) (n-div am bm) (values (if (if as bs (not bs)) q (z-neg q)) (if as (z-neg r) r)))))) (defun n-div (a b) "Unsigned integer division" ;; We still need to scale. (let* ((na (n-length a)) (nb (n-length b)) (q +z-zero+)) (when (= nb 0) (error "Division by zero")) (when (>= na nb) (let* ((b* (ith-digit b (1- nb)))) (do ((i na (- i 1))) ((< i nb)) (do* ((guess (min (1- +base+) (floor (+ (* +base+ (ith-digit a i)) (ith-digit a (1- i))) b*)) (- guess 1)) (r (z-sub a (z-scale-by-k**base (- i nb) (z-scale guess b 0))) (z-add r (z-scale-by-k**base (- i nb) b)))) ((not (z-minusp r)) (setq a r) (setq q (z-cons guess q))))))) (values q a))) (defun n-length (a) "Length of unsigned integer in digits" (do ((r 0) (i 0 (+ i 1)) (q a (cdr q))) ((null q) r) (unless (zerop (car q)) (setq r (+ i 1))))) (defun ith-digit (x i) "Return the ith digit of /x/. Zero-based and LSB-first." (or (nth i x) 0)) ;;;; -- Tests --------------------------------------------------------------------------------- (defun integer-z (x) "Integer to bignum" (cond ((< x 0) (z-neg (integer-z (- x)))) ((= x 0) (list 0)) (t (multiple-value-bind (hi low) (floor x +base+) (z-cons low (integer-z hi)))))) (defun z-integer (x) "Integer from bignum" (if (z-minusp x) (- (z-integer (z-neg x))) (+ (car x) (if (cdr x) (* +base+ (z-integer (cdr x))) 0)))) (defun test-op (num-op z-op &optional (a-min -100) (a-max 100) (b-min a-min) (b-max a-max)) (do ((a a-min (+ a 1))) ((> a a-max)) (do ((b b-min (+ b 1))) ((> b b-max)) (let* ((r (ignore-errors (funcall num-op a b))) (az (integer-z a)) (bz (integer-z b)) (rz (ignore-errors (funcall z-op az bz))) (rz* (and r (integer-z r))) (r* (and rz (z-integer rz)))) (unless (and (eql r* r) (equal rz* rz)) (error "Fail for A = ~D, B = ~D, want ~D, have ~D" a b r r*)))))) (defun run-tests () (dolist (k '(100 1000)) (test-op #'+ #'z-add (- k) k (- k) k) (test-op #'- #'z-sub (- k) k (- k) k) (test-op #'* #'z-mul (- k) k (- k) k) (test-op (lambda (x y) (nth-value 0 (truncate x y))) (lambda (x y) (nth-value 0 (z-div x y))) (- k) k (- k) k) (test-op (lambda (x y) (nth-value 1 (truncate x y))) (lambda (x y) (nth-value 1 (z-div x y))) (- k) k (- k) k))) (defun random-tests (n min max &optional (bmin min) (bmax max)) (let ((ops (list (list #'+ #'z-add) (list #'- #'z-sub) (list #'* #'z-mul) (list #'truncate #'z-div) (list (lambda (x y) (nth-value 1 (truncate x y))) (lambda (x y) (nth-value 1 (z-div x y))))))) (dotimes (i n) (let ((a (+ min (random (- max min)))) (b (+ bmin (random (- bmax bmin)))) (op (elt ops (random (length ops))))) (test-op (car op) (cadr op) a a b b)))))