基本数论
说一下一些基本数论算法的实现和简洁证明。
持续更新中。
\(exgcd\)
目的：求解方程 \(ax+by=(a,b)\) 的一组特解。
做法：
由于他叫exgcd，我们就按求 \(\gcd\) 方法试一下。 \[ \begin{aligned} 不妨设 (a,b)=1,a>b\\ a\cdot x+b\cdot y=&1\\ b\cdot y+(a\bmod b+\lfloor\frac{a}{b}\rfloor b)\cdot x=&1\\ 提取b\\ b\cdot(y+\lfloor \frac{a}{b}\rfloor x)+(a\bmod b)x=&1\\ 问题化归，求解\ bx_1+(a\bmod b)y_1=&1\\ \end{aligned} \] \(\mathcal{Code}\)1
2
3
4
5inline void exgcd(int a,int b,int &x,int &y)
{
if(b==0){x=1,y=0;return ;}
exgcd(b,a%b,y,x); y-=(a/b)*x;
}为啥能求逆元？ \[ \begin{aligned} ax+by&=1\\ ax+by&=1\pmod{b}\\ ax&=1\pmod{b} \end{aligned} \] 显然 \(x\) 是 \(a\) 在 \(\bmod b\) 意义下的逆元
\(excrt\)
目的：合并若干线性同余方程组。 \[ \begin{cases} x\equiv a_1\pmod{p_1}\\ x\equiv a_2\pmod{p_2}\\ \ \ \ \ \vdots\\ x\equiv a_n\pmod{p_n} \end{cases} \] 问题本质等价，合并两个方程。
\(x\equiv a_1\pmod{p_1},x\equiv a_2\pmod{p_2}\)
若存在解则解一定能表达成 \(x=x_0+k\cdot lcm(p_1,p_2)\),\(x_0\) 是一个特解。
显然 \(p_1\mid lcm(p_1,p_2),p_2\mid lcm(p_1,p_2)\)。
现在目的找到这样的一个特解 \(k\)。 \[ k=k_1\cdot p_1+a_1,k=k_2\cdot p_2+a_2\Rightarrow k_1\cdot p_1+a_1=k_2\cdot p_2+b_2\\ \Rightarrow k_1\cdot p_1-k_2\cdot p_2=b_2-b_1 \] 可以用 \(exgcd\) 求一组 \(k_1,k_2\) 特解。
·
\(Miller\ rabin\)
首先费马小定理判断一波,\(p\in \text{prime} \Leftrightarrow \forall (a,p)=1,a^{p-1}\equiv1\pmod{p}\)。
如果 \(a^{2}\equiv 1\pmod{p}\),如果 \(p\) 是质数，那么必定有 \(p\mid(a-1)\),或者 \(p\mid(a+1)\)。
否则 \(p\) 是合数。
\(a^{p-1}\equiv 1\pmod{p}\)。
此时如果 \(p\)是质数，\(a^{\frac{p-1}{2}}\equiv 1\ or\ -1\pmod{p}\)
如果是 \(-1\) 的时候，继续不下去，只得暂时假装他是素数。
否则继续检查 \(a^{\frac{p-1}{2}}\equiv 1\pmod{p}\)
本质上我们就是通过看 \(a^b\equiv 1\pmod{p}\) 这种形式来看 \(p\) 是否为质数。
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25inline bool mr(int x,int p)
{
int b=x-1;
while(!(b&1)) b>>=1;
int k=qpow(p,b,x);
while(b<=x-1)
{
int muls=k*k%x;
if(muls==1&&(k!=1&&k!=x-1)) return 0;
k=muls;b<<=1;
}
if(k!=1) return 0;
return 1;
}
int p[]={2,3,5,7,11,13,15,17,19,23};
inline bool judge(int x)
{
if(!(x&1)) return 0;
for(int i=0;i<=8;i++)
{
if(x==p[i]) return 1;
if(mr(x,p[i])==0) return 0;
}
return 1;
}当 \(p\) 的值域是
long long范围时，乘法可能会需要快速乘。
找阶&原根。
这个东西说实话只要真正了解定义就能如同行云流水般解决了。
阶
最小的使得 \(a^x\equiv 1\pmod{p}\) 的正整数 \(x\),被称为 \(a\bmod p\) 的阶。
由于 \(a^{\phi(p)}\equiv 1\pmod{p}\),所以 \(ord_p(a)\mid\phi(p)\)。
证明很显然，如果不整除，可以用类似辗转相除的方法推出，\(\exists \ x<ord_p(a),a^x\equiv 1\pmod{p}\)。
求法可以把 \(\phi(p)\) 的质因子一个一个消去，来看是否有 \(a^{\frac{\phi(p)}{k}}\equiv 1\pmod{p}\)。
复杂度 \(O(pollard-\rho)+\log(n)\)。
其实，如果我们知道阶的大小，也能很好知道阶为这个数的数有多少个。
若阶为 \(d\)。
那么 \(ord_{n}(x)=d\) 的 \(x\) 一定有 \(\phi(d)\) 个。
考虑我们随便找一个 \(ord_n(x)=d\),那么 \(x,x^2,cdots,x^{d-1},x^{d}\) 这样的一个集合。
显然我们应该在指数上与 \(d\) 互质的数。
所以 \(ord_x(x^k)\) 当且仅当 \((k,d)=1\)，此时也就是 \([1,d]\) 与 \(d\) 互质的数，也就是 \(\phi(d)\) 个。
原根
满足 \(ord_p(a)=\phi(p)\) 的 \(a\) 被叫做 \(p\) 的一个原根。
原根性质：若 \(g\) 为 \(p\) 其中一个原根，那么 \(1,g,g^2\cdots g^{\phi(p)-1}\) 互不相同。
证明可以用反证法。
由于 如果 \((k,p)=1\),有 \(0,k,2k,\cdots,(p-1)k\bmod p\) 一定互不相同。
也是反证，如果 \(xk\equiv yk\pmod{p}\),那么 \((x-y)k\equiv 0\pmod{p}\),与 \((k,p)=1\) 矛盾。
所有 \(p\) 有 \(\phi(\phi(p))\) 个原根。
原根有很多，随机找可以在 \(O(\frac{p}{\phi(\phi(p))})\) 的期望次找到。
设 \(\phi(p)=\prod p_i^{\alpha_i}\)。
那么判断一个数 \(x\) 是不是原根就看 是否 \(\forall \ k,x^\frac{\phi(p)}{p_k}\not \equiv 1\pmod{p}\)。
复杂度 \(O(pollard-\rho)+\frac{p}{\phi(\phi(p))}\log(n)\)。
\(pollard-\rho\)
在 \(\mathcal(O(\sqrt p)\)，\(p\) 为 \(n\) 最小质因子的复杂度，找到 \(n\) 的一个因数。
考虑对于最小质因子 \(p\)，我们随机枚举，枚举到 \(p\) 的倍数的概率为 \(\frac{1}{p}\)。
根据生日悖论，我们随机枚举，发生碰撞的概率则为 \(\frac{1}{\sqrt p}\)。
似乎问题马上就解决了。
但是考虑我们并不知道 \(p\) 的具体值，判断碰撞也只能使用 \(\gcd(\mid a-b\mid,n)\) 这么来判断。
那该怎么办？
\(pollard\) 发现了一种伪随机数列 \(a_{i+1}=a_i^2+c\bmod m\)。
可以证明该数列的形状为一条链加上一个环。
如果我们认为该数列是纯随机的那，环和链的长度均为期望 \(\sqrt m\) 级别的，根据生日悖论。
现在我们的任务就是给定 \(c,a_0\) 的条件下找到这个环。
但是，注意这样一点，我们并不能得到 \(a\) 的值，因为我们甚至都得不到 \(p\) 的值。
我们接下来做的都是用 \(\bmod n\) 的数列，反应 \(\bmod p\) 的数列。
找环的过程可以使用 \(a_i,a_{2i}\) 看什么呢时候冲突。
设链长为 \(l\)，环长为 \(m\)。
显然当 \(i=km>l\) 时，\(a_i=a_{l+(i-l)}=a_{l+(i-l+km)}=a_{i+km}=a_{2i}\)
也就是说，碰撞一定会产生。
在考虑 \(\bmod{p}\) 的数列，环和链长均为期望 \(\sqrt p\) 级别的，并且一定会发生碰撞，那么找到 \(p\) 的次数，一定是期望 \(\sqrt p\) 次。
那么一个复杂度大概为 \(\mathcal{O(\log_2 n\log n\sqrt p)}\) 的 \(pollard-\rho\) 算法出现了。
考虑在 \(a_i,a_{2i}\) 找环过程，我们是对每一组 \(\mid a_i-a_{2i}\mid\) 都和 \(n\) 求 \(\gcd\) 了。
这让我们有将所有 \(\mid a_i-a_{2i}\mid\) 乘起来一起做 \(\gcd\)。
考虑 \(\gcd(a,b)=\gcd(a,b\bmod a)\)，我们可以让 \(\mid a_{i}-a_{2i}\mid\) 按 \(\bmod n\) 乘积。
那就可以分块，每 \(128\) 分成一块，看 \(\gcd\)，这样我们可以认为一次找因数的复杂度被均摊成 \(\mathcal{O(\sqrt p)}\)。
如果 \(\sqrt p\) 比较小时，\(128\) 会太大，所以可以考虑倍增。
复杂度可以认为是 \(\mathcal{O(\log_2 n\sqrt p)}\)。由于每次找到因数后 \(n\) 的最小因数会改变，所以复杂度可能会更优，但是具体怎么证明我不是很会，在具体应用时你可以认为前面那个 \(\log_2\) 因子会十分的小，甚至忽略。
\(\mathcal{O(\log_2 n \sqrt p)}\) 是期望复杂度上界。
实现，边界。
- 如果发现 \(a_i=a_{2i}\) 立马停止，重新换 \(a_0,c\) 重新 \(pollard-\rho\)。
- 如果 \(n=4\) 那么 \(a_i=a_{2i}\) 据说不可避免，可以特判掉。
- 如果 \(n\) 为质数，立马停止。
大概代码：
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
41inline ll nxt(ll a,ll c,ll mod){
ll tmp=mul(a,a,mod)+c;
return tmp>=mod?tmp-mod:tmp;
}
inline ll gcd(ll a,ll b){
if(!b) return a;
return gcd(b,a%b);
}
inline ll pollard(ll n){
if(n==4) return 2;
ll a1=1ll*rand()*rand()%n,c=1ll*rand()*rand()%n;
ll a2=a1;
a1=nxt(a1,c,n),a2=nxt(nxt(a2,c,n),c,n);
int lim=1;
while(a1!=a2){
ll v=1;
for(int i=1;i<=lim&&a1!=a2;i++,a1=nxt(a1,c,n),a2=nxt(nxt(a2,c,n),c,n)){
ll backv=v;
v=mul(v,abs(a1-a2),n);
if(!v) return gcd(abs(a1-a2),n);
}
ll g=gcd(v,n);
if(g>1) return g;
if(lim<128) lim<<=1;
}
return pollard(n);
}
vector<ll> f;
inline void factor(ll n){
if(prime(n)){
f.push_back(n);
return ;
}
ll fac=pollard(n);
if(prime(fac)){
f.push_back(fac);
return factor(n/fac);
}
factor(fac);
factor(n/fac);
}