zoukankan      html  css  js  c++  java
  • Miller-Rabin素性测试

    有时候我们想快速的知道一个数是不是素数,而这个数又特别的大导致 $O(sqrt n)$ 的算法也难以通过,这时候我们可以对其进行 Miller-Rabin 素数测试,可以很大概率测出其是否为素数。

    两个理论基础

    (1)费马小定理:当 $p$ 为质数,有 $a^{p-1}equiv 1(mod p)$,反过来不一定成立,也就是说,如果 $a, p$ 互质,且 $a^{p-1}equiv 1(mod p)$,不能推出 $p$ 是质数,比如 $Carmichael$ 数

    (2)二次探测:如果 $p$ 是素数,$0 < x < p$,则方程 $x^2 equiv 1(mod p)$ 的解为 $x=1$ 或 $x=p-1$.

    证一下二次探测定理:

    易知 $ x^2-1 equiv 0(mod p)$

    $ecause(x+1)(x-1) equiv  0 (mod p)$

    $ecause$ $p$ 为素数,只能分解为

    $ecause x+1 equiv 0 (mod p)$ 或 $x-1 equiv 0 (mod p)$

    根据 $x$ 的取值范围,$x=1$ 或 $x=p-1$.

    Fermat 素性测试

    只根据费马小定理,易得一种检测质数的思路:

    它的基本思想就是多次从 $left [ 2, n-1 ight ]$ 中选取基 $a$,并检验是否每次都有 $a^{n-1} equiv 1 (mod n)$

    bool millerRabin(int n) {
      if (n < 3) return n == 2;
      // test_time 为测试次数,建议设为不小于 8
      // 的整数以保证正确率,但也不宜过大,否则会影响效率
      int  test_time = 20;
      for (int i = 1; i <= test_time; ++i) {
        int a = rand() % (n - 2) + 2;
        if (qpow(a, n - 1, n) != 1)
        {
            printf("a:%d
    ", a);
            return false;
        }
      }
      return true;
    }

    遗憾的是,由于费马小定理的逆定理并不成立 ,即满足  $a^{n-1} equiv 1 (mod n)$,$n$ 也不一定是素数。

    卡迈克尔数

    上面的做法中随机地选择 $a$,很大程度上降低了犯错概率,但是仍有一类数,上面的做法并不能准确地判断。

    对于合数,如果对于所有与 $n$ 互素的正整数 $a$,都有同余式 $a^{n-1} equiv 1 (mod n)$ 成立,则合数 $n$ 为卡迈克尔数(Carmichael Number),又称费马伪素数。

    而且我们知道,若 $n$ 为卡迈克尔数,则 $2^n-1$ 也是卡迈克尔数,从而卡迈克尔数的个数是无穷的。

    在1~100000000范围内的整数中,只有255个Carmichael数。前3个Carmichael数是561,1105,1729。

    Miller-Rabin测试

    根据卡迈克尔数的性质,可知其一定不是 $p^e$.

    不妨将费马小定理和二次探测定理结合起来使用:

    1. 对偶数和0、1、2直接判断
    2. 设要测试的数为 $n$,设$m$、$k$,满足 $n-1 = m2^k$(使其中 $m$ 为奇数)
    3. 我们先算出 $a^t$,然后不断地平方且进行二次探测(进行 $k$ 次)
    4. 最后根据费马小定理探测一次
    5. 多次取不同的 $a$ 进行上面3步

    可以证明Miller-Rabbin算法单次出错的概率小于等于 $frac{1}{4}$,若重复 $k$ 次,则出错概率降低至 ${(frac{1}{4})}^k$。这是一个很保守的估计,实际效果要好得多。

    #include<bits/stdc++.h>
    using namespace std;
    
    typedef long long int ll;
    
    ll mod_mul(ll a, ll b, ll mod)
    {
        ll res = 0;
        while (b)
        {
            if (b & 1)
                res = (res + a) % mod;
            a = (a + a) % mod;
            b >>= 1;
        }
        return res;
    }
    
    ll mod_pow(ll a, ll n, ll mod)
    {
        ll res = 1;
        while (n)
        {
            if (n & 1)
                res = mod_mul(res, a, mod);
            a = mod_mul(a, a, mod);
            n >>= 1;
        }
        return res;
    }
    
    // Miller-Rabin随机算法检测n是否为素数
    bool Miller_Rabin(ll n)
    {
        if (n == 2)
            return true;
        if (n < 2 || !(n & 1))
            return false;
        ll m = n - 1, k = 0;
        while (!(m & 1))
        {
            k++;
            m >>= 1;
        }
        for (int i = 1; i <= 10; i++)  // 20为Miller-Rabin测试的迭代次数
        {
            ll a = rand() % (n - 1) + 1;
            ll x = mod_pow(a, m, n);
            ll y;
            for (int j = 1; j <= k; j++)
            {
                y = mod_mul(x, x, n);
                if (y == 1 && x != 1 && x != n - 1)
                    return false;
                x = y;
            }
            if (y != 1)
                return false;
        }
        return true;
    }
    
    int main()
    {
        ll n;
        while(scanf("%lld", &n) == 1)
        {
            printf("%d
    ", Miller_Rabin(n));
        }
    }
    View Code

    注:

    (1)时间复杂度 $O(Slog^3)$($S$ 为测试次数)

    (2)long long相乘可能会爆long long,所以使用快速乘

    (3)当 $a$ 取遍30内的素数,$int$ 都不会出错

    (4)记住几个素数19260817、23456789、2092876369、11111111111111111111111(23个1)

     (5)还可以去 hihocodeCoder 1287 测一下板子

    参考链接:

    1. 百度百科——卡迈克尔数

    2. OI WIKI——素数

    3. https://blog.csdn.net/forever_dreams/article/details/82314237?tdsourcetag=s_pctim_aiomsg

    4. https://blog.csdn.net/ECNU_LZJ/article/details/72675595

    5. http://www.matrix67.com/blog/archives/234

  • 相关阅读:
    D
    hdu2376 Average distance (树形dp)
    hdu2376 Average distance (树形dp)
    选拔赛——旅游
    选拔赛——旅游
    cf 990c(思维+括号匹配)
    cf 990c(思维+括号匹配)
    Garland CodeForces
    Garland CodeForces
    Sherlock and his girlfriend CodeForces
  • 原文地址:https://www.cnblogs.com/lfri/p/11272199.html
Copyright © 2011-2022 走看看