Mercurial > emacs

comparison lisp/calc/calc.el @ 81573:d5640ed7c397

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
(math-bignum-digit-length,math-bignum-digit-size,math-small-integer-size): New constants. (math-normalize,math-bignum-big,math-make-float,math-div10-bignum) (math-scale-left,math-scale-left-bignum,math-scale-right) (math-scale-right-bignum,math-scale-rounding,math-add,math-add-bignum) (math-sub-bignum,math-sub,math-mul,math-mul-bignum,math-mul-bignum-digit) (math-idivmod,math-quotient,math-div-bignum,math-div-bignum-digit) (math-div-bignum-part,math-format-bignum-decimal,math-read-bignum): Use math-bignum-digit-length, math-bignum-digit-size and math-small-integer-size.
author Jay Belanger <jay.p.belanger@gmail.com>
date Sat, 23 Jun 2007 04:05:29 +0000
parents 041563cafc04
children 38726b7738d9
comparison
equal deleted inserted replaced
81572:0991efe3cafa 81573:d5640ed7c397
2281 (calcDigit-nondigit)))) 2281 (calcDigit-nondigit))))
2282 2282
2283 2283
2284 2284
2285 2285
2286 2286 (defconst math-bignum-digit-length 3
2287 "The length of a \"digit\" in Calc bignums.
2288 If a big integer is of the form (bigpos N0 N1 ...), this is the
2289 length of the allowable Emacs integers N0, N1,...
2290 The value of 2*10^(2*MATH-BIGNUM-DIGIT-LENGTH) must be less than the
2291 largest Emacs integer.")
2292
2293 (defconst math-bignum-digit-size (expt 10 math-bignum-digit-length)
2294 "An upper bound for the size of the \"digit\"s in Calc bignums.")
2295
2296 (defconst math-small-integer-size (expt 10 (* 2 math-bignum-digit-length))
2297 "An upper bound for the size of \"small integer\"s in Calc.")
2287 2298
2288 2299
2289 ;;;; Arithmetic routines. 2300 ;;;; Arithmetic routines.
2290 ;;; 2301 ;;;
2291 ;;; An object as manipulated by one of these routines may take any of the 2302 ;;; An object as manipulated by one of these routines may take any of the
2292 ;;; following forms: 2303 ;;; following forms:
2293 ;;; 2304 ;;;
2294 ;;; integer An integer. For normalized numbers, this format 2305 ;;; integer An integer. For normalized numbers, this format
2295 ;;; is used only for -999999 ... 999999. 2306 ;;; is used only for
2307 ;;; negative math-small-integer-size + 1 to
2308 ;;; math-small-integer-size - 1
2296 ;;; 2309 ;;;
2297 ;;; (bigpos N0 N1 N2 ...) A big positive integer, N0 + N1*1000 + N2*10^6 ... 2310 ;;; (bigpos N0 N1 N2 ...) A big positive integer,
2298 ;;; (bigneg N0 N1 N2 ...) A big negative integer, - N0 - N1*1000 ... 2311 ;;; N0 + N1*math-bignum-digit-size
2299 ;;; Each digit N is in the range 0 ... 999. 2312 ;;; + N2*(math-bignum-digit-size)^2 ...
2313 ;;; (bigneg N0 N1 N2 ...) A big negative integer,
2314 ;;; - N0 - N1*math-bignum-digit-size ...
2315 ;;; Each digit N is in the range
2316 ;;; 0 ... math-bignum-digit-size -1.
2300 ;;; Normalized, always at least three N present, 2317 ;;; Normalized, always at least three N present,
2301 ;;; and the most significant N is nonzero. 2318 ;;; and the most significant N is nonzero.
2302 ;;; 2319 ;;;
2303 ;;; (frac NUM DEN) A fraction. NUM and DEN are small or big integers. 2320 ;;; (frac NUM DEN) A fraction. NUM and DEN are small or big integers.
2304 ;;; Normalized, DEN > 1. 2321 ;;; Normalized, DEN > 1.
2384 (defvar math-normalize-a) 2401 (defvar math-normalize-a)
2385 (defun math-normalize (math-normalize-a) 2402 (defun math-normalize (math-normalize-a)
2386 (cond 2403 (cond
2387 ((not (consp math-normalize-a)) 2404 ((not (consp math-normalize-a))
2388 (if (integerp math-normalize-a) 2405 (if (integerp math-normalize-a)
2389 (if (or (>= math-normalize-a 1000000) (<= math-normalize-a -1000000)) 2406 (if (or (>= math-normalize-a math-small-integer-size)
2407 (<= math-normalize-a (- math-small-integer-size)))
2390 (math-bignum math-normalize-a) 2408 (math-bignum math-normalize-a)
2391 math-normalize-a) 2409 math-normalize-a)
2392 math-normalize-a)) 2410 math-normalize-a))
2393 ((eq (car math-normalize-a) 'bigpos) 2411 ((eq (car math-normalize-a) 'bigpos)
2394 (if (eq (nth (1- (length math-normalize-a)) math-normalize-a) 0) 2412 (if (eq (nth (1- (length math-normalize-a)) math-normalize-a) 0)
2399 (setcdr last nil))) 2417 (setcdr last nil)))
2400 (if (cdr (cdr (cdr math-normalize-a))) 2418 (if (cdr (cdr (cdr math-normalize-a)))
2401 math-normalize-a 2419 math-normalize-a
2402 (cond 2420 (cond
2403 ((cdr (cdr math-normalize-a)) (+ (nth 1 math-normalize-a) 2421 ((cdr (cdr math-normalize-a)) (+ (nth 1 math-normalize-a)
2404 (* (nth 2 math-normalize-a) 1000))) 2422 (* (nth 2 math-normalize-a)
2423 math-bignum-digit-size)))
2405 ((cdr math-normalize-a) (nth 1 math-normalize-a)) 2424 ((cdr math-normalize-a) (nth 1 math-normalize-a))
2406 (t 0)))) 2425 (t 0))))
2407 ((eq (car math-normalize-a) 'bigneg) 2426 ((eq (car math-normalize-a) 'bigneg)
2408 (if (eq (nth (1- (length math-normalize-a)) math-normalize-a) 0) 2427 (if (eq (nth (1- (length math-normalize-a)) math-normalize-a) 0)
2409 (let* ((last (setq math-normalize-a (copy-sequence math-normalize-a))) 2428 (let* ((last (setq math-normalize-a (copy-sequence math-normalize-a)))
2413 (setcdr last nil))) 2432 (setcdr last nil)))
2414 (if (cdr (cdr (cdr math-normalize-a))) 2433 (if (cdr (cdr (cdr math-normalize-a)))
2415 math-normalize-a 2434 math-normalize-a
2416 (cond 2435 (cond
2417 ((cdr (cdr math-normalize-a)) (- (+ (nth 1 math-normalize-a) 2436 ((cdr (cdr math-normalize-a)) (- (+ (nth 1 math-normalize-a)
2418 (* (nth 2 math-normalize-a) 1000)))) 2437 (* (nth 2 math-normalize-a)
2438 math-bignum-digit-size))))
2419 ((cdr math-normalize-a) (- (nth 1 math-normalize-a))) 2439 ((cdr math-normalize-a) (- (nth 1 math-normalize-a)))
2420 (t 0)))) 2440 (t 0))))
2421 ((eq (car math-normalize-a) 'float) 2441 ((eq (car math-normalize-a) 'float)
2422 (math-make-float (math-normalize (nth 1 math-normalize-a)) 2442 (math-make-float (math-normalize (nth 1 math-normalize-a))
2423 (nth 2 math-normalize-a))) 2443 (nth 2 math-normalize-a)))
2533 (cons 'bigneg (math-bignum-big (- a))))) 2553 (cons 'bigneg (math-bignum-big (- a)))))
2534 2554
2535 (defun math-bignum-big (a) ; [L s] 2555 (defun math-bignum-big (a) ; [L s]
2536 (if (= a 0) 2556 (if (= a 0)
2537 nil 2557 nil
2538 (cons (% a 1000) (math-bignum-big (/ a 1000))))) 2558 (cons (% a math-bignum-digit-size)
2559 (math-bignum-big (/ a math-bignum-digit-size)))))
2539 2560
2540 2561
2541 ;;; Build a normalized floating-point number. [F I S] 2562 ;;; Build a normalized floating-point number. [F I S]
2542 (defun math-make-float (mant exp) 2563 (defun math-make-float (mant exp)
2543 (if (eq mant 0) 2564 (if (eq mant 0)
2550 (let ((digs (cdr mant))) 2571 (let ((digs (cdr mant)))
2551 (if (= (% (car digs) 10) 0) 2572 (if (= (% (car digs) 10) 0)
2552 (progn 2573 (progn
2553 (while (= (car digs) 0) 2574 (while (= (car digs) 0)
2554 (setq digs (cdr digs) 2575 (setq digs (cdr digs)
2555 exp (+ exp 3))) 2576 exp (+ exp math-bignum-digit-length)))
2556 (while (= (% (car digs) 10) 0) 2577 (while (= (% (car digs) 10) 0)
2557 (setq digs (math-div10-bignum digs) 2578 (setq digs (math-div10-bignum digs)
2558 exp (1+ exp))) 2579 exp (1+ exp)))
2559 (setq mant (math-normalize (cons (car mant) digs)))))) 2580 (setq mant (math-normalize (cons (car mant) digs))))))
2560 (while (= (% mant 10) 0) 2581 (while (= (% mant 10) 0)
2568 (signal 'math-overflow nil) 2589 (signal 'math-overflow nil)
2569 (list 'float mant exp))))) 2590 (list 'float mant exp)))))
2570 2591
2571 (defun math-div10-bignum (a) ; [l l] 2592 (defun math-div10-bignum (a) ; [l l]
2572 (if (cdr a) 2593 (if (cdr a)
2573 (cons (+ (/ (car a) 10) (* (% (nth 1 a) 10) 100)) 2594 (cons (+ (/ (car a) 10) (* (% (nth 1 a) 10)
2595 (expt 10 (1- math-bignum-digit-length))))
2574 (math-div10-bignum (cdr a))) 2596 (math-div10-bignum (cdr a)))
2575 (list (/ (car a) 10)))) 2597 (list (/ (car a) 10))))
2576 2598
2577 ;;; Coerce A to be a float. [F N; V V] [Public] 2599 ;;; Coerce A to be a float. [F N; V V] [Public]
2578 (defun math-float (a) 2600 (defun math-float (a)
2599 (defun math-numdigs (a) 2621 (defun math-numdigs (a)
2600 (if (consp a) 2622 (if (consp a)
2601 (if (cdr a) 2623 (if (cdr a)
2602 (let* ((len (1- (length a))) 2624 (let* ((len (1- (length a)))
2603 (top (nth len a))) 2625 (top (nth len a)))
2604 (+ (* len 3) (cond ((>= top 100) 0) ((>= top 10) -1) (t -2)))) 2626 (+ (* (1- len) math-bignum-digit-length) (math-numdigs top)))
2605 0) 2627 0)
2606 (cond ((>= a 100) (+ (math-numdigs (/ a 1000)) 3)) 2628 (cond ((>= a 100) (+ (math-numdigs (/ a 1000)) 3))
2607 ((>= a 10) 2) 2629 ((>= a 10) 2)
2608 ((>= a 1) 1) 2630 ((>= a 1) 1)
2609 ((= a 0) 0) 2631 ((= a 0) 0)
2620 (defun math-scale-left (a n) ; [I I S] 2642 (defun math-scale-left (a n) ; [I I S]
2621 (if (= n 0) 2643 (if (= n 0)
2622 a 2644 a
2623 (if (consp a) 2645 (if (consp a)
2624 (cons (car a) (math-scale-left-bignum (cdr a) n)) 2646 (cons (car a) (math-scale-left-bignum (cdr a) n))
2625 (if (>= n 3) 2647 (if (>= n math-bignum-digit-length)
2626 (if (or (>= a 1000) (<= a -1000)) 2648 (if (or (>= a math-bignum-digit-size)
2649 (<= a (- math-bignum-digit-size)))
2627 (math-scale-left (math-bignum a) n) 2650 (math-scale-left (math-bignum a) n)
2628 (math-scale-left (* a 1000) (- n 3))) 2651 (math-scale-left (* a math-bignum-digit-size)
2629 (if (= n 2) 2652 (- n math-bignum-digit-length)))
2630 (if (or (>= a 10000) (<= a -10000)) 2653 (let ((sz (expt 10 (- (* 2 math-bignum-digit-length) n))))
2631 (math-scale-left (math-bignum a) 2) 2654 (if (or (>= a sz) (<= a (- sz)))
2632 (* a 100)) 2655 (math-scale-left (math-bignum a) n)
2633 (if (or (>= a 100000) (<= a -100000)) 2656 (* a (expt 10 n))))))))
2634 (math-scale-left (math-bignum a) 1)
2635 (* a 10)))))))
2636 2657
2637 (defun math-scale-left-bignum (a n) 2658 (defun math-scale-left-bignum (a n)
2638 (if (>= n 3) 2659 (if (>= n math-bignum-digit-length)
2639 (while (>= (setq a (cons 0 a) 2660 (while (>= (setq a (cons 0 a)
2640 n (- n 3)) 3))) 2661 n (- n math-bignum-digit-length))
2662 math-bignum-digit-length)))
2641 (if (> n 0) 2663 (if (> n 0)
2642 (math-mul-bignum-digit a (if (= n 2) 100 10) 0) 2664 (math-mul-bignum-digit a (expt 10 n) 0)
2643 a)) 2665 a))
2644 2666
2645 (defun math-scale-right (a n) ; [i i S] 2667 (defun math-scale-right (a n) ; [i i S]
2646 (if (= n 0) 2668 (if (= n 0)
2647 a 2669 a
2649 (cons (car a) (math-scale-right-bignum (cdr a) n)) 2671 (cons (car a) (math-scale-right-bignum (cdr a) n))
2650 (if (<= a 0) 2672 (if (<= a 0)
2651 (if (= a 0) 2673 (if (= a 0)
2652 0 2674 0
2653 (- (math-scale-right (- a) n))) 2675 (- (math-scale-right (- a) n)))
2654 (if (>= n 3) 2676 (if (>= n math-bignum-digit-length)
2655 (while (and (> (setq a (/ a 1000)) 0) 2677 (while (and (> (setq a (/ a math-bignum-digit-size)) 0)
2656 (>= (setq n (- n 3)) 3)))) 2678 (>= (setq n (- n math-bignum-digit-length))
2657 (if (= n 2) 2679 math-bignum-digit-length))))
2658 (/ a 100) 2680 (if (> n 0)
2659 (if (= n 1) 2681 (/ a (expt 10 n))
2660 (/ a 10) 2682 a)))))
2661 a))))))
2662 2683
2663 (defun math-scale-right-bignum (a n) ; [L L S; l l S] 2684 (defun math-scale-right-bignum (a n) ; [L L S; l l S]
2664 (if (>= n 3) 2685 (if (>= n math-bignum-digit-length)
2665 (setq a (nthcdr (/ n 3) a) 2686 (setq a (nthcdr (/ n math-bignum-digit-length) a)
2666 n (% n 3))) 2687 n (% n math-bignum-digit-length)))
2667 (if (> n 0) 2688 (if (> n 0)
2668 (cdr (math-mul-bignum-digit a (if (= n 2) 10 100) 0)) 2689 (cdr (math-mul-bignum-digit a (expt 10 (- math-bignum-digit-length n)) 0))
2669 a)) 2690 a))
2670 2691
2671 ;;; Multiply (with rounding) the integer A by 10^N. [I i S] 2692 ;;; Multiply (with rounding) the integer A by 10^N. [I i S]
2672 (defun math-scale-rounding (a n) 2693 (defun math-scale-rounding (a n)
2673 (cond ((>= n 0) 2694 (cond ((>= n 0)
2674 (math-scale-left a n)) 2695 (math-scale-left a n))
2675 ((consp a) 2696 ((consp a)
2676 (math-normalize 2697 (math-normalize
2677 (cons (car a) 2698 (cons (car a)
2678 (let ((val (if (< n -3) 2699 (let ((val (if (< n (- math-bignum-digit-length))
2679 (math-scale-right-bignum (cdr a) (- -3 n)) 2700 (math-scale-right-bignum
2680 (if (= n -2) 2701 (cdr a)
2681 (math-mul-bignum-digit (cdr a) 10 0) 2702 (- (- math-bignum-digit-length) n))
2682 (if (= n -1) 2703 (if (< n 0)
2683 (math-mul-bignum-digit (cdr a) 100 0) 2704 (math-mul-bignum-digit
2684 (cdr a)))))) ; n = -3 2705 (cdr a)
2685 (if (and val (>= (car val) 500)) 2706 (expt 10 (+ math-bignum-digit-length n)) 0)
2707 (cdr a))))) ; n = -math-bignum-digit-length
2708 (if (and val (>= (car val) (/ math-bignum-digit-size 2)))
2686 (if (cdr val) 2709 (if (cdr val)
2687 (if (eq (car (cdr val)) 999) 2710 (if (eq (car (cdr val)) (1- math-bignum-digit-size))
2688 (math-add-bignum (cdr val) '(1)) 2711 (math-add-bignum (cdr val) '(1))
2689 (cons (1+ (car (cdr val))) (cdr (cdr val)))) 2712 (cons (1+ (car (cdr val))) (cdr (cdr val))))
2690 '(1)) 2713 '(1))
2691 (cdr val)))))) 2714 (cdr val))))))
2692 (t 2715 (t
2701 (defun math-add (a b) 2724 (defun math-add (a b)
2702 (or 2725 (or
2703 (and (not (or (consp a) (consp b))) 2726 (and (not (or (consp a) (consp b)))
2704 (progn 2727 (progn
2705 (setq a (+ a b)) 2728 (setq a (+ a b))
2706 (if (or (<= a -1000000) (>= a 1000000)) 2729 (if (or (<= a (- math-small-integer-size)) (>= a math-small-integer-size))
2707 (math-bignum a) 2730 (math-bignum a)
2708 a))) 2731 a)))
2709 (and (Math-zerop a) (not (eq (car-safe a) 'mod)) 2732 (and (Math-zerop a) (not (eq (car-safe a) 'mod))
2710 (if (and (math-floatp a) (Math-ratp b)) (math-float b) b)) 2733 (if (and (math-floatp a) (Math-ratp b)) (math-float b) b))
2711 (and (Math-zerop b) (not (eq (car-safe b) 'mod)) 2734 (and (Math-zerop b) (not (eq (car-safe b) 'mod))
2750 (if a 2773 (if a
2751 (if b 2774 (if b
2752 (let* ((a (copy-sequence a)) (aa a) (carry nil) sum) 2775 (let* ((a (copy-sequence a)) (aa a) (carry nil) sum)
2753 (while (and aa b) 2776 (while (and aa b)
2754 (if carry 2777 (if carry
2755 (if (< (setq sum (+ (car aa) (car b))) 999) 2778 (if (< (setq sum (+ (car aa) (car b)))
2779 (1- math-bignum-digit-size))
2756 (progn 2780 (progn
2757 (setcar aa (1+ sum)) 2781 (setcar aa (1+ sum))
2758 (setq carry nil)) 2782 (setq carry nil))
2759 (setcar aa (+ sum -999))) 2783 (setcar aa (+ sum -999)))
2760 (if (< (setq sum (+ (car aa) (car b))) 1000) 2784 (if (< (setq sum (+ (car aa) (car b))) math-bignum-digit-size)
2761 (setcar aa sum) 2785 (setcar aa sum)
2762 (setcar aa (+ sum -1000)) 2786 (setcar aa (- sum math-bignum-digit-size))
2763 (setq carry t))) 2787 (setq carry t)))
2764 (setq aa (cdr aa) 2788 (setq aa (cdr aa)
2765 b (cdr b))) 2789 b (cdr b)))
2766 (if carry 2790 (if carry
2767 (if b 2791 (if b
2788 (if borrow 2812 (if borrow
2789 (if (>= (setq diff (- (car aa) (car b))) 1) 2813 (if (>= (setq diff (- (car aa) (car b))) 1)
2790 (progn 2814 (progn
2791 (setcar aa (1- diff)) 2815 (setcar aa (1- diff))
2792 (setq borrow nil)) 2816 (setq borrow nil))
2793 (setcar aa (+ diff 999))) 2817 (setcar aa (+ diff (1- math-bignum-digit-size))))
2794 (if (>= (setq diff (- (car aa) (car b))) 0) 2818 (if (>= (setq diff (- (car aa) (car b))) 0)
2795 (setcar aa diff) 2819 (setcar aa diff)
2796 (setcar aa (+ diff 1000)) 2820 (setcar aa (+ diff math-bignum-digit-size))
2797 (setq borrow t))) 2821 (setq borrow t)))
2798 (setq aa (cdr aa) 2822 (setq aa (cdr aa)
2799 b (cdr b))) 2823 b (cdr b)))
2800 (if borrow 2824 (if borrow
2801 (progn 2825 (progn
2802 (while (eq (car aa) 0) 2826 (while (eq (car aa) 0)
2803 (setcar aa 999) 2827 (setcar aa (1- math-bignum-digit-size))
2804 (setq aa (cdr aa))) 2828 (setq aa (cdr aa)))
2805 (if aa 2829 (if aa
2806 (progn 2830 (progn
2807 (setcar aa (1- (car aa))) 2831 (setcar aa (1- (car aa)))
2808 a) 2832 a)
2838 ;;; Compute the difference of A and B. [O O O] [Public] 2862 ;;; Compute the difference of A and B. [O O O] [Public]
2839 (defun math-sub (a b) 2863 (defun math-sub (a b)
2840 (if (or (consp a) (consp b)) 2864 (if (or (consp a) (consp b))
2841 (math-add a (math-neg b)) 2865 (math-add a (math-neg b))
2842 (setq a (- a b)) 2866 (setq a (- a b))
2843 (if (or (<= a -1000000) (>= a 1000000)) 2867 (if (or (<= a (- math-small-integer-size)) (>= a math-small-integer-size))
2844 (math-bignum a) 2868 (math-bignum a)
2845 a))) 2869 a)))
2846 2870
2847 (defun math-sub-float (a b) ; [F F F] 2871 (defun math-sub-float (a b) ; [F F F]
2848 (let ((ediff (- (nth 2 a) (nth 2 b)))) 2872 (let ((ediff (- (nth 2 a) (nth 2 b))))
2865 2889
2866 ;;; Compute the product of A and B. [O O O] [Public] 2890 ;;; Compute the product of A and B. [O O O] [Public]
2867 (defun math-mul (a b) 2891 (defun math-mul (a b)
2868 (or 2892 (or
2869 (and (not (consp a)) (not (consp b)) 2893 (and (not (consp a)) (not (consp b))
2870 (< a 1000) (> a -1000) (< b 1000) (> b -1000) 2894 (< a math-bignum-digit-size) (> a (- math-bignum-digit-size))
2895 (< b math-bignum-digit-size) (> b (- math-bignum-digit-size))
2871 (* a b)) 2896 (* a b))
2872 (and (Math-zerop a) (not (eq (car-safe b) 'mod)) 2897 (and (Math-zerop a) (not (eq (car-safe b) 'mod))
2873 (if (Math-scalarp b) 2898 (if (Math-scalarp b)
2874 (if (and (math-floatp b) (Math-ratp a)) (math-float a) a) 2899 (if (and (math-floatp b) (Math-ratp a)) (math-float a) a)
2875 (require 'calc-ext) 2900 (require 'calc-ext)
2934 d (car b) 2959 d (car b)
2935 c 0 2960 c 0
2936 aa a) 2961 aa a)
2937 (while (progn 2962 (while (progn
2938 (setcar ss (% (setq prod (+ (+ (car ss) (* (car aa) d)) 2963 (setcar ss (% (setq prod (+ (+ (car ss) (* (car aa) d))
2939 c)) 1000)) 2964 c)) math-bignum-digit-size))
2940 (setq aa (cdr aa))) 2965 (setq aa (cdr aa)))
2941 (setq c (/ prod 1000) 2966 (setq c (/ prod math-bignum-digit-size)
2942 ss (or (cdr ss) (setcdr ss (list 0))))) 2967 ss (or (cdr ss) (setcdr ss (list 0)))))
2943 (if (>= prod 1000) 2968 (if (>= prod math-bignum-digit-size)
2944 (if (cdr ss) 2969 (if (cdr ss)
2945 (setcar (cdr ss) (+ (/ prod 1000) (car (cdr ss)))) 2970 (setcar (cdr ss) (+ (/ prod math-bignum-digit-size) (car (cdr ss))))
2946 (setcdr ss (list (/ prod 1000)))))) 2971 (setcdr ss (list (/ prod math-bignum-digit-size))))))
2947 sum))) 2972 sum)))
2948 2973
2949 ;;; Multiply digit list A by digit D. [L L D D; l l D D] 2974 ;;; Multiply digit list A by digit D. [L L D D; l l D D]
2950 (defun math-mul-bignum-digit (a d c) 2975 (defun math-mul-bignum-digit (a d c)
2951 (if a 2976 (if a
2952 (if (<= d 1) 2977 (if (<= d 1)
2953 (and (= d 1) a) 2978 (and (= d 1) a)
2954 (let* ((a (copy-sequence a)) (aa a) prod) 2979 (let* ((a (copy-sequence a)) (aa a) prod)
2955 (while (progn 2980 (while (progn
2956 (setcar aa (% (setq prod (+ (* (car aa) d) c)) 1000)) 2981 (setcar aa
2982 (% (setq prod (+ (* (car aa) d) c))
2983 math-bignum-digit-size))
2957 (cdr aa)) 2984 (cdr aa))
2958 (setq aa (cdr aa) 2985 (setq aa (cdr aa)
2959 c (/ prod 1000))) 2986 c (/ prod math-bignum-digit-size)))
2960 (if (>= prod 1000) 2987 (if (>= prod math-bignum-digit-size)
2961 (setcdr aa (list (/ prod 1000)))) 2988 (setcdr aa (list (/ prod math-bignum-digit-size))))
2962 a)) 2989 a))
2963 (and (> c 0) 2990 (and (> c 0)
2964 (list c)))) 2991 (list c))))
2965 2992
2966 2993
2969 ;;; if A or B is negative. B must be nonzero. [I.I I I] [Public] 2996 ;;; if A or B is negative. B must be nonzero. [I.I I I] [Public]
2970 (defun math-idivmod (a b) 2997 (defun math-idivmod (a b)
2971 (if (eq b 0) 2998 (if (eq b 0)
2972 (math-reject-arg a "*Division by zero")) 2999 (math-reject-arg a "*Division by zero"))
2973 (if (or (consp a) (consp b)) 3000 (if (or (consp a) (consp b))
2974 (if (and (natnump b) (< b 1000)) 3001 (if (and (natnump b) (< b math-bignum-digit-size))
2975 (let ((res (math-div-bignum-digit (cdr a) b))) 3002 (let ((res (math-div-bignum-digit (cdr a) b)))
2976 (cons 3003 (cons
2977 (math-normalize (cons (car a) (car res))) 3004 (math-normalize (cons (car a) (car res)))
2978 (cdr res))) 3005 (cdr res)))
2979 (or (consp a) (setq a (math-bignum a))) 3006 (or (consp a) (setq a (math-bignum a)))
2988 (defun math-quotient (a b) ; [I I I] [Public] 3015 (defun math-quotient (a b) ; [I I I] [Public]
2989 (if (and (not (consp a)) (not (consp b))) 3016 (if (and (not (consp a)) (not (consp b)))
2990 (if (= b 0) 3017 (if (= b 0)
2991 (math-reject-arg a "*Division by zero") 3018 (math-reject-arg a "*Division by zero")
2992 (/ a b)) 3019 (/ a b))
2993 (if (and (natnump b) (< b 1000)) 3020 (if (and (natnump b) (< b math-bignum-digit-size))
2994 (if (= b 0) 3021 (if (= b 0)
2995 (math-reject-arg a "*Division by zero") 3022 (math-reject-arg a "*Division by zero")
2996 (math-normalize (cons (car a) 3023 (math-normalize (cons (car a)
2997 (car (math-div-bignum-digit (cdr a) b))))) 3024 (car (math-div-bignum-digit (cdr a) b)))))
2998 (or (consp a) (setq a (math-bignum a))) 3025 (or (consp a) (setq a (math-bignum a)))
2999 (or (consp b) (setq b (math-bignum b))) 3026 (or (consp b) (setq b (math-bignum b)))
3000 (let* ((alen (1- (length a))) 3027 (let* ((alen (1- (length a)))
3001 (blen (1- (length b))) 3028 (blen (1- (length b)))
3002 (d (/ 1000 (1+ (nth (1- blen) (cdr b))))) 3029 (d (/ math-bignum-digit-size (1+ (nth (1- blen) (cdr b)))))
3003 (res (math-div-bignum-big (math-mul-bignum-digit (cdr a) d 0) 3030 (res (math-div-bignum-big (math-mul-bignum-digit (cdr a) d 0)
3004 (math-mul-bignum-digit (cdr b) d 0) 3031 (math-mul-bignum-digit (cdr b) d 0)
3005 alen blen))) 3032 alen blen)))
3006 (math-normalize (cons (if (eq (car a) (car b)) 'bigpos 'bigneg) 3033 (math-normalize (cons (if (eq (car a) (car b)) 'bigpos 'bigneg)
3007 (car res))))))) 3034 (car res)))))))
3011 ;;; The following division algorithm is borrowed from Knuth vol. II, sec. 4.3.1 3038 ;;; The following division algorithm is borrowed from Knuth vol. II, sec. 4.3.1
3012 (defun math-div-bignum (a b) 3039 (defun math-div-bignum (a b)
3013 (if (cdr b) 3040 (if (cdr b)
3014 (let* ((alen (length a)) 3041 (let* ((alen (length a))
3015 (blen (length b)) 3042 (blen (length b))
3016 (d (/ 1000 (1+ (nth (1- blen) b)))) 3043 (d (/ math-bignum-digit-size (1+ (nth (1- blen) b))))
3017 (res (math-div-bignum-big (math-mul-bignum-digit a d 0) 3044 (res (math-div-bignum-big (math-mul-bignum-digit a d 0)
3018 (math-mul-bignum-digit b d 0) 3045 (math-mul-bignum-digit b d 0)
3019 alen blen))) 3046 alen blen)))
3020 (if (= d 1) 3047 (if (= d 1)
3021 res 3048 res
3026 3053
3027 ;;; Divide a bignum digit list by a digit. [l.D l D] 3054 ;;; Divide a bignum digit list by a digit. [l.D l D]
3028 (defun math-div-bignum-digit (a b) 3055 (defun math-div-bignum-digit (a b)
3029 (if a 3056 (if a
3030 (let* ((res (math-div-bignum-digit (cdr a) b)) 3057 (let* ((res (math-div-bignum-digit (cdr a) b))
3031 (num (+ (* (cdr res) 1000) (car a)))) 3058 (num (+ (* (cdr res) math-bignum-digit-size) (car a))))
3032 (cons 3059 (cons
3033 (cons (/ num b) (car res)) 3060 (cons (/ num b) (car res))
3034 (% num b))) 3061 (% num b)))
3035 '(nil . 0))) 3062 '(nil . 0)))
3036 3063
3042 (res2 (math-div-bignum-part num b blen))) 3069 (res2 (math-div-bignum-part num b blen)))
3043 (cons 3070 (cons
3044 (cons (car res2) (car res)) 3071 (cons (car res2) (car res))
3045 (cdr res2))))) 3072 (cdr res2)))))
3046 3073
3047 (defun math-div-bignum-part (a b blen) ; a < b*1000 [D.l l L] 3074 (defun math-div-bignum-part (a b blen) ; a < b*math-bignum-digit-size [D.l l L]
3048 (let* ((num (+ (* (or (nth blen a) 0) 1000) (or (nth (1- blen) a) 0))) 3075 (let* ((num (+ (* (or (nth blen a) 0) math-bignum-digit-size)
3076 (or (nth (1- blen) a) 0)))
3049 (den (nth (1- blen) b)) 3077 (den (nth (1- blen) b))
3050 (guess (min (/ num den) 999))) 3078 (guess (min (/ num den) (1- math-bignum-digit-size))))
3051 (math-div-bignum-try a b (math-mul-bignum-digit b guess 0) guess))) 3079 (math-div-bignum-try a b (math-mul-bignum-digit b guess 0) guess)))
3052 3080
3053 (defun math-div-bignum-try (a b c guess) ; [D.l l l D] 3081 (defun math-div-bignum-try (a b c guess) ; [D.l l l D]
3054 (let ((rem (math-sub-bignum a c))) 3082 (let ((rem (math-sub-bignum a c)))
3055 (if (eq rem 'neg) 3083 (if (eq rem 'neg)
3356 3384
3357 (defun math-format-bignum-decimal (a) ; [X L] 3385 (defun math-format-bignum-decimal (a) ; [X L]
3358 (if a 3386 (if a
3359 (let ((s "")) 3387 (let ((s ""))
3360 (while (cdr (cdr a)) 3388 (while (cdr (cdr a))
3361 (setq s (concat (format "%06d" (+ (* (nth 1 a) 1000) (car a))) s) 3389 (setq s (concat
3390 (format
3391 (concat "%0"
3392 (number-to-string (* 2 math-bignum-digit-length))
3393 "d")
3394 (+ (* (nth 1 a) math-bignum-digit-size) (car a))) s)
3362 a (cdr (cdr a)))) 3395 a (cdr (cdr a))))
3363 (concat (int-to-string (+ (* (or (nth 1 a) 0) 1000) (car a))) s)) 3396 (concat (int-to-string
3397 (+ (* (or (nth 1 a) 0) math-bignum-digit-size) (car a))) s))
3364 "0")) 3398 "0"))
3365 3399
3366 3400
3367 3401
3368 ;;; Parse a simple number in string form. [N X] [Public] 3402 ;;; Parse a simple number in string form. [N X] [Public]
3445 (if (match-beginning n) 3479 (if (match-beginning n)
3446 (substring s (match-beginning n) (match-end n)) 3480 (substring s (match-beginning n) (match-end n))
3447 "")) 3481 ""))
3448 3482
3449 (defun math-read-bignum (s) ; [l X] 3483 (defun math-read-bignum (s) ; [l X]
3450 (if (> (length s) 3) 3484 (if (> (length s) math-bignum-digit-length)
3451 (cons (string-to-number (substring s -3)) 3485 (cons (string-to-number (substring s (- math-bignum-digit-length)))
3452 (math-read-bignum (substring s 0 -3))) 3486 (math-read-bignum (substring s 0 (- math-bignum-digit-length))))
3453 (list (string-to-number s)))) 3487 (list (string-to-number s))))
3454 3488
3455 3489
3456 (defconst math-tex-ignore-words 3490 (defconst math-tex-ignore-words
3457 '( ("\\hbox") ("\\mbox") ("\\text") ("\\left") ("\\right") 3491 '( ("\\hbox") ("\\mbox") ("\\text") ("\\left") ("\\right")