标签:题解 unsigned lin href srand 讲解 检验 表示 定义
传送门:http://www.lydsy.com/JudgeOnline/problem.php?id=4802
【题解】
参考:http://www.matrix67.com/blog/archives/234
Millar-Rabin质数检验方法:
根据费马小定理,如果p是素数,a<p,那么有a^(p-1) mod p = 1。
直观想法我们直接取若干个a,如果都有一个不满足,那么p就是合数。
遗憾的是,存在Carmichael数:你无论取多少个a,有一个不满足,算我输。
比如:561 = 11*51就是一个Carmichael数。
那么就很江了啊。。我们需要改进算法。
首先有:如果p是素数,x是小于p的正整数,且x^2 mod p = 1,那么要么x=1,要么x=p-1
(这个废话,x=p-1模意义下等于x=-1)
然后我们可以展示下341满足2^340 mod 341 = 1,却不是素数(341=31*11)的原因:
2^340 mod 341 = 1
2^170 mod 341 = 1
2^85 mod 341 = 32
(32这个数很py啊怎么不等于340也不等于1啊。。这明显有交易嘛)
那么就能说明这个数不是素数。
如果是素数,一定是从p-1变到1,或是把所有2的次幂去除完,本来就等于1(这样平方完就一直是1了)
所以要么把所有2的次幂去除完,本来就等于1,要么存在某一个次幂=p-1(这样就正常多了)
这就是Millar-Rabin素数验证的二次探测。
应该来说Millar-Rabin算法也是挺好写的
其中mul(a,b,c)表示a*b%c(因为a*b会爆longlong,所以用快速加)
namespace Millar_Rabin { const int Prime[14] = {0, 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41}; const int PN = 13; inline bool witness(int pr, ll res, int times, ll n) { ll p = pwr2((ll)pr, res, n); if(p == 1) return 0; for (int i=0; i<times; ++i) { if(p == n-1) return false; if(p == 1) return false; p = mul(p, p, n); } return true; } inline bool main(ll n) { for (int i=1; i<=PN; ++i) { if(n == Prime[i]) return 1; if(n % Prime[i] == 0) return 0; } ll p = n-1; int times = 0; while(!(p&1)) { ++times; p >>= 1; } for (int i=1; i<=PN; ++i) if(witness(Prime[i], p, times, n)) return false; return true; } }
然后我们会检验素数了,现在要质因数分解。
好了下一个是Pollard-Rho算法:
如果现在拆分的是n:Pollard-Rho(n)
主要流程:Millar-Rabin判断是否质数,是返回,否就试图找出其中一个因子d,然后递归做Pollard-Rho(d)和Pollard-Rho(n/d)。
我猜你会说废话这谁都会。问题在于:试图找出其中一个因子d
参考:https://wenku.baidu.com/view/3db5c7a6ad51f01dc381f156.html?from=search
参考文章讲的非常详细了。。我就不细讲了qwq
所以这题就是分解因数,按照欧拉函数定义式求解即可。
# include <stdio.h> # include <string.h> # include <iostream> # include <algorithm> // # include <bits/stdc++.h> using namespace std; typedef long long ll; typedef long double ld; typedef unsigned long long ull; const int M = 5e5 + 10; const int mod = 1e9+7; # define RG register # define ST static inline ll mul(ll a, ll b, ll mod) { ll ret = 0; a %= mod, b %= mod; while(b) { if(b&1) { ret = ret + a; if(ret >= mod) ret -= mod; } a <<= 1; if(a >= mod) a -= mod; b >>= 1; } return ret; } inline ll pwr2(ll a, ll b, ll mod) { ll ret = 1; a %= mod; while(b) { if(b&1) ret = mul(ret, a, mod); a = mul(a, a, mod); b >>= 1; } return ret; } inline ll gcd(ll a, ll b) { return b==0 ? a : gcd(b, a%b); } namespace Millar_Rabin { const int Prime[14] = {0, 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41}; const int PN = 13; inline bool witness(int pr, ll res, int times, ll n) { ll p = pwr2((ll)pr, res, n); if(p == 1) return 0; for (int i=0; i<times; ++i) { if(p == n-1) return false; if(p == 1) return false; p = mul(p, p, n); } return true; } inline bool main(ll n) { for (int i=1; i<=PN; ++i) { if(n == Prime[i]) return 1; if(n % Prime[i] == 0) return 0; } ll p = n-1; int times = 0; while(!(p&1)) { ++times; p >>= 1; } for (int i=1; i<=PN; ++i) if(witness(Prime[i], p, times, n)) return false; return true; } } namespace PollardRho { const int N = 110; ll q[N]; int qn; inline void PR(ll n) { if(Millar_Rabin::main(n)) { q[++qn] = n; return ; } ll a, b, c, del; while(1) { c = rand() % n; a = b = rand() % n; b = (mul(b, b, n) + c) % n; while(a != b) { del = a-b; del = gcd(abs(del), n); if(del > 1 && del < n) { PR(del); PR(n/del); return ; } a = (mul(a, a, n) + c) % n; b = (mul(b, b, n) + c) % n; b = (mul(b, b, n) + c) % n; } } } inline ll getphi(ll n) { if(n == 1) return 1ll; sort(q+1, q+qn+1); ll res = q[1] - 1; for (int i=2; i<=qn; ++i) { if(q[i] != q[i-1]) res = res * (q[i] - 1); else res = res * q[i]; } return res; } inline void main(ll n) { qn = 0; PR(n); cout << getphi(n) << endl; } } int main() { srand(19260817); ll n; cin >> n; if(n == 1) { puts("1"); return 0; } PollardRho::main(n); return 0; }
bzoj4802 欧拉函数(附Millar-Rabin和Pollard-Rho讲解)
标签:题解 unsigned lin href srand 讲解 检验 表示 定义
原文地址:http://www.cnblogs.com/galaxies/p/bzoj4802.html