搜索大质数的算法PHP版本【算法导论实践】
搜索大质数的算法PHP版本【算法导论实践】
产生大质数是现代公钥加密的基础（一般为RSA）。
在这里，我们需要搜索一个大质数（二进制，可能长达上千位！）必然会存在着很多问题。
1.质数有多少？
根据公式n/lgn，在n趋于无穷的情况下，这个式子的结果=数字n包含的质数数量。
虽然说是约等于，但是误差还是很小的（在原书中著有，10^9数量级的n误差为6%左右）
那么假设我要找一个长度为512位的质数的概率是多少呢？
[512/lg(2^512)-511/lg(2^511)]/(2^512-2^511) 必须注意，这个对数以2为底数。
原书中有歧义的是，原书中给出的概率应该包含0~512位的概率，所以这里给予修正。
如果你计算过的话，会发现概率不低，这就为我们搜索质数提供了可能。
我们可以找一个位数足够的随机数进行素数测试，测试几百个就应该会出现一个素数。
如果我们只测试奇数，那么上面得出的概率还会倍增！
2.产生一个随机数
对于2^n这样一个数量级的数字而言，用普通的算法来产生已经很困难了（在PHP中，整型的变量上限为2^31，这取决于cpu的寻址能力），对此我们可以用巧妙的方法迅速的产生一个大的随机数。
这个函数产生的随机数必定为一个奇数。
2.a^b mod n
解决了随机数的问题，我们接下来要解决的问题就是素数的测试。
通常我们所用的算法是
对待测数字n进行试除，一直到sqrt(n)的时候为止，若没有能使之模数字为0的数字，就说它是一个质数。
然后对于一个巨大的数字而言，这工作简直太可怕了，所以我们要另外寻找方法。
原书给出的通常算法是Miller-Rabin。在这里我不做介绍，这里给出它的前身。
首先给出求a^b mod n的方法。
对于巨大的数字求a^b 再对其模n几乎是不可能的（试想(2^n)^n)的增长率），但是我们由其公式可以得知
a*a mod n=a mod n * a mod n
这就对我们的运算产生了极大的方便。于是便有了反复平方法求解的思想。
3.判断质数
判断一个数是否为一个质数是很容易的。根据费马定理（page545）的逆命题（几乎和原命题同真）
我们可以得到一个测试素数的算法：
通过这个函数我们就可以判断这个数究竟是不是质数了。
必须注意！这个方法的判断是可能出现错误的。比如是合数，却被判断为了质数，是质数却被判断为合数。但错误的概率随着n增加，将变得微乎其微（Miller-Rabin正是对此的改进，它通过选用不同的基数来对其进行判断，但十分遗憾，对于所有的基数都有Carmichael数符合条件，也就是说凡是Carmilchael数都会被当做质数输出。幸运的是这种数字是非常稀有的（随着搜索的数字的数量级变大，概率逐渐降低，如果你找到的数字是这样一个数，那么恭喜你，你可以中500万了！）。
4.搜索质数
最后奉上搜索质数的主程序
*赠品：
输出的结果居然是二进制数字！这简直太难看了吧？（真正应用的时候应该没人看）
我们可以通过简单的转换把它转化为其他进制的数字。以下为几个进制转换的小程序。
在这里，我们需要搜索一个大质数（二进制，可能长达上千位！）必然会存在着很多问题。
1.质数有多少？
根据公式n/lgn，在n趋于无穷的情况下，这个式子的结果=数字n包含的质数数量。
虽然说是约等于，但是误差还是很小的（在原书中著有，10^9数量级的n误差为6%左右）
那么假设我要找一个长度为512位的质数的概率是多少呢？
[512/lg(2^512)-511/lg(2^511)]/(2^512-2^511) 必须注意，这个对数以2为底数。
原书中有歧义的是，原书中给出的概率应该包含0~512位的概率，所以这里给予修正。
如果你计算过的话，会发现概率不低，这就为我们搜索质数提供了可能。
我们可以找一个位数足够的随机数进行素数测试，测试几百个就应该会出现一个素数。
如果我们只测试奇数，那么上面得出的概率还会倍增！
2.产生一个随机数
对于2^n这样一个数量级的数字而言，用普通的算法来产生已经很困难了（在PHP中，整型的变量上限为2^31，这取决于cpu的寻址能力），对此我们可以用巧妙的方法迅速的产生一个大的随机数。
复制内容到剪贴板
这里提出要求，我们要产生的是一个二进制的随机数，那么$digit务必要求为16的倍数。PHP代码:
function random_int($digit)
{
$digit_pro=floor($digit/16);
for($i=1;$i <= $digit_pro;$i++)
{
$random.=decbin(mt_rand(0,65536));
}
return "1".decbin(mt_rand(0,16384)).$random."1";
}
这个函数产生的随机数必定为一个奇数。
2.a^b mod n
解决了随机数的问题，我们接下来要解决的问题就是素数的测试。
通常我们所用的算法是
对待测数字n进行试除，一直到sqrt(n)的时候为止，若没有能使之模数字为0的数字，就说它是一个质数。
然后对于一个巨大的数字而言，这工作简直太可怕了，所以我们要另外寻找方法。
原书给出的通常算法是Miller-Rabin。在这里我不做介绍，这里给出它的前身。
首先给出求a^b mod n的方法。
对于巨大的数字求a^b 再对其模n几乎是不可能的（试想(2^n)^n)的增长率），但是我们由其公式可以得知
a*a mod n=a mod n * a mod n
这就对我们的运算产生了极大的方便。于是便有了反复平方法求解的思想。
复制内容到剪贴板
上面给出了算法，里面的bc函数代表PHP的高精度数学运算函数，尽管说这样做避免了恐怖的运算，但是数字仍然很大，以至于我们不能直接对其进行普通数学运算。PHP代码:
//m-e
function modular_exponentiation($a,$b,$n)
{
$d=1;
for($k=strlen($b);$k > 0;$k--)
{
$d=bcmod(bcmul($d,$d),$n);
if(substr($b,0-$k,1)==1)
{
$d=bcmod(bcmul($d,$a),$n);
}
}
return $d;
}
3.判断质数
判断一个数是否为一个质数是很容易的。根据费马定理（page545）的逆命题（几乎和原命题同真）
我们可以得到一个测试素数的算法：
复制内容到剪贴板
这个过程用到了上面一节所写的函数，这个函数所进行的运算是整个算法中最为复杂的。PHP代码:
function pseudoprime($n)
{
$dec=bindec_pro($n);
if(bcmod(modular_exponentiation(2,decbin_pro(bcsub($dec,1)),$dec),$dec)!=1)
{
return false;
}
else
{
return true;
}
}
通过这个函数我们就可以判断这个数究竟是不是质数了。
必须注意！这个方法的判断是可能出现错误的。比如是合数，却被判断为了质数，是质数却被判断为合数。但错误的概率随着n增加，将变得微乎其微（Miller-Rabin正是对此的改进，它通过选用不同的基数来对其进行判断，但十分遗憾，对于所有的基数都有Carmichael数符合条件，也就是说凡是Carmilchael数都会被当做质数输出。幸运的是这种数字是非常稀有的（随着搜索的数字的数量级变大，概率逐渐降低，如果你找到的数字是这样一个数，那么恭喜你，你可以中500万了！）。
4.搜索质数
最后奉上搜索质数的主程序
复制内容到剪贴板
对于256位及以下的质数的搜索，该算法能在30s内找到（我尝试搜索了一个4096位的质数，在运气十分好的情况下（只测试了92个随机数），仍然用去了超过7000s的时间，这取决于CPU的运算能力以及语言的高效性）PHP代码:
//查找素数
function searchprime($digit=256)
{
$i=0;
while(true)
{
$num=random_int($digit);
if(pseudoprime($num))
{
return $num;
}
}
}
*赠品：
输出的结果居然是二进制数字！这简直太难看了吧？（真正应用的时候应该没人看）
我们可以通过简单的转换把它转化为其他进制的数字。以下为几个进制转换的小程序。
复制内容到剪贴板
PHP代码:
//10->2
function decbin_pro($num)
{
$b=0;
$i=0;
do
{
$q=bcdiv($num,2);
$r=bcmod($num,2);
$b=$r.$b;
$num=$q;
$i++;
}while($q!=0);
return $b;
}
//2->10
function bindec_pro($num)
{
$b=0;
$t=str_split($num);
$len=strlen($num)-1;
for($l=$len;$l >= 0;$l--)
{
$b=bcadd(bcmul($t[$len-$l],bcpow(2,$l)),$b);
}
return $b;
}
//2->16
function binhex_pro($num)
{
$arr=array(
"0000"=>"0",
"0001"=>"1",
"0010"=>"2",
"0011"=>"3",
"0100"=>"4",
"0101"=>"5",
"0110"=>"6",
"0111"=>"7",
"1000"=>"8",
"1001"=>"9",
"1010"=>"A",
"1011"=>"B",
"1100"=>"C",
"1101"=>"D",
"1110"=>"E",
"1111"=>"F");
$len=strlen($num);
$d=8-$len%8;
for($i=1;$i <= $d;$i++){
$num='0'.$num;
}
$t=str_split($num,8);
$p=ceil($len/8)-1;
for($i=0;$i <= $p;$i++)
{
$k=str_split($t[$i],4);
$b.=$arr[$k[0]].$arr[$k[1]]." ";
}
return $b;
}
