首页 > 代码库 > 素数检测杂谈

素数检测杂谈

2024-08-30 20:42:42   222人阅读

最开始我们学到的是这种朴素的写法：

 

(define (prime? n)
  (cond
   [(or (= n 2) (= n 3)) #t]
   [(even? n) #f]
   [else (prime-test-loop n)]))

(define (prime-test-loop n)
  (let ((top (ceiling (sqrt n))))
    (let iter ((start 3))
      (cond [(> start top) #t]
            [(= (mod n start) 0) #f]
            [else (iter (+ start 2))]))))

 

虽然很笨拙，但其实这也是经过“优化”的，首先，除了 2 以外的偶数就不用判断了。其次，试除迭代从 3 起步，每次 +2，同样避开了偶数。最后，循环的结束条件是 (sqrt n). 

后来又看到一种生成素数序列的方法，大体思想是维护一个已知的素数表，每一个新的数字 n 都用已知素数表里的数去试除。如果都不能整除，说明 n 是素数，将 n 添加到已知的素数表里，进入下一轮迭代。原文是用 C 实现的，事先搞一个大数组，每一轮迭代都遍历数组，将已知素数的整数倍所对应的数组下标都划掉，最后留在数组中的就是素数。下面是我写的 Scheme 版：

 

(define (prime-list n)
  (let iter ((start 3) (primes ‘(2)))
    (cond [(> start n) primes]
          [(andmap (lambda (p) (not (= (mod start p) 0)))
                   primes)
           (iter (+ start 2) (cons start primes))]
          [else (iter (+ start 2) primes)])))

 

后来知道有一种基于概率的检测方法叫费马检查，SICP 上就有一个 Scheme 实现：

(define (square n) (* n n))

(define (expmod base exp m)
  (cond ((= exp 0) 1)
        ((even? exp)
         (mod (square (expmod base (/ exp 2) m))
              m))
        (else
         (mod (* base (expmod base (- exp 1) m))
              m))))

(define (fermat-test n)
  (define (try-it a)
    (= (expmod a n n) a))
  (try-it (+ 1 (random (- n 1)))))

(define (fast-prime? n times)
  (cond ((= times 0) #t)
        ((fermat-test n) (fast-prime? n (- times 1)))
        (else #f)))

这个算法很不可靠，因为它所依据的理论是个伪命题，费马小定理断言：如果 n 是一个素数，那么小于 n 的任意正整数 a 的 n 次方再对 n 取模，结果依然为 a。

它的逆命题是，a 是任意一个小于 n 的正整数,  如果 a ^ n % n == a，那么 n 就是一个素数。很遗憾，这个逆命题不成立，因为有一类数叫做 Carmichael 数，它符合  a ^ n % n == a 这个特性，但却不是素数。直接用这段程序检查一个数是不是素数有很大的机会出错。这是已知 Carmichael 数的一部分：

561, 1105, 1729, 2465, 2821, 6601, 8911, 10585, 15841, 29341, 41041, 46657, 52633, 62745, 63973, 75361, 101101, 115921, 126217, 162401, 172081, 188461, 252601, 278545, 294409, 314821, 334153, 340561, 399001, 410041, 449065, 488881, 512461

可以看到，最小的数是561。如果你准备写一个程序，生成一个素数表，然而它迭代了500多次以后就开始胡言乱语，那玩笑就有点开大了。SICP 提到一个费马检测的变形，叫做 Miller-Rabin 检测，它能够避免 Carmichael 数的愚弄。

这个检测依据的断言是：如果 n 是一个素数，a 是任意一个小于 n 的正整数，那么 a ^ (n-1) % n == 1 % n。

要用 米勒-拉宾 检测算法测试一个数字是不是素数，我们应当选择一个小于 n 的随机正整数 a, 然后利用 expmod 函数计算 a ^ (n-1) % n . 但是当执行到 expmod 里的 square 那一步时，我们检测是不是找到了一个 (1 % n) 的非平凡的平方根，换言之，看一看是不是有这样一个数字，它既不是 1, 也不是 n-1，但是它的平方等于 1 % n。很容易证明，如果这样一个 1 的非平凡的平方根存在，那么 n 肯定不是素数。同样可以证明，如果 n 是一个并非素数的奇数，在小于 n 的正整数中，至少有一半在用这种方法计算 a ^ (n-1) 时能够找到这种非凡平方根。

中文版的 SICP 看不太明白，我自己翻译英文题目，还是觉得看不太明白。但大概意思明白了，如果一个小于 n 的正整数 a，它符合 a ^ (n-1) & n == 1，但是 a 既不是 1, 也不是 n-1，那么这个 n 肯定不是素数。对一个伪素数来说，只做一次测试，正确与错误的概率对半分，基本和算命差不多。如果做两次测试，那就很有可能发现它的真面目。重复的次数越多，结果越可信。如果重复多次依然通过了测试，那个这个数 n 大概真的是个素数了。这就是它能避免 Carmicheal 数愚弄的原因。

 

SICP 留了一个作业 1.28 ，要求实现 Miller-Rabin 检测，下面就是实现：

(define (expmod base exp m)
  (cond ((= exp 0) 1)
        ((even? exp)
         (check-nontrivial-sqrt (expmod base (/ exp 2) m) m)) ;; look here
        (else
         (mod (* base (expmod base (- exp 1) m)) m))))

(define (check-nontrivial-sqrt n m)
  (let ((x (mod (square n) m)))
    (if (and (not (= n 1))
             (not (= n (- m 1)))
             (= x 1))
        0
        x)))

(define (miller-rabin-test n)
  (define (try-it a)
    (= (expmod a (- n 1) n) 1))
  (try-it (+ 1 (random (- n 1)))))

(define (fast-prime? n times)
  (cond ((= times 0) #t)
        ((miller-rabin-test n) (fast-prime? n (- times 1)))
        (else #f)))

;; 搞了个包装函数，偶数就不必判断了吧
(define (prime? n)
  (cond ((= n 2) #t)
        ((even? n) #f)
        (else (fast-prime? n 20))))

好吧，写到是写出来了，但是测试发现，它依然胡言乱语。。。

等有空再调试一下。

素数检测杂谈


联系
我们