;; This LISP program computes the continued fraction for pi using ;; W. Gosper's technique. Update formulas: ;; r+2+p p q m is leading term ratio as (quo rem) ;; r+2 r s l is last term ratio as (quo rem) ;; ar+bq a b ;; cr+dq c d -S.McD 5/4/83 ;; (defun gcd(a b) ;; gcd(a,b) ;; (prog (rem) ;; loop (setq rem (remainder a b)) ;; (if (zerop rem) (return b)) ;; (setq a b)(setq b rem)(go loop))) #lang racket (define-syntax mdefine (syntax-rules () ((_ id ...) (begin (define id 'undefined) ...)))) (mdefine tmp a b c d e f g h i j k l m n o p q r s t u v w x y z) (define (input) (set! tmp a) (set! a (+ (* a r) (* b q))) (set! b tmp) (set! tmp c) (set! c (+ (* c r) (* d q))) (set! d tmp) (set! tmp r) (set! r (+ r 2)) (set! s tmp) (set! tmp p) (set! p (+ p r)) (set! o q) (set! q tmp) (set! tmp (gcd d o)) (when (> tmp 1) (set! tmp (gcd b tmp)) (when (> tmp 1) (set! tmp (gcd c tmp)) (when (> tmp 1) (set! tmp (gcd a tmp))))) (when (> tmp 1) (set! a (*quo a tmp)) (set! b (*quo b tmp)) (set! c (*quo c tmp)) (set! d (*quo d tmp))) (when (positive? d) (set! l (divide b d)))) (define (*quo a b) (/ a b)) (define (divide a b) (let-values (((q r) (quotient/remainder a b))) (list q r))) (define patom write) (define (output) (set! m (car m)) (set! tmp c) (set! c (- a (* c m))) (set! a tmp) (set! tmp d) (set! d (- b (* d m))) (set! b tmp) (set! l (divide b d)) (set! n (sub1 n)) (write m) (newline)) (define (decide) (let/ec return (when (zero? c) (return 0)) (when (zero? d) (return 0)) (set! m (divide a c)) (when (equal? (car m) (car l)) (return 1)) (when (equal? (car m) (sub1 (car l))) (when (zero? (cadr l)) (return 1))) (set! l m) (return 0))) (define (init) (set! p 1) (set! q 1) (set! a 0) (set! b 4) (set! r 1) (set! c 1) (set! d 0) (set! l (divide a c))) (define (cfpi nn) (set! n nn) ; compute first nn coeffs (do () ((zero? n)) (if (positive? (decide)) (output) (input))) (patom -99999) (newline)) (init) (cfpi 65536) (exit)