数论
一.扩展欧几里德（exgcd）
设 a>b。
1，显然当 b=0，gcd（a，b）=a。此时 x=1，y=0；
2，ab!=0 时
设 ax1+by1=gcd(a,b);
bx2+(a mod b)y2=gcd(b,a mod b);
根据欧几里德有 gcd(a,b)=gcd(b,a mod b);
则:ax1+by1=bx2+(a mod b)y2;
即:ax1+by1=bx2+(a-(a/b)*b)y2=ay2+bx2-(a/b)*by2=ay2+b(x2-(a/b)*y2);
根据恒等定理得：x1=y2; y1=x2-(a/b)*y2;
这样我们就得到了求解 x1,y1 的方法：x1，y1 的值基于 x2，y2.
上面的思想是以递归定义的，因为 gcd 不断的递归求解一定会有个时候 b=0，所以递归可以结束。
对于一般的不定式 a*x+b*y==c；有整数解得充分必要条件是（c % gcd（a，b）==0）。
若已经求出 a*x+b*y==gcd(a,b)的解。则原不定式的解是x1=x*(c/gcd(a,b)),y1=y*(c/gcd(a,b))。
那么原不定式的所有解xi=x1+b/gcd(a,b)*t; yi=y1-a/gcd(a,b)*t;其中t为任意常整数。
当gcd(a,b)==1时，形象化理解就是当你的x要+b时，你的y要-a，所以，把x变成非负整数的方法是：x=(x%b+b)%b。
void exgcd(int a,int b,int &x,int &y) { if(b==0){x=1;y=0;} else { exgcd(b,a%b); int t=x;x=y;y=t-a/b*y; } }
二.乘法逆元（模数为p）
求a^(-1)（mod p意义下）。
以下仅讨论gcd(a,p)==1的情况，否则不存在逆元。
在 mod p的意义下我们把x的乘法逆元写作 x^(-1)。
乘法逆元有如下的性质：x*x^(-1)≡1(mod p)。
应用：x/y≡x*y^(-1)(mod p)
1.费马小定理。
欧拉定理：当gcd(a,p)==1时，a^(φ(p))≡1(mod p)
a^φ(p)≡1(mod p)
a*a^(φ(p)-1)≡1(mod p)
∴a^(φ(p)-1)==a^(-1)。
用快速幂即可出解。
通常用在p为质数的情况：当p为质数时，φ(p)=p-1。
a^(p-2)==a^(-1)
2.exgcd。
求ax≡1(mod p)的x。
原方程可化为：ax+py=1。
由于gcd(a,p)==1，所以用exgcd求出x即可，最后把x弄成非负整数。
通常用此计算逆元。
3.线性求逆元（inv[i]表示i的逆元）。
要求p是质数。我们要求inv[1~p-1]。
设当前已知inv[1~i-1]。记a=p/i，b=p%i。则p=ai+b。
ai+b≡0(mod p)
方程两边同乘i^(-1)*b^(-1)
(ai+b)*i^(-1)*b^(-1)≡0(mod p)
a*b^(-1)+i^(-1)≡0(mod p)
i^(-1)≡-a*b^(-1)(mod p)
∴inv[i]=(-p/i*inv[p%i])%p=((p-p/i)*inv[p%i])%p。
初始化：inv[1]=1,inv[0]=1。
#include<cstdio> long long inv[10000000]; int main() { int n,p;scanf("%d%d",&n,&p);inv[0]=inv[1]=1; for(int i=2;i<=n;i++)inv[i]=((long long)(p-p/i)*inv[p%i])%(long long)p; for(int i=1;i<=n;i++)printf("%d\n",inv[i]); return 0; }
通过这种方法，我们还可以用递归的方式用O(log p)的时间求某个数的逆元。
即我们要查inv[i]，就dfs调用inv[p%i]。
三.BSGS
以下仅讨论p为质数且gcd(a,p)==1的情况
我们要求a^x≡b(mod p)（给出a,b）的最小非负整数的x
我们先考虑暴力枚举x。那么我们枚举的x应该何时停止？
由费马小定理得：a^(p−1)≡1(mod p)
a^x/a^(p−1)^m≡a^x(mod p)（m为非负整数）
a^(x-m*(p-1))≡a^x(mod p)
∵x-m*(p-1)==x%(p-1)
∴a^(x%(p-1))≡a^x(mod p)
所以，我们的x只要从0枚举到p即可。
设x=ki-j。
则a^(ki-j)≡b(mod p)
a^ki≡b*a^j(mod p)
同样的，我们的ki-j只要枚举到p即可。
那么我们令k=√p（上取整），把i从1枚举到k，j从0枚举到k，这样我们的ki-j就只会枚举到p了。
可是这样还是O(k^2)=O(p)的啊。
以下的b*a^j与a^ki均膜p。
我们先枚举j，把所有枚举出来的b*a^j压入hash表，如果有两个j，b*a^j的值是一样的，那么我们取大的一个j，因为这样ki-j最小。
然后再枚举i，对于每个a^ki，我们查表，如果查到一个j满足a^ki==b*a^j，那么我们停止搜索，此时的ki-j即为答案。
时间复杂度（可以使用hash表来代替map，降低时间复杂度，map好写）为O(√plogp)
typedef long long ll; map<ll,ll>mp; ll mod; ll ksm(ll a,ll p) { ll ans=1; while(p) { if(p&1)ans=(ans*a)%mod; a=(a*a)%mod;p>>=1; } return ans; } ll BSGS(ll a,ll b) { if(a%mod==0)return -1; mp.clear(); ll m=(ll)sqrt((double)mod)+1;ll am=ksm(a,m); ll now=1; for(int j=0;j<=m;j++){mp[(b*now)%mod]=j;now=(now*a)%mod;}now=1; for(int i=1;i<=m;i++){now=(now*am)%mod;if(mp[now])return i*m-mp[now];} return -1; }
四.中国剩余定理（CRT）
我们要求一个数x，满足以下这些同余方程（要求p1,p2,...,pn互质）。
x≡a1(mod p1)
x≡a2(mod p2)
...
x≡an(mod pn)
我们设x1是第一个同余方程的一个解，x2是第二个同余方程的一个解...xn是第n个同余方程的一个解
当x2~xn都是p1的倍数，x1,x3~xn都是x2的倍数......x1~xn-1都是xn的倍数时
我们的答案x=x1+x2+...+xn，此时的x可以符合每一个同余方程。
我们不妨假设此时的a1<p1,a2<p2...an<pn（否则可以mod p）
x1 mod p1==a1且是p2*p3*...*pn的公倍数。
x2 mod p2==a2且是p1*p3*...*pn的公倍数。
...
xn mod pn==an且是p1*p2*...*pn-1的公倍数。
我们知道：如果a mod b==1，那么ac mod b==c（c为正整数且c<b）
x1 mod p1==1且是p2*p3*...*pn的公倍数。x1*=a1。
x2 mod p2==1且是p1*p3*...*pn的公倍数。x2*=a2。
...
xn mod pn==1且是p1*p2*...*pn-1的公倍数。xn*=an。
记M=p1*p2*...*pn。
显然，当我们的xi是M/pi的倍数时，我们满足了第2个条件。
我们要找到一个ki，使得M/pi *ki≡1(mod pi)，这样就能满足第1个条件。
显然，ki是M/pi的逆元（mod pi意义下）！
于是，xi=M/pi * (M/pi)^(-1) * ai。x=x1+x2+...+xn。
记Mi=M/pi。则x=Σai*Mi*Mi^(-1)。
这只是一个解，其余解可以表示为Σai*Mi*Mi^(-1) +kM（k为整数）（因为我们不管是加上还是减去M，所得的x仍会符合以上的同余方程）。
所以最小非负整数解x'=Σai*Mi*Mi^(-1) %M
扩展中国剩余定理（ExCRT）
考虑合并
x≡a1(mod p1)
x≡a2(mod p2)
x=a1+k1p1=a2+k2p2
k1p1-k2p2=a2-a1
考虑Xp1+Yp2=gcd(p1,p2)
k1=X(a2-a1)/gcd(p1,p2)，k2=-Y(a2-a1)/gcd(p1,p2)
因为x=a1+k1p1 >=0 ，所以X为满足Xp1+Yp2=gcd(p1,p2)最小的>=0的数
这样a1+k1p1=a2+k2p2
所以，当x=a1+k1p1=a2+k2p2时可以满足上面两个方程，且x=a1+k1p1 + K*lcm(p1,p2)的x也都可以满足上面两个方程（K为任意自然数）
所以x≡a1+k1p1(mod lcm(p1,p2))的x，就能满足上面两个方程，且不满足该条件的x，不能满足上面两个方程
