「Pollard-Rho」学习笔记

发表于 | 分类于 | 8

Pollard-Rho 是用来分解大数质因子的一个玩意。

Miller-Rabin 素性判断

众所周知，费马告诉我们对于一个素数 p，有 ap−1≡1(modp)。
然鹅它的逆命题是不成立的，即 ap−1≡1(modp) 不是 p 的充分条件。

不过，大部分合数都可以被这个式子筛掉，所以我们可以基于这个式子进行下一步改进。
然后有二次探测定理：对于素数 p，如果有 a2≡1(modp)，那么有 a≡1(modp) 或 a≡−1(modp)。

简单地证明一下：
x2=(x+1)(x−1)+1。
所以 (x+1)(x−1)≡0(modp)。
要么 x+1≡0(modp)，要么 x−1≡0(modp)。
然后易得那两个式子。

考虑到费马小定理和二次探测定理的形式很像，并且 p−1 一定为偶数，所以我们可以一直运用二次探测定理。
对于不同的 a，我们把上述算法称之为 以 a 为底的 Miller-Rabin 测试。

一般底数仍然是随机选取，但当待测数不太大时，选择测试底数就有一些技巧了。
比如，如果被测数小于 4759123141，那么只需要测试三个底数 2,7,61 就足够了。
当然，你测试的越多，正确的范围肯定也越大。如果你每次都用前 7 个素数 (2,3,5,7,11,13,17) 进行测试，所有不超过 341550071728320 的数都是正确的。
如果选用 2,3,7,61,24251 作为底数，那么 1016 内唯一的强伪素数为 46856248255981。
这样的一些结论使得 Miller-Rabin 算法在OI中非常实用。
——Matrix67

根据本人测试，对于上面那个 46856248255981，增加底数 19∗∗∗∗17 即可。
（大雾）

代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
//洛谷 3383.「模板」线性筛素数
#include <cstdio>
using namespace std;
int n,m;
inline int fpow(int a,int b,int mod)
{
    int ret = 1;
    for(;b;b >>= 1)
        (b & 1) && (ret = (long long)ret * a % mod),a = (long long)a * a % mod;
    return ret;
}
inline bool witness(int x,int p)
{
    if(fpow(x,p - 1,p) != 1)
        return 0;
    int k = p - 1;
    while(!(k & 1))
    {
        k >>= 1;
        int t = fpow(x,k,p);
        if((t ^ 1) && (t ^ p - 1))
            return 0;
        if(t == p - 1)
            return 1;
    }
    return 1;
}
inline bool Miller_Rabin(int n)
{
    const int base[] = {2,3,7,61,24251,19260817};
    if(n == 1)
        return 0;
    for(register int i = 0;i < 6;++i)
        if(n == base[i])
            return 1;
        else if(!witness(base[i],n) || !(n % base[i]))
            return 0;
    return 1;
}
int main()
{
    scanf("%d%d",&n,&m);
    while(m--)
    {
        scanf("%d",&n);
        puts(Miller_Rabin(n) ? "Yes" : "No");
    }
}

Pollard-Rho

考虑对于一个合数 n，它的一个质因子为 p，如果有 a≡b(modp) 但 a≢b(modn)，那么 |a−b| 便是 n 的一个非平凡因子。
所以一个比较暴力的方法就是一直随机两个数并计算它们的差与 n 的 gcd 然后递归做。

然鹅这样子考虑到随机数的效率和质量，并没有什么卵用……
我们考虑手写随机数，设我们生成的第 i 个随机数是 xi，并选取一个常数 c，有 xi≡x2i−1+c(modn)。
但是接下来又有一个问题，怎样提高准确率？
显然地，x 数列会出现循环。如果 xi≡xj(modd) (i<j)，那么数列在模 d 意义下是周期的，其长度整除 j−i。
那么对于 j−i 最小的 >i 的倍数 s，有 xs≡x2s(modd)。
所以我们可以对于 xk 与 x2k 进行这个过程。

又有一个问题，如果判环？
哈希？太麻烦了！
考虑在一个环上，甲、乙从同一个点出发，甲经过一条边的时间乙能经过两条边，甲、乙向同一个方向走。
那么当乙追上甲时，甲一定至少走了一圈。
这个方式恰好可以与上面的过程结合起来。

在实际应用中，一定要注意先用线性筛出空间允许的范围内的素数并试除，最后再上 Pollard-Rho。
因为它只能找到部分较大的质因子。

代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
//洛谷 4718.「模板」Pollard-Rho 算法
#include <cstdio>
#include <cmath>
#include <algorithm>
using namespace std;
long long ans;
inline long long fmul(long long a,long long b,long long mod)
{
    return (a * b - (long long)((long double)a / mod * b) * mod + mod) % mod;
}
inline long long fpow(long long a,long long b,long long mod)
{
    long long ret = 1;
    for(;b;b >>= 1)
        (b & 1) && (ret = fmul(ret,a,mod)),a = fmul(a,a,mod);
    return ret;
}
inline bool witness(long long x,long long p)
{
    if(fpow(x,p - 1,p) ^ 1)
        return 0;
    long long k = p - 1;
    while(!(k & 1))
    {
        k >>= 1;
        long long t = fpow(x,k,p);
        if((t ^ 1) && (t ^ p - 1))
            return 0;
        if(t == p - 1)
            return 1;
    }
    return 1;
}
inline bool Miller_Rabin(long long n)
{
    if(n == 1)
        return 0;
    const int base[] = {2,3,7,61,24251,19260817};
    for(register int i = 0;i < 6;++i)
        if(n == base[i])
            return 1;
        else if(!witness(base[i],n))
            return 0;
    return 1;
}
void Pollard_Rho(long long n)
{
    if(n <= ans)
        return ;
    if(Miller_Rabin(n))
    {
        ans = max(ans,n);
        return ;
    }
    const long long c = 20060312;
    const long long k = 127;
    long long t1 = 20070119;
    long long t2 = (fmul(t1,t1,n) + c) % n;
    long long i = 0,p = 1;
    while(t1 ^ t2)
    {
        ++i,p = fmul(p,abs(t1 - t2),n);
        if(!p)
        {
            long long d = __gcd(abs(t1 - t2),n);
            if((d ^ 1) && (d ^ n))
            {
                Pollard_Rho(d),Pollard_Rho(n / d);
                return ;
            }
        }
        if(!(i % k))
        {
            long long d = __gcd(p,n);
            p = 1;
            if((d ^ 1) && (d ^ n))
            {
                Pollard_Rho(d),Pollard_Rho(n / d);
                return ;
            }
        }
        t1 = (fmul(t1,t1,n) + c) % n,t2 = (fmul(t2,t2,n) + c) % n,t2 = (fmul(t2,t2,n) + c) % n;
    }
    long long d = __gcd(p,n);
    if((d ^ 1) && (d ^ n))
        Pollard_Rho(d),Pollard_Rho(n / d);
}
int T;
long long n;
int main()
{
    scanf("%d",&T);
    while(T--)
    {
        ans = 0;
        scanf("%lld",&n);
        Pollard_Rho(n);
        ans == n ? puts("Prime") : printf("%lld\n",ans);
    }
}

Related Issues not found

Please contact @Alpha1022 to initialize the comment