Miller_Rabin素性测试学习小结

引入——威尔逊定理，费马小定理

相信大家都知道这个威尔逊定理定理吧，定理内容就是对于一个素数p，$(p-1)!\equiv -1\pmod p$,
这个东西的证明可以看我前面写的文章。
费马小定理：对于一个素数p，且a不整除p,那么
$a^{p-1}\equiv 1 \pmod p$,但是这定理的逆定理却不成立。
而威尔逊定理是判断一个数是否是素数的充分必要条件。
http://blog.csdn.net/ganjingxian/article/details/76268706，这里就不再赘述。利用这个东西我们可以推出另一个定理。

二次探测定理

对于一个素数n，方程$x^2\equiv 1 \pmod p$模p的根有两个，一个是n-1，另一个是1。我们把这两个根称为以n为模的平凡平方根。
结论：如果一个数存在1的非平凡平方根，则n是合数。

Miller_Rabin素性测试算法思路

①选取多个基数进行测试，（也就是说找多几个数进行测试）。
②寻找模n为1的非平凡平方根,做法如下：
先令$n-1=2^u\times t$,其中$u\ge 1,t为奇数$，
$a^{n-1}=a^{2^u \times t}=(a^t)^{2^u}$
先算出$a^t$，当然是用快速幂，然后再平方u次，如果在平方的过程中发现当前这个$(a^t)^{2^k}$模n等于1，那么n就不是素数。
③算完之后，看看算出来的这个值是不是1，若不是则这个数不是素数。

举个例子

举个例子，
$n=341，a=2，n-1=340＝2^2\times 85$,$x=2^{85}\equiv 32\pmod {341}$，$32^2\equiv 1\pmod {341}$,那么341肯定不是质数。

特别声明

这个测试的正确率并不是百分之百的，经过研究表明，它测试s次的出错概率最多为$2^{-s}$，但是一般来说这个算法的正确性很高，时间复杂度也不高。

代码

使用pascal的同志，对不起，我写的程序太烂了，测试23333333333，就215了。

var
        s,n:int64;
function mull(a,b,n:int64):int64;
var
        tmp:int64;
begin
        int64(tmp):=int64(double(a*b/n+1e-6))*n;
        exit(a*b-tmp);
end;
function f(a,b,n:int64):int64;
var
        t,y:int64;
begin
        t:=1;
        y:=a;
        while b<>0 do
        begin
                if(b and 1)=1 then
                t:=t*y mod n;
                y:=y*y mod n;
                b:=b shr 1;
        end;
        exit(t);
end;
function judge(n,a:int64):boolean;
var
        u,t,x,next,i:int64;
begin
        u:=0;
        t:=n-1;
        while t mod 2=0 do
        begin
                inc(u);
                t:=t div 2;
        end;
        x:=f(a,t,n);
        i:=1;
        while i<=u do
        begin
                next:=x*x mod n;
                if ((next=1)and(x<>1)and(x<>n-1)) then
                begin
                        exit(true);
                end;
                x:=next;
                inc(i);
        end;
        if x<>1 then exit(true);
        exit(false);
end;
function miller(n,s:int64):boolean;
var
        i:longint;
        a:int64;
begin
        if n=2 then
                exit(true);
        if (n=1) or (n mod 2=0) then exit(false);
        i:=1;
        while i<=s do
        begin
                a:=random(n-2)+2;
                if judge(n,a) then
                        exit(false);
                inc(i);
        end;
        exit(true);
end;
begin
        randomize;
        readln(n,s);
        if miller(n,s) then
                writeln('yes')
        else
                writeln('no');
end.

c++

#include<cstdio>
#include<ctime>
#include<cstdlib>
#define fo (i,a,b) for (int i=(a);i<=(b);i++)
#define fd (i,b,a) for (int i=(b);i>=(a);i--);
typedef  long long ll;
typedef  long double db;
ll n,s;
ll mult(ll a,ll b,ll n)
{
    ll temp=((ll)((db)a*b/n+1e-6)*n);
    return a*b-temp;
}
ll pow(ll a,ll b,ll n)
{
    ll t=1,y=a;
    while (b)
    {
        if (b&1)
            t=mult(t,y,n);
        y=mult(y,t,n);
        b=b>>1;
    }
    return t;
}
bool judge(ll n,ll a)
{
    ll u=0,t=n-1;
    while(t%2==0){u++;t=t/2;}
    ll x=pow(a,t,n);
    for(int i=1;i<=u;i++)
    {
        ll next=mult(x,x,n);
        if ((next==1)&&(x!=1)&&(x!=n-1))
        return true;
        x=next;
    }
    if (x!=1)return true;
    return false;
}
bool Miller_Rabin(ll n,int s)
{
    if (n==2)return true;
    if (n<2 || n % 2==0)return false;
    for(int i=1;i<=s;i++)
    {
        ll a=rand()%(n-2)+2;
        if (judge(n,a))return false;
    }
    return true;
}
int main()
{
    srand(time(0));
    scanf("%lld%d",&n,&s);
    if (Miller_Rabin(n,s))
        printf("yes");
    else
        printf("no");
    return 0;
}

c++就不怕会像pascal一样215。我觉得应该是我写的代码太丑了。

最后补充

上面的代码中有一个黑科技，就是模的时候将它转化一下，变成mult（pascal中为mull），原理：$a \pmod b=a-[\frac{a}{b}]\times b$,这个东西的证明十分简单，证法如下：
设a=bk+r。

$$\begin{align} a\pmod b & = r \\ & = a-[\frac{bk+r}{b}]\times b \\ & = a-bk \\ & = r \end{align}$$

the end

由于我的水平有限，难免会有些写错的地方，希望大家批评指正，多多包容，thank you for your patience.