; rsa backdoor encryption

(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 (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?
  (let ((q11 (make-vector 11 #f))
        (q63 (make-vector 63 #f))
        (q64 (make-vector 64 #f))
        (q65 (make-vector 65 #f)))
    (do ((k 0 (+ k 1))) ((< 5 k))
      (vector-set! q11 (modulo (* k k) 11) #t))
    (do ((k 0 (+ k 1))) ((< 31 k))
      (vector-set! q63 (modulo (* k k) 63) #t))
    (do ((k 0 (+ k 1))) ((< 31 k))
      (vector-set! q64 (modulo (* k k) 64) #t))
    (do ((k 0 (+ k 1))) ((< 32 k))
      (vector-set! q65 (modulo (* k k) 65) #t))
    (lambda (n)
      (if (not (vector-ref q64 (modulo n 64))) #f
        (let ((r (modulo n 45045)))
          (if (not (vector-ref q63 (modulo r 63))) #f
            (if (not (vector-ref q65 (modulo r 65))) #f
              (if (not (vector-ref q11 (modulo r 11))) #f
                (let ((q (isqrt n)))
                  (if (= (* q q) n) q #f))))))))))

(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 (euclid x y)
  (let loop ((a 1) (b 0) (g x) (u 0) (v 1) (w y))
    (if (zero? w) (values a b g)
      (let ((q (quotient g w)))
        (loop u v w (- a (* q u)) (- b (* q v)) (- g (* q w)))))))

(define (inverse x m)
  (if (not (= (gcd x m) 1))
      (error 'inverse "divisor must be coprime to modulus")
      (call-with-values
        (lambda () (euclid x m))
        (lambda (a b g) (modulo a m)))))

(define (primes n) ; list of primes not exceeding n
  (let* ((len (quotient (- n 1) 2)) (bits (make-vector len #t)))
    (let loop ((i 0) (p 3) (ps (list 2)))
      (cond ((< n (* p p))
              (do ((i i (+ i 1)) (p p (+ p 2))
                   (ps ps (if (vector-ref bits i) (cons p ps) ps)))
                  ((= i len) (reverse ps))))
            ((vector-ref bits i)
              (do ((j (+ (* 2 i i) (* 6 i) 3) (+ j p)))
                  ((<= len j) (loop (+ i 1) (+ p 2) (cons p ps)))
                (vector-set! bits j #f)))
            (else (loop (+ i 1) (+ p 2) ps))))))

(define prime? ; strong pseudoprime to prime bases less than 100
  (let* ((ps (primes 100)) (p100 (apply * ps)))
    (lambda (n)
      (define (expm b e m)
        (let loop ((b b) (e e) (x 1))
          (if (zero? e) x
            (loop (modulo (* b b) m) (quotient e 2)
                  (if (odd? e) (modulo (* b x) m) x)))))
      (define (spsp? n a) ; #t if n is a strong pseudoprime base a
        (do ((d (- n 1) (/ d 2)) (s 0 (+ s 1)))
            ((odd? d) (if (= (expm a d n) 1) #t
              (do ((r 0 (+ r 1)))
                  ((or (= (expm a (* d (expt 2 r)) n) (- n 1)) (= r s))
                    (< r s)))))))
      (if (< n 2) #f (if (< 1 (gcd n p100)) (if (member n ps) #t #f)
        (do ((ps ps (cdr ps)))
            ((or (null? ps) (not (spsp? n (car ps)))) (null? ps))))))))

(define (next-prime n) ; smallest prime larger than n
  (define (wheel n)
    (vector-ref (vector 1 6 5 4 3 2 1 4 3 2 1 2 1 4
      3 2 1 2 1 4 3 2 1 6 5 4 3 2 1 2) (modulo n 30)))
  (if (< n 2) 2 (if (< n 3) 3 (if (< n 5) 5
    (let loop ((p (+ n (wheel n))))
      (if (prime? p) p (loop (+ p (wheel p)))))))))

(define (prev-prime n) ; largest prime smaller than n
  (define (wheel n)
    (vector-ref (vector 1 2 1 2 3 4 5 6 1 2 3 4 1 2
      1 2 3 4 1 2 1 2 3 4 1 2 3 4 5 6) (modulo n 30)))
  (if (<= n 2) #f (if (<= n 3) 2 (if (<= n 5) 3 (if (<= n 7) 5
    (let loop ((p (- n (wheel n))))
      (if (prime? p) p (loop (- p (wheel p))))))))))

(define random-prime ; random prime with n digits in base b
  (let* (;(seed (time-second (current-time))) ; chez
        ;(seed (current-seconds)) ; racket/chicken
        (seed (current-time)) ; guile
        ;(seed (inexact->exact (round (time->seconds (current-time))))) ; gambit
         (rand ; knuth linear congruential method
           (let* ((a 69069) (c 1234567) (m 4294967296))
             (lambda () (set! seed (modulo (+ (* a seed) c) m))
                        (/ seed m))))
         (randint (lambda (lo hi) (+ lo (floor (* (rand) (- hi lo)))))))
    (lambda (n . base)
      (let ((b (if (pair? base) (car base) 2)))
        (let loop ((p (randint 1 b)) (n (- n 1)))
          (if (zero? n) (prev-prime p)
            (loop (+ (* p b) (randint 0 b)) (- n 1))))))))

(define (factors n) ; fermat
  (if (even? n) (cons 2 (factors (/ n 2)))
    (let ((s (+ (isqrt n) 1)))
      (let loop ((u (+ s s 1)) (v 1) (r (- (* s s) n)))
        (cond ((positive? r) (loop u (+ v 2) (- r v)))
              ((negative? r) (loop (+ u 2) v (+ r u)))
              ((= (- u v) 2) (list (/ (+ u v -2) 2)))
              (else (append (factors (/ (+ u v -2) 2))
                            (factors (/ (- u v) 2)))))))))

(define (encipher str)
  (undigits
    (map (lambda (c) (- (char->integer c) 65))
         (string->list str))
    26))

(define (decipher x)
  (list->string
    (map (lambda (x) (integer->char (+ x 65)))
         (digits x 26))))

(define a 11314087)

(define (backdoor x-prime e)
  (let loop ((k x-prime))
    (let ((x (+ (* a k) x-prime)))
      (if (or (< 1 (gcd e (- x 1)))
              (not (= (modulo x 4) 3))
              (not (prime? x)))
          (loop (+ k 1))
          x))))

(define (keygen t*2 e)
  (let* ((t (quotient t*2 2))
         (len (quotient t 4))
         (p-prime (random-prime len 2))
         (q-prime (random-prime len 2))
         (p (backdoor p-prime e))
         (q (backdoor q-prime e))
         (n (* p q))
         (d (inverse e (* (- p 1) (- q 1)))))
    (values n d)))

(define (crypt text modulus key)
  (expm text key modulus))

(define t*2 64)
(define e 65537)
(define n #f)
(define d #f)
(call-with-values
  (lambda () (keygen t*2 e))
  (lambda (nn dd)
    (set! n nn)
    (set! d dd)))
    
(display n) (newline)
(display d) (newline)

(display (crypt (encipher "PRAXIS") n e)) (newline)
(display (decipher (crypt (crypt (encipher "PRAXIS") n e) n d))) (newline)

(display (factors (modulo n a))) (newline)
(display (backdoor 157 e)) (newline)
(display (backdoor 137 e)) (newline)
(display (* 2375958427 1787625883)) (newline)
(display (inverse e (* 2375958426 1787625882))) (newline)
(display (decipher (crypt 3852399687424634890
                          4247324781037166041
                          646136198260976393))) (newline)