Guest User

factor

a guest
Feb 8th, 2019
118
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. (define start
  2.     ;;; This program is a tribute to the test program 'factor.nj', which severely penetrated my implementation of the garbage collector for the Ninja VM.
  3.     (lambda ()
  4.         ;(small-primes 100)
  5.         (display "upper limit: ")
  6.         (let ((lim (read)))
  7.             (compute 1 lim)
  8.         )
  9.     )
  10. )
  11. (define small-primes-result '())
  12. (define (small-primes limit . rest) ; get a list of small primes up to limit
  13.     (if (< limit 0)
  14.         (begin
  15.             (set! small-primes-result '()) ;flush primes with a negative arg
  16.             ()
  17.         )
  18.         (if (null? small-primes-result)
  19.             (let (( i (if (null? rest) 0 (car rest)) ))
  20.                 (if (>= i limit)
  21.                     ()
  22.                     (if (= i 0)
  23.                         (begin
  24.                             (set! small-primes-result (cons 2 (cons 3 (small-primes limit 5))))
  25.                             small-primes-result
  26.                         )
  27.                         (if (is-prime i)
  28.                             (cons i (small-primes limit (+ i 2)))
  29.                             (small-primes limit (+ i 2))
  30.                         )
  31.                     )
  32.                 )
  33.             )
  34.             small-primes-result
  35.         )
  36.     )
  37. )
  38. (define (is-prime n . rest)
  39.     (let (( i (if (null? rest) 0 (car rest)) ))
  40.         (if (< i 3)
  41.             (is-prime n 3) ;set start value -> i=3
  42.             (if (<= (square i) n)
  43.                 (if (= (modulo n i) 0)
  44.                     #f
  45.                     (is-prime n (+ i 2))
  46.                 )
  47.                 #t
  48.             )
  49.         )
  50.     )
  51. )
  52. (define (compute n limit)
  53.     (if (= n limit)
  54.         #t
  55.         (begin
  56.             (display "10^")
  57.             (display n)
  58.             (display "+1 = ")
  59.             (let ((m (+ (expt 10 n) 1)))
  60.                 (fmt m)
  61.                 (let find-fac ((factors (factorize m)))
  62.                     (fmt "result: " 0)
  63.                     (fmt factors)
  64.                     ; check if all primes -> if yes, check product. if no, factorize remaining composites
  65.                     (fmt "check: product = " 0)
  66.                     (let ((prod (reduce factors (lambda (x y) (* x y)) 1)))
  67.                         (fmt prod)
  68.                         (if (not (= m prod))
  69.                             #f
  70.                             (compute (+ n 1) limit)
  71.                         )
  72.                     )
  73.                 )
  74.             )
  75.         )
  76.     )
  77. )
  78. (define (all-primes l)
  79.     (if (null? l)
  80.         #t
  81.         (if (is-prime (car l))
  82.             (all-primes (cdr l))
  83.             #f
  84.         )
  85.     )
  86. )
  87. (define (factorize n) ;factorize n using trial division
  88.     (if (< n 2)
  89.         ()
  90.         (let ((p (small-prime-factor n)))
  91.             (if (> p -1)
  92.                 (if (< (square p) n)
  93.                     (begin
  94.                         (fmt "detected small prime factor " 0)
  95.                         (fmt p)
  96.                         (cons p (factorize (/ n p)))
  97.                     )
  98.                     (begin
  99.                         (fmt "the remaining factor " 0)
  100.                         (fmt n 0)
  101.                         (fmt " is a small prime")
  102.                         (cons n ())
  103.                     )
  104.                 )
  105.                 (begin
  106.                     (fmt "the remaining factor " 0)
  107.                     (fmt n 0)
  108.                     (fmt " doesn't have any prime factors < " 0)
  109.                     (fmt SMALL_PRIMES_LIM)
  110.                     (if #t  ;;; disable composite testing for now
  111.                         (cons n ())
  112.                         (if (is-composite n)
  113.                             (begin
  114.                                 (fmt "but is definitely composite")
  115.                                 (let ((f1 (find-factor n)))
  116.                                     (if (= f1 -1)
  117.                                         (begin
  118.                                             (fmt "cannot factorize " 0)
  119.                                             (fmt n 0)
  120.                                             (fmt ", giving up")
  121.                                             (cons n ())
  122.                                         )
  123.                                         (let ((f2 (/ n f1)))
  124.                                             (fmt "this number can be split into " 0)
  125.                                             (fmt f1 0)
  126.                                             (fmt " and " 0)
  127.                                             (fmt f2)
  128.                                             (cons (factorize f1) (cons (factorize f2)))
  129.                                         )
  130.                                     )
  131.                                 )
  132.                             )
  133.                             (begin
  134.                                 (fmt "and is very probably prime")
  135.                                 ;if prove-prime(x) then
  136.                                 ;   primality has been proven
  137.                                 (cons n ())
  138.                             )
  139.                         )
  140.                     )
  141.                 )
  142.             )
  143.         )
  144.     )
  145. )
  146. (define SMALL_PRIMES_LIM 2000)
  147. (define (small-prime-factor n)
  148.     ; find a small prime factor of n (using prime list plist)
  149.     ; -1: not found
  150.     (let cont ((plist (small-primes SMALL_PRIMES_LIM)))
  151.         (if (< (length plist) 1)
  152.             -1
  153.             (if (= (modulo n (car plist)) 0)
  154.                 (car plist)
  155.                 (cont (cdr plist))
  156.             )
  157.         )
  158.     )
  159. )
  160. (define (is-composite n)     ;Rabin-Miller test
  161.     (let cont1 ((q (- n 1)) (t 0))
  162.         (if (= (modulo q 2) 0)
  163.             (cont1 (/ q 2) (+ t 1))
  164.             'todo ;...
  165.         )
  166.     )
  167. )
  168. ;procedure pollard(n):
  169. ;   x<-2, y<-2, d<-1
  170. ;   while d=1 do
  171. ;       x<-g(x)=(x^2-1) mod n
  172. ;       y<-g(g(y))=((y^2-1) mod n)^2-1 mod n)
  173. ;       d<-gcd(|x-y|,n)
  174. ;   if d=n then
  175. ;       return FAIL
  176. ;   else
  177. ;       return n
  178.  
  179. (define (find-factor n)      ;Pollard rho method
  180.  'todo 
  181. )
  182.  
  183. (start)
  184.  
  185. ;module tests
  186. (define (test-is-prime)
  187.     (fmt (is-prime 71 0))
  188. )
  189. (define (test-small-primes)
  190.     (fmt (small-primes 100))
  191. )
  192. (define (test-factorize)
  193.     (assert-eq (factorize 144) '(2 2 2 2 3 3))
  194. )
  195. (define (test-compute)
  196.     (compute 1 40)
  197. )
RAW Paste Data