Miller-Rabin 方案实现不可预测的输出

Question

我是 Scheme 的新手。我已经尝试并使用 PLT Scheme 实现了 Rabin-Miller 算法的概率变体。我知道这是概率性的，但我大部分时间都得到错误的结果。我用 C 实现了同样的东西，而且效果很好（从来没有失败过）。我在调试时得到了预期的输出，但是当我运行时，它几乎总是返回错误的结果。我使用了Wikipedia的算法。

(define expmod( lambda(b e m)
                 ;(define result 1)
                 (define r 1)
                 (let loop()
                   (if (bitwise-and e 1)
                       (set! r (remainder (* r b) m)))
                   (set! e (arithmetic-shift e -1))
                   (set! b (remainder (* b b) m))
                   (if (> e 0)
                       (loop)))r))

(define rab_mil( lambda(n k)
                  (call/cc (lambda(breakout)
                  (define s 0)
                  (define d 0)
                  (define a 0)
                  (define n1 (- n 1))
                  (define x 0)          
                  (let loop((count 0))
                    (if (=(remainder n1 2) 0)
                        (begin
                          (set! count (+ count 1))
                          (set! s count)
                          (set! n1 (/ n1 2))
                          (loop count))
                        (set! d n1)))
                  (let loop((count k))
                    (set! a (random (- n 3)))
                    (set! a (+ a 2))
                    (set! x (expmod a d n))
                    (set! count (- count 1))

                    (if (or (= x 1) (= x (- n 1)))
                        (begin
                          (if (> count 0)(loop count))))
                    (let innerloop((r 0))
                      (set! r (+ r 1))
                      (if (< r (- s 1)) (innerloop  r))
                      (set! x (expmod x 2 n))
                      (if (= x 1)
                          (begin
                          (breakout #f)))
                      (if (= x (- n 1)) 
                          (if (> count 0)(loop count)))
                      )
                    (if (= x (- s 1)) 
                        (breakout #f))(if (> count 0) (loop count)))#t))))

另外，我在 Scheme 中的编程方式是否正确？ （我不确定我使用 call/cc 的循环部分的中断。我在某个网站上找到了它，从那以后一直在使用它。）

提前致谢。

Answer 1

一般来说，您正在以一种过于“命令式”的方式进行编程；一个更优雅的 expmod 将是

(define (expmod b e m)
  (define (emod b e)
    (case ((= e 1) (remainder b m))
          ((= (remainder e 2) 1)
           (remainder (* b (emod b (- e 1))) m)
          (else (emod (remainder (* b b) m) (/ e 2)))))))
  (emod b e))

避免使用 set!并且只是递归地实现规则

b^1 == b (mod m)     
b^k == b b^(k-1) (mod m) [k odd]
b^(2k) == (b^2)^k (mod m)

类似地，rab_mil 的东西是以非常非方案的方式编程的。这是一个替代实现。请注意，循环没有“中断”，也没有 call/cc；相反，突破是作为尾递归调用实现的，它真正对应于 Scheme 中的“goto”：

(define (rab_mil n k)
  ;; calculate the number 2 appears as factor of 'n'
  (define (twos-powers n)
     (if (= (remainder n 2) 0)
         (+ 1 (twos-powers (/ n 2)))
         0))
  ;; factor n to 2^s * d where d is odd:
  (let* ((s (twos-powers n 0))
         (d (/ n (expt 2 s))))
    ;; outer loop
    (define (loop k)
      (define (next) (loop (- k 1)))
      (if (= k 0) 'probably-prime
          (let* ((a (+ 2 (random (- n 2))))
                 (x (expmod a d n)))
            (if (or (= x 1) (= x (- n 1)))
                (next)
                (inner x next))))))
    ;; inner loop
    (define (inner x next)
      (define (i r x)
        (if (= r s) (next)
            (let ((x (expmod x 2 n)))
              (case ((= x 1) 'composite)
                    ((= x (- n 1)) (next))
                    (else (i (+ 1 r) x))))
      (i 1 x))
    ;; run the algorithm
    (loop k)))