【问题标题】:Bug In Miller-Rabin ImplementationMiller-Rabin 实现中的错误
【发布时间】:2012-08-03 17:36:14
【问题描述】:

我正在实施Wikipedia's Miller-Rabin algorithm,但似乎并没有得到非常恰当的结果。 7, 11, 19, 23 等被报道为复合材料。事实上,当 k>12 时,即使 5 也显示为复合。我已经阅读了 Miller-Rabin 背后的数学,但不太了解它,并且盲目地依赖算法。关于我哪里出错的任何线索?

这是我的代码:

#include<stdio.h>
#include<math.h>

int modpow(int b, int e, int m) {
    long result = 1;

    while (e > 0) {
        if ((e & 1) == 1) {
            result = (result * b) % m;
        }
        b = (b * b) % m;
        e >>= 1;
    }

    return result;
}

int isPrime(long n,int k){
        int a,s,d,r,i,x,loop;
        if(n<2)return 0;
        if(n==2||n==3)return 1;
        if(n%2==0)return 0;

        d=n-1;
        s=0;
        while(d&1==0){
                d>>=1;
                s++;
        }


        for(i=0;i<k;i++){
                loop=0;
                a=(int)(rand()*(n-1))+1;
                x=modpow(a,d,n);
                if(x==1 || x==n-1){
                        continue;

                }
                for(r=1;r<=s;r++){
                        x=modpow(x,2,n);
                        if(x==1)return 0;
                        if(x==n-1){
                                loop=1;
                                break;
                        }
                }
                if(!loop)return 0;

        }
        return 1; 

}

int main(){
        int i,k;
        scanf("%d",&k);
        for(i=5;i<100;i+=2){
                printf("%d : %d\n",i,isPrime(i,k));
        }
        return 0;
}

【问题讨论】:

  • 可惜我现在不在校园里,否则我会让米勒自己看一看。
  • 你确定你没有读错输出吗?您究竟期望输出是什么以及您得到什么?
  • 我看到至少一个问题;我会很快发布答案。
  • @DennisMeng :太棒了。但我想用实现错误来打扰他已经太多了。 :)
  • @another.anon.coward :我应该得到所有(可能)质数的 1 和复合数的 0。这就是我得到的:ideone.com/bYDjV

标签: c primes


【解决方案1】:

如果碱基与候选者不互质,则强费马检查总是返回“不是可能的素数”。

你的错误

a=(int)(rand()*(n-1))+1;

对于素数 p,当且仅当 rand() 的结果具有以下形式时,底数不是 pp 的倍数)的互质数

k*p + 1

对于小素数,即使迭代很少,也几乎可以保证发生这种情况。

基数应介于 2 和 n/2 之间(选择大于 n/2 的基数是不必要的,因为当且仅当 n - a 是一个时,a 是复合性的见证),所以你想要类似

a = rand() % (n/2 - 2) + 2;

如果您不介意随机数生成中偏向小余数的模偏差,或者

a = rand() /(RAND_MAX + 1.0) * (n/2 - 2) + 2;

如果您想在整个可能范围内分配偏差。

【讨论】:

  • 谢谢!解决了这个问题。使用您建议的最后一种形式。由于我会将 k 保持在相当低的水平,因此最好在整个范围内进行分布。
  • @Adarsh a 也以第一种形式分布在整个范围内。这只是由以下事实引入的偏差,通常人们感兴趣的数字范围不是RAND_MAX + 1 的除数。 rand()N = RAND_MAX + 1 可能的结果。如果您想要m 可能的数字之一,则必须将N 可能的数字尽可能均匀地分布在m 存储桶上。除非mN 的除数,否则某些桶必须包含比其他桶多一个数字。第一个公式使存储桶 0 到 RAND_MAX % m(含)成为多一个数字的存储桶,另一种方法
  • 以难以预测的方式在整个范围内分配更完整的桶。如果m 相对于N 较小,则通常可以忽略分布的偏斜,但如果m 较大,则意味着一些桶包含一个数字,而另一些桶包含两个,或者如果@,则甚至一些1 和其他0 987654346@。在后一种情况下,应该采取措施消除偏差(对于m &lt; NN = k*m + r,只能接受小于k*mrand() 的结果,如果结果是,则再次调用rand()被拒绝。
  • 感谢您的解释。在我可以在这里应用任何概率分布之前,我需要复习我的数论。目前,我只有一个模糊的想法。我可以清楚地看到使用 rand()/(1.0 + RAND_MAX) 的不同之处——显然更统一,但你评论的最后一点需要我做一些修改。我会尽快解决的。谢谢。 :)
猜你喜欢
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 2011-01-15
  • 1970-01-01
  • 2016-03-19
  • 2012-06-29
  • 2019-10-04
相关资源
最近更新 更多