Project:
 ```1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 ``` ```; baillie-wagstaff pseudoprimality test (define (isqrt n) (if (not (and (positive? n) (integer? n))) (error 'isqrt "must be positive integer") (let loop ((x n)) (let ((y (quotient (+ x (quotient n x)) 2))) (if (< y x) (loop y) x))))) (define (square? n) (let ((n2 (isqrt n))) (= (* n2 n2) n))) (define (expm b e m) (define (m* x y) (modulo (* x y) m)) (cond ((zero? e) 1) ((even? e) (expm (m* b b) (/ e 2) m)) (else (m* b (expm (m* b b) (/ (- e 1) 2) m))))) (define (jacobi a m) (if (not (integer? a)) (error 'jacobi "must be integer") (if (not (and (integer? m) (positive? m) (odd? m))) (error 'jacobi "modulus must be odd positive integer") (let loop1 ((a (modulo a m)) (m m) (t 1)) (if (zero? a) (if (= m 1) t 0) (let ((z (if (member (modulo m 8) (list 3 5)) -1 1))) (let loop2 ((a a) (t t)) (if (even? a) (loop2 (/ a 2) (* t z)) (loop1 (modulo m a) a (if (and (= (modulo a 4) 3) (= (modulo m 4) 3)) (- t) t)))))))))) (define (primes n) (let ((sieve (make-vector n #t))) (let loop ((p 2) (ps (list))) (cond ((= n p) (reverse ps)) ((vector-ref sieve p) (do ((i p (+ i p))) ((<= n i)) (vector-set! sieve i #f)) (loop (+ p 1) (cons p ps))) (else (loop (+ p 1) ps)))))) (define (strong-pseudoprime? n a) (let loop ((r 0) (s (- n 1))) (if (even? s) (loop (+ r 1) (/ s 2)) (if (= (expm a s n) 1) #t (let loop ((r r) (s s)) (cond ((zero? r) #f) ((= (expm a s n) (- n 1)) #t) (else (loop (- r 1) (* s 2))))))))) (define (selfridge n) (let loop ((d-abs 5) (sign 1)) (let ((d (* d-abs sign))) (cond ((< 1 (gcd d n)) (values d 0 0)) ((= (jacobi d n) -1) (values d 1 (/ (- 1 d) 4))) (else (loop (+ d-abs 2) (- sign))))))) (define (chain n u v u2 v2 d q m) (let loop ((u u) (v v) (u2 u2) (v2 v2) (q q) (qkd q) (m m)) (if (zero? m) (values u v qkd) (let* ((u2 (modulo (* u2 v2) n)) (v2 (modulo (- (* v2 v2) (* 2 q)) n)) (q (modulo (* q q) n)) (qkd (modulo (* qkd q) n))) (if (odd? m) (let* ((t1 (* u2 v)) (t2 (* u v2)) (t3 (* v2 v)) (t4 (* u2 u d)) (u (+ t1 t2)) (v (+ t3 t4))) (loop (modulo (quotient (if (odd? u) (+ u n) u) 2) n) (modulo (quotient (if (odd? v) (+ v n) v) 2) n) u2 v2 q qkd (quotient m 2))) (loop u v u2 v2 q qkd (quotient m 2))))))) (define (standard-lucas-pseudoprime? n) (call-with-values (lambda () (selfridge n)) (lambda (d p q) (if (zero? p) (= n d) (call-with-values (lambda () (chain n 0 2 1 p d q (/ (+ n 1) 2))) (lambda (u v qkd) (zero? u))))))) (define (powers-of-two n) (let loop ((s 0) (n n)) (if (odd? n) (values s n) (loop (+ s 1) (/ n 2))))) (define (strong-lucas-pseudoprime? n) (call-with-values (lambda () (selfridge n)) (lambda (d p q) (if (zero? p) (= n d) (call-with-values (lambda () (powers-of-two (+ n 1))) (lambda (s t) (call-with-values (lambda () (chain n 1 p 1 p d q (quotient t 2))) (lambda (u v qkd) (if (or (zero? u) (zero? v)) #t (let loop ((r 1) (v v) (qkd qkd)) (if (= r s) #f (let* ((v (modulo (- (* v v) (* 2 qkd)) n)) (qkd (modulo (* qkd qkd) n))) (if (zero? v) #t (loop (+ r 1) v qkd)))))))))))))) (define prime? (let ((ps (primes 100))) (lambda (n) (and (integer? n) (< 1 n) (not (square? n)) (let loop ((ps ps)) (cond ((null? ps) #t) ((zero? (modulo n (car ps))) (= n (car ps))) (else (loop (cdr ps))))) (strong-pseudoprime? n 2) (strong-pseudoprime? n 3) (strong-lucas-pseudoprime? n))))) (do ((ns (list 83 89 111 323 5777 3825123056546413051 (- (expt 2 89) 1)) (cdr ns))) ((null? ns)) (let ((n (car ns))) (display n) (display " ") (display (strong-pseudoprime? n 2)) (display " ") (display (strong-pseudoprime? n 3)) (display " ") (display (standard-lucas-pseudoprime? n)) (display " ") (display (strong-lucas-pseudoprime? n)) (display " ") (display (prime? n)) (newline))) ```
 ```1 2 3 4 5 6 7 ``` ```83 #t #t #t #t #t 89 #t #t #t #t #t 111 #f #f #f #f #f 323 #f #f #t #f #f 5777 #f #f #t #t #f 3825123056546413051 #t #t #f #f #f 618970019642690137449562111 #t #t #t #t #t ```