python 费马检测

jokerchu

　　费马小定理：
　　如果 n 是一个素数，a 是小于 n 的任意正整数，那么 a 的 n 次方与 a 模 n 同余。（两个数称为是模 n 的同余，如果它们除以 n 的余数相同。数 a 除以 n 的余数称为 a 取模 n 的余数，或简称为 a 取模 n）。
如果 n 不是素数，那么，一般而言，大部分的 a < n 都将满足上面的关系。这就引出了下面这个检查素数的算法：
对于给定的整数 n，随机任取一个 a < n 并计算出 an 取模 n 的余数。如果得到的结果不等于 a，那么 n 就肯定不是素数。如果它就是 a，那么 n 是素数的机会就很大。现在再另取一个随机的 a 并采用同样的方式检查。如果它满足上述等式，那么我们就能对 n 是素数有更大的信心了。通过检查越来越多的 a 值，我们就可以不断增加对有关结果的信心。这一算法称为费马检查。
读起来理解不直观？那么我这么总结下吧。
假如a是一个整数，p是一个素数，那么 ap = a (mod p)
如果a不是p的倍数，这个定理也可以写成 ap-1 = 1 (mod p)
举个例子，67是一个素数，则266 mod 67 = 1。

import random
def expmod(base,exp,m):
if exp==0:
return 1
elif exp%2==0:
return expmod(base,exp/2,m)**2 % m
else:
return expmod(base,exp-1,m)*base % m
def fermat_test(n):
a=random.randint(1,n-1)
if a==expmod(a,n,n):
return 1
else:
return 0
def fermat_prime(n,time):
if time==0:
return 1
elif fermat_test(n)==1:
return fermat_prime(n,time-1)
else:
return 0
#print fermat_test(11)
print fermat_prime(99,100)
　　对于为啥不直接用 a^exp 来求幂值 解释如下：
　　【转】http://www.blogjava.net/killme2008/archive/2007/05/11/116813.html
　　直接定义(expmod base exp m)为base^exp基于m的模(modulo)
　　(define (expmod base exp m)
　　(remainder (fast-expt base exp) m))
　　在理想情况下是正确的，在语义上与原定义是等价的，甚至递归层数
　　都是一样的，而且更加直观。
　　但在实践中，这样的定义是不可行的，这也是为什么我们要采用原文中的定义
　　的原因：因为base^exp很可能是一个非常大的数。比如，当我们取第二个
　　测试例子中的n=1000000时，base是[1,999999]区间中的任意
　　随机数，它的平均取值为50000，而exp=1000000。50000^1000000是一个大得
　　惊人的数，无论是fast-expt的层层函数调用计算，还是用remainder对其取模，
　　计算量都是很大的。
　　而相反，原始的expmod函数
　　(define (expmod base exp m)
　　(cond ((= exp 0) 1)
　　((even? exp)
　　(remainder (square (expmod base (/ exp 2) m))
　　m))
　　(else
　　(remainder (* base (expmod base (- exp 1) m))
　　m))))
　　通过分解(a*b mod n) 为 ((a mod n) * (b mod n) mod n)，使得每层递归
　　调用的计算参数和返回值总是小于n (x mod n < n)，即便是计算过程中出现
　　的最大数(a mod n) * (b mod n) 也依然是要小于n^2, 相对于前者n^n的数
　　量级，实在是小得多，这样就避免了大数的计算问题。
　　比如经测试，在使用新的expmod的情况下，1009的验证需要1.2ms的时间，
　　1000003的验证需要 93680ms，远高于原来的expmod函数。
　　这也同时印证了我们在1.24题中的分析，同样的操作在不同字长的数字上的
　　代价是不同的。1000000的验证时间现在是1000的8000多倍，远不再是2倍左右
　　的关系。在1.24中的，因为expmod算法的控制，操作数最多是n^2级的，字长
　　所引起的差距并不明显，只在10^6-10^12间产生了一点点阶跃；而这里的算法
　　中, 操作数可以达到n^n级，当n=1000时，1000^1000的字长大约在10000位(二
　　进制数)左右，而n=1000000时1000000^1000000的字长则达到在20000000位(二
　　进制数)左右，这字长的巨大差距导致了我们在1.24中已经发现的问题更加明显。