#lang racket

(provide (contract-out
          (prime? (-> integer? integer? boolean?))))

(define (mod-op op a b n)
  (modulo (op 
           (modulo a n)
           (modulo b n))
          n))

(define (get-s-d x)
  (let find-s-d ((s 0)
                 (d x))
    (if (even? d)
        (find-s-d (add1 s) (arithmetic-shift d -1))
        (values s d))))

(define (mod-exp a x n)
  (cond ((zero? x) 1)
        ((odd? x)
         (mod-op * (mod-exp (mod-op * a a n) (quotient x 2) n) a n))
        ((even? x)
         (mod-exp (mod-op * a a n) (quotient x 2) n))))

(define (prime? x k)
  (let-values (((s d) (get-s-d (sub1 x))))
    (let loop ((i 0))
      (cond ((>= i k)
             #t)
            (else 
             (let ((test (mod-exp (add1 (random (- x 3))) d x)))
               (if (or (= test 1) (= test (sub1 x)))
                   (loop (add1 i))
                   (let carmichael-test ((r 1)
                                         (ctest (modulo (* test test) x)))
                     (cond ((>= r s)
                            #f)
                           ((= ctest 1)
                            #f)
                           ((= ctest (sub1 x))
                            (loop (add1 i)))
                           (else
                            (carmichael-test (add1 r) (modulo (* ctest ctest) x))))))))))))