zoukankan      html  css  js  c++  java
  • Miller_Rabin大质数检验

    质数检验有不少算法,一般使用的质数检验复杂度是(O(sqrt{n}))

    又如线性筛可以在(O(n))的时间内求出所有1~n的质数

    但是,当n非常大,连(O(sqrt{n}))的复杂度也难以接受时,上述算法便不能满足要求

    这篇blog记录了一些关于Miller_Rabin算法的内容


    大家都知道的著名的费马小定理:

    [a^{p-1}equiv1pmod p ]

    其中(a,p)互质

    我们猜想,任意选取(a),如果一个数(p)满足以上式子,那么它就很有可能是一个质数

    但是这个猜想很容易找到反例:(a=2)(p=341,561, 645, 1105)时,费马小定理逆命题不成立,人们把(a^{p-1}equiv1pmod p)的合数称为以(a)为底的“伪素数”

    ([1,10^9])中,质数一共有(50847534)个,而满足(2^{p-1}equiv1pmod p)的合数(p)(5597)个,算法出错的概率太高了

    一个想法是,同时使用多个a来进行判断,例如同时检验(a=2,3)

    ([1,10^9])中,同时以(2,3)为底的伪素数只有(1272)个,不到以2为底的(frac{1}{4})

    选取的(a)越多时,算法越准确,这便是费马素性检验

    (以下内容引自Matrix67的博客)

    人们自然会想,如果考虑了所有小于n的底数a,出错的概率是否就可以降到0呢?没想到的是,居然就有这样的合数,它可以通过所有a的测试(这个说法不准确,详见我在地核楼层的回复)。Carmichael第一个发现这样极端的伪素数,他把它们称作Carmichael数。你一定会以为这样的数一定很大。错。第一个Carmichael数小得惊人,仅仅是一个三位数,561。前10亿个自然数中Carmichael数也有600个之多。Carmichael数的存在说明,我们还需要继续加强素性判断的算法。

    由此观之,费马素性检验仍不能满足我们的要求,我们需要改进素性检验算法


    引理:

    对于(forall a<p)(p)是质数

    (a^2equiv1pmod p)

    那么有(a=1)(a=p-1)

    证明:

    [a^2equiv1pmod p ]

    [Rightarrow a^2-1=(a+1)*(a-1)equiv0pmod p ]

    Miller_Rabin素数测试的方法是,对于待检验数(x),不断地提取(x-1)中的因子(2),把(x-1)表示成(2^r*d)的形式,那么我们需要计算的东西变成了(a^{2^r*d}mod x)

    于是,如果(x)是质数,(a^{2^r*d}mod x)要么等于(1),要么等于(x-1)

    如果(a^{2^{r-1}*d}mod x =1),定理同样适用于(a^{2^{r-2}*d}mod x)

    不断地开方,直到存在一个(iin[0,r)),使得(a^{2^i*d}mod x=x-1)

    至此,费马素性检验被强化为如下形式:

    尽可能提取因子(2),将待检验数(x)表示为(2^r*d)的形式,其中(d)是一个奇数
    如果(a^dmod n=1),或者找到一个(iin[0,r)),使得(a^{2^i*d}mod x=x-1),那么(x)极有可能是一个质数,否则(x)就不是一个质数

    上文提到极有可能,是因为Miller_Rabin同样是不确定算法,我们称能够通过以(a)为底的Miller_Rabin检验的合数称为以(a)为底的“强伪素数”,第一个以(2)为底的强伪素数为(2047),而第一个以(2)(3)为底的强伪素数则为1373653

    如果选用(2,3,7,61)(24251)作为(a),那么(10^{16})内唯一的强伪素数为(46space856space248space255space981),正确率可以接受


    Miller_Rabin核心代码

    这里,我采用先将n-1的所有因子2提取出来,再不断地平方回到原来的n-1的枚举方法

    typedef long long ll;
    const ll a[10]={0,2,3,7,61,24251};
    
    ll fastpow(ll a,ll b,ll p)
    {
    	ll re=1,base=a;
    	while(b)
    	{
    		if(b&1)
    			re=re*base%p;
    		base=base*base%p;
    		b>>=1;
    	}
    	return re;
    }
    
    bool Miller_Rabin(ll x,ll a)
    {
    	if(x==a)
    		return true;
    	if(x%a==0 || x<2)
    		return false;
    	ll d=x-1,r=0;
    	while(!(d&1))
    	{
    		d>>=1;
    		++r;
    	}
    	ll k=fastpow(a,d,x);
    	if(k==1 || k==x-1)
    		return true;
    	for(ll i=1;i<r;++i)
    	{
    		k=k*k%x;
    		if(k==x-1)
    			return true;
    	}
    	return false;
    }
    
    bool test(ll x)
    {
    	for(int i=1;i<=5;++i)
    		if(!Miller_Rabin(x,a[i]))
    			return false;
    	return true;
    }
    

    题外话

    由于水平及时间有限,博客中可能存在较多疏漏,希望读者能够指出,非常感谢

    撰写过程中,很多内容参考了matrix67的博客,在此表示感谢

    参考出处:http://www.matrix67.com/blog/archives/234

  • 相关阅读:
    maven的安装教程
    webstorm的中文教程和技巧分享
    WebStorm
    grunt配置任务
    grunt快速入门
    CSS简介
    浅介HTML DOM
    【转】计算机是如何启动的?
    【转】深入理解C++中public、protected及private用法
    【转】VS2013动态库文件的创建及其使用详解
  • 原文地址:https://www.cnblogs.com/lizbaka/p/Miller_Rabin.html
Copyright © 2011-2022 走看看