codepad
[
create a new paste
]
login
|
about
Language:
C
C++
D
Haskell
Lua
OCaml
PHP
Perl
Plain Text
Python
Ruby
Scheme
Tcl
; home primes (define (digits n . args) (let ((b (if (null? args) 10 (car args)))) (let loop ((n n) (d '())) (if (zero? n) d (loop (quotient n b) (cons (modulo n b) d)))))) (define (undigits ds . args) (let ((b (if (null? args) 10 (car args)))) (let loop ((ds ds) (n 0)) (if (null? ds) n (loop (cdr ds) (+ (* n b) (car ds))))))) (define sort #f) (define merge #f) (let () (define dosort (lambda (pred? ls n) (if (= n 1) (list (car ls)) (let ((i (quotient n 2))) (domerge pred? (dosort pred? ls i) (dosort pred? (list-tail ls i) (- n i))))))) (define domerge (lambda (pred? l1 l2) (cond ((null? l1) l2) ((null? l2) l1) ((pred? (car l2) (car l1)) (cons (car l2) (domerge pred? l1 (cdr l2)))) (else (cons (car l1) (domerge pred? (cdr l1) l2)))))) (set! sort (lambda (pred? l) (if (null? l) l (dosort pred? l (length l))))) (set! merge (lambda (pred? l1 l2) (domerge pred? l1 l2)))) (define (prime? n) (letrec ( (expm (lambda (b e m) (let ((times (lambda (x y) (modulo (* x y) m)))) (cond ((zero? e) 1) ((even? e) (expm (times b b) (/ e 2) m)) (else (times b (expm (times b b) (/ (- e 1) 2) m))))))) (digits (lambda (n) (let loop ((n n) (ds '())) (if (zero? n) ds (loop (quotient n 2) (cons (modulo n 2) ds)))))) (isqrt (lambda (n) (let loop ((x n) (y (quotient (+ n 1) 2))) (if (<= 0 (- y x) 1) x (loop y (quotient (+ y (quotient n y)) 2)))))) (square? (lambda (n) (let ((n2 (isqrt n))) (= n (* n2 n2))))) (jacobi (lambda (a n) (let loop ((a a) (n n)) (cond ((= a 0) 0) ((= a 1) 1) ((= a 2) (case (modulo n 8) ((1 7) 1) ((3 5) -1))) ((even? a) (* (loop 2 n) (loop (/ a 2) n))) ((< n a) (loop (modulo a n) n)) ((and (= (modulo a 4) 3) (= (modulo n 4) 3)) (- (loop n a))) (else (loop n a)))))) (inverse (lambda (x m) (let loop ((a 1) (b 0) (g x) (u 0) (v 1) (w m)) (if (zero? w) (modulo a m) (let ((q (quotient g w))) (loop u v w (- a (* q u)) (- b (* q v)) (- g (* q w)))))))) (strong-pseudo-prime? (lambda (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)))))))))) (chain (lambda (m f g u v) (let loop ((ms (digits m)) (u u) (v v)) (cond ((null? ms) (values u v)) ((zero? (car ms)) (loop (cdr ms) (f u) (g u v))) (else (loop (cdr ms) (g u v) (f v))))))) (lucas-pseudo-prime? (lambda (n) (let loop ((a 11) (b 7)) (let ((d (- (* a a) (* 4 b)))) (cond ((square? d) (loop (+ a 2) (+ b 1))) ((not (= (gcd n (* 2 a b d)) 1)) (loop (+ a 2) (+ b 2))) (else (let* ((x1 (modulo (- (* a a (inverse b n)) 2) n)) (m (quotient (- n (jacobi d n)) 2)) (f (lambda (u) (modulo (- (* u u) 2) n))) (g (lambda (u v) (modulo (- (* u v) x1) n)))) (let-values (((xm xm1) (chain m f g 2 x1))) (zero? (modulo (- (* x1 xm) (* 2 xm1)) n))))))))))) (if (not (integer? n)) (error 'prime? "must be integer") (if (< n 2) #f (if (even? n) (= n 2) (if (zero? (modulo n 3)) (= n 3) (and (strong-pseudo-prime? n 2) (strong-pseudo-prime? n 3) (lucas-pseudo-prime? n)))))))) (define (factors n) (letrec ( (last-pair (lambda (xs) (if (pair? (cdr xs)) (last-pair (cdr xs)) xs))) (cycle (lambda xs (set-cdr! (last-pair xs) xs) xs)) (cons* (lambda (x . xs) (if (null? xs) x (cons x (apply cons* (car xs) (cdr xs))))))) (if (not (integer? n)) (error 'factors "must be integer") (if (member n '(-1 0 1)) (cdr (assoc n '((-1 -1) (0 0) (1 1)))) (if (negative? n) (cons -1 (factors (- n))) (let wheel ((n n) (p 2) (ps '()) (ws (cons* 1 2 2 (cycle 4 2 4 2 4 6 2 6)))) (if (= n 1) (reverse ps) (if (< n (* p p)) (reverse (cons n ps)) (if (zero? (modulo n p)) (wheel (/ n p) p (cons p ps) ws) (if (< p 1000) (wheel n (+ p (car ws)) ps (cdr ws)) (let rho ((ps ps) (cs (list n)) (a 1) (m 10000000)) (letrec ( (f (lambda (y) (modulo (+ (* y y) a) (car cs)))) (g (lambda (p x y) (modulo (* p (abs (- x y))) (car cs))))) (if (null? cs) (sort < ps) (let ((n (car cs))) (if (= n 1) (rho ps (cdr cs) a m) (if (prime? n) (rho (cons n ps) (cdr cs) a m) (let loop ((x 2) (y (f 2)) (j 1) (q 2) (p 1)) (if (= j m) (cons (sort < ps) cs) (if (= j q) (let ((t (f y))) (loop y t (+ j 1) (* q 2) (g p y t))) (if (positive? (modulo j 100)) (loop x (f y) (+ j 1) q (g p x y)) (let ((d (gcd p n))) (if (= d 1) (loop x (f y) (+ j 1) q 1) (if (= d n) (rho ps cs (+ a 1) (- m j)) (rho ps (cons* d (/ n d) (cdr cs)) a (- m j))))))))))))))))))))))))) (define (home-prime n) (let* ((fs (factors n)) (hp (undigits (apply append (map digits fs))))) (if (prime? hp) hp (home-prime hp)))) (define (home-factors n) (let loop ((n n) (hps '())) (if (prime? n) (reverse (cons n hps)) (let* ((fs (factors n)) (hp (undigits (apply append (map digits fs))))) (loop hp (cons fs hps)))))) (display (home-prime 99)) (newline) (display (home-factors 8)) (newline) (define (td-least n) (let loop ((p 3)) (cond ((< 1000000 p) #f) ((zero? (modulo n p)) p) (else (loop (+ p 2)))))) (define (euclid-mullin) (let loop ((n 1) (ems '(2))) (display n) (display " ") (display (car ems)) (newline) (let* ((next (+ (apply * ems) 1)) (lpf (td-least next))) (if lpf (loop (+ n 1) (cons lpf ems)) (let ((fs (factors next))) (cond ((= n 16) (loop 17 (cons 30693651606209 ems))) ((= n 27) (loop 28 (cons 643679794963466223081509857 ems))) ((number? (car fs)) (loop (+ n 1) (cons (car fs) ems))) ((pair? (car fs)) (loop (+ n 1) (cons (caar fs) ems))) (else (error 'euclid-mullin "can't complete factorization")))))))) (euclid-mullin)
Private
[
?
]
Run code
Submit