gcd
朴素欧几里得（辗转相除法）
- 这一节我们默认：\(a,b\in\mathbb{Z}\)
对于求解 \(\gcd(a,b)\) 需要用朴素欧几里得定理。
\[\gcd(a,b)=\gcd(b,a\bmod b) \]gcd
-
\(\gcd\) 为 \(\text{greatest common divisor}\) 的缩写，即最大公约数。（或称最大公因数）
-
对于 \(\gcd\) 函数，有以下几个性质：
- \(\gcd(a,b)=\gcd(b,a)\)
- \(\gcd(a,-b)=\gcd(a,b)\)
- \(\gcd(a,a)=|a|\)
- \(\gcd(a,0)=|a|\)
- \(\gcd(a,1)=1\)
- 若有常数 \(k\in \mathbb{Z}\)，则：
\(\gcd(ka,kb)=k\times\gcd(a,b)\)
\(\gcd(a+kb,b)=\gcd(a,b)\)
对于定理的证明：
令 \(c=\gcd(a,b)\)，则存在 \(a=mc,b=nc\)。
令 \(r=a\bmod b\)，则存在\(k\)，使\(r=a-kb=mc-knc=c\times(m-kn)\)。
\(\therefore \gcd(b,a\bmod b)=\gcd(b,r)=\gcd(nc,c\times(m-kn))=c\times gcd(n,m-kn)\)。
\(\therefore c\) 为 \(a\) 和 \(b\) 的公约数。
令 \(d=\gcd(n,m-kn)\)，则存在 \(x,y\) 使 \(n=xd,m-kn=yd\)
\(\therefore m=yd+kn=yd+kxd=d\times (y+kx)\)
\(\therefore a=mc=d\times(y+kx)c=dc\times(y+kx)\)，\(b=nc=xdc=dc\times x\)
\(\therefore \gcd(a,b)=\gcd(\red{dc}\times(y+kx),\red{dc}\times x)=\red{dc}\)
\(\because \gcd(a,b)=c\)
联立得：\(dc=c\)
\(\therefore \gcd(b,a\bmod b)=c\)
联立得：\(\gcd(a,b)=c=\gcd(b,a\bmod b)\)
- 代码实现：
有两种写法，个人推荐第一种。时间复杂度为 \(O(\log n)\)。
inline int gcd(int a,int b){
int temp;
while(b){//边界条件为b=0;
temp=b;
b=a%b;
a=temp;
}
return a;
}
inline int gcd(int a,int b){return b?gcd(b,a%b):a;}
扩展欧几里得(exgcd)
- 对于三个整数 \(a,b,c\)，求解 \(ax+by=c\) 的整数解。
首先要判断是否有解，有解的充分条件是 \(\gcd(a,b)\mid c\)（裴蜀定理）。
然后判定是否有解后，我们需要在这个基础上求一组解 \((x,y)\)
若 \(a<0\) 或 \(b<0\) 用上述的性质化为正数即可。
先将式子两边同除 \(\gcd(a,b)\)，所以式子转化为了：
\(ax+by=1(a\bot b)\)
用前面的朴素欧几里得推式子：
\[\begin{aligned}ax+by&=\gcd(a,b)\\&=\gcd(b,a\bmod b)\\&=bx+(a\bmod b)y\\&=bx+(a-\left\lfloor\dfrac{a}{b}\right\rfloor b)y\\&=ay+b(x-\left\lfloor\dfrac{a}{b}\right\rfloor y)\end{aligned} \]此时 \(x\Leftarrow y,y\Leftarrow ay+b(x-\left\lfloor\dfrac{a}{b}\right\rfloor y)\)。
所以可以递归求解。
边界条件和朴素欧几里得一样，当 \(b=0\) 时，有 \(a=1\) 且 \(ax+by=1\)，所以 \(x=1,y=0\)。
这样做完的话我们用 \(O(\log n)\) 的时间就会得到一组 \((x,y)\) 的特殊解。
不定方程通解公式
有方程 \(ax+by=c\) 且已知 \(x=x_0,y=y_0\)。
\(\therefore a(x-x_0)+b(y-y_0)=0\)
\(\therefore a(x-x_0)=b(y_0-y)\)
设 \(a\bot b\)，则 \(x-x_0\) 必含因子 \(b\)，\(y_0-y\) 必含因子 \(a\)
\(\therefore x-x_0=kb,y_0-y=ka(k\in \mathbb Z)\)
\(\therefore \begin{cases}x=x_0+kb\\y=y_0-ka\end{cases}(k\in \mathbb Z)\)
解系拓展
有时题目要求的不是 \((x,y)\) 这一解，而是某种特殊解，这时需将其拓展成解系。
根据不定方程通解公式可得：
\[\begin{cases}x=x_0+kb\\y=y_0-ka\end{cases}(k\in \mathbb Z) \]然后根据题面找出特殊解即可。
- 注意：
-
最后都要乘\(c\)。
-
记得令 \(x=(x\bmod b+b)\bmod b\) 就行了。防止爆 \(\text {int}\) 或其他玄学错误
- 代码实现
inline void Exgcd(int a,int b,int &x,int &y){
if(b)Exgcd(b,a%b,y,x),y-=a/b*x;//同样的递归条件
else x=1,y=0;
}
调用的话，\(x=y',y=x'\)，只修改 \(y\) 的值。
更相减损术
- 更相减损术是出自《九章算术》的一种求最大公约数的算法，它原本是为约分而设计的，但它适用于任何需要求最大公约数的场合。
原文：
可半者半之，不可半者，副置分母、子之数，以少减多，更相减损，求其等也。以等数约之。
译文(来自百度)：（如果需要对分数进行约分，那么）可以折半的话，就折半（也就是用2来约分）。如果不可以折半的话，那么就比较分母和分子的大小，用大数减去小数，互相减来减去，一直到减数与差相等为止，用这个相等的数字来约分。
- 代码实现：
inline int gcd(int a,int b){
while(a!=b)
a>b?a-=b:b-=a;
return a;
}
时间复杂度为 \(O(\log n)\)，但是如果遇到一个数很大，另一个数比较小的情况，可能要进行很多次减法才能达到一次除法的效果，从而使得算法的时间复杂度退化为 \(O(n)\)，因此几乎不会用到此算法。
乘法逆元
- 众所周知，取模有可加、减、乘性：
但是，取模是没有可除性的。模意义下的除法如何运算呢？这里引进乘法逆元。
- 定义：若\(a\times x\equiv 1\pmod p\)，则称 \(x\) 为 \(a\) 在模 \(p\) 意义下的乘法逆元，记为 \(a^{-1}\pmod p\)。
求逆元
- 我们可以从 \(a\times x\equiv 1\pmod p\) 这个式子入手。
不难发现，它与 \(\exists y\in\mathbb{Z},ax+by=1\) 等价。
- 所以它就变成一个关于 \((a,p)\) 的不定方程，可用上一节讲到的\(\text{exgcd}\)解决。
并且可以由裴蜀定理得到：在 \(a\) 在模 \(p\) 意义下的乘法逆元当且仅当 \(\gcd(a,p)=1\)，\(a\bot b\)。
- 代码实现：
inline int solve(int a,int p){
int x,y;
exgcd(a,p,x,y);//求不定方程ax+py=1的一组解;
return x;//求出的x即为a在模p意义下的乘法逆元;
}
- 当然还可以用费马小定理，即欧拉定理的一种特殊情况求解，这里我们后面再讲。
中国剩余定理(CRT)
-
缩写为 \(\text{CRT}\)，全称为 \(\text{Chinese remainder theorem}\)
-
一般地，中国剩余定理可以求解以下这类同余方程组：
- 其中，\(m_1,m_2,\cdots m_n\) 两两互质，则对于任意地整数 \(a_1,a_2,\cdots a_n\) 上面这个方程组有解。
算法流程
- 举个例子：\(n=3,a=\{2,3,2\},m=\{3,5,7\}\)。
则方程组为：
\[\begin{cases}x\equiv 2\pmod{3}\\x\equiv 3\pmod{5}\\x\equiv 2\pmod{7}\end{cases} \]- 设 \(M=m_1\times m_2\times \cdots \times m=\prod\limits_{i=1}^na_i\)。
- 设 \(M_i=\dfrac{M}{m_i}\)。
- 设 \(t_i=m_i^{-1}\) 为 \(M_i\) 模 \(m_i\) 的逆元。（因为 $ M_i\bot m_i$，所以逆元必存在）
- \(x\) 的通解为 \(a_1 t_1 M_1+a_2 t_2 M_2+\cdots+a_n t_n M_n+kM=\sum\limits^n_{i=1}a_i t_i M_i+kM(k\in \mathbb{Z})\)
- 故通解为：\(x=23+k\times 105(k\in \mathbb{Z})\)
例题：P1495 【模板】中国剩余定理（CRT）/ 曹冲养猪
-
读题可知，题目中的 \(a_i,b_i\) 分别对应上述的 \(m_i,a_i\)。直接使用中国剩余定理解决即可。
-
code:
#include<bits/stdc++.h>
using namespace std;
#define LL long long
inline LL read(){
LL s=0;int f=1;char ch=getchar();
while(ch<'0' or ch>'9'){if(ch=='-')f=-1;ch=getchar();}
while(ch>='0' and ch<='9'){s=(s<<1)+(s<<3)+(ch^48);ch=getchar();}
return f*s;
}
inline void write(LL num){
if(num<0)putchar('-'),num=-num;
if(num>9)write(num/10);
putchar(num%10+48);
}
const int N=15;
LL n,x,y,MM,m[N],a[N],t[N],M[N];
inline LL exgcd(LL a,LL b,LL &x,LL &y){
if(b==0){
x=1,y=0;
return a;
}
LL d=exgcd(b,a%b,y,x);
y-=a/b*x;
return d;
}
int main(){
n=read();
MM=1,x=0;
for(int i=1;i<=n;++i){
m[i]=read(),a[i]=read();
MM*=m[i];//M=m1*m2*...mn
}
for(int i=1;i<=n;++i){
M[i]=MM/m[i];//Mi=M/mi
exgcd(M[i],m[i],t[i],y);//ti=Mi^-1(mod mi)
t[i]=(m[i]+t[i]%m[i])%m[i];//保证0<=ti<mi
for(int j=1;j<=a[i];++j)//防止爆long long，用多次加法模拟乘法
x=(x+M[i]*t[i])%MM;
}
write(x);
return 0;
}
扩展中国剩余定理（exCRT）
-
因为上述的 \(m_1,m_2,\cdots,m_n\) 两两互质，会使最后构造出来的 \(M_i\) 特别的美好，只用把和式 \(\sum\limits_{1}^n M_i t_i\) 输出就行了，因为其 \(\bmod\ m_1\) 余 \(a_1\)，\(\bmod\ m_2\) 余 \(a_2\)，$\cdots\ \(，\)\bmod\ m_n$ 余 \(a_n\)。但求 \(M_i\) 需要应用到逆元，因此 \(m_1,m_2,\cdots,m_n\) 必须两两互质。
-
那么当它们两两不互质时，又给如何解决呢？
先从最简单的情况看：假设方程组只有一个方程，即：
\[x\equiv a\pmod m \]-
这个式子等价于：\(k\times a+m=x(x\in\mathbb{Z})\)，这就是这个同于方程的解系。
-
再到两个方程，此时很明显无法直接求解。我们是否可以构造一个新的方程，使其解系于原方程组的解系等价呢？即：
答案是肯定的。
算法流程
- 左边那一个方程组等价于：
移项，即可得到：
\[k_1 a_1-k_2 a_2=m_2-m_1 \]因为 \(a_1,a_2,m_2-m_1\) 都是已知的，所以这是的典型的不定方程。
- 因此可以用裴蜀定理判断是否有解，用 \(\text{exgcd}\) 求解系 \((k_1,k_2)\)。
- 接下来考虑如何求出一组 \((a_1, a_2)\)。
令 \(d=\gcd(a_1,a_2),p_1=\dfrac{a_1}{d},p_2=\dfrac{a_2}{d}\)。
显然 \(p_1\bot p_2\)，所以 \(a_1=p_1 d,a_2=p_2 d\)。
- 因此原式化为：
- 右边那一串东西是整数，因为当且仅当 \(d|(m_2-m_1)\) 才会有解。左边满足了 \(p_1\bot p_2\)，因此求出整个解系是很容易的。现在假设我们拿\(\text{exgcd}\)求出了下面这个方程的解：\((\lambda_1, \lambda_2)\)
- 显然，这是个非常标准的 \(\text{exgcd}\)，因此可用这个把 \((\lambda_1,\lambda_2)\) 求出来后，就可以将 \((k_1,k_2)\) 整出来：
-
至此，我们成功地构造出了一个特殊解 \(x\)，满足\(\begin{cases}x\equiv a_1\pmod{m_1}\\x\equiv a_2\pmod{m_2}\end{cases}\)。
-
如何求出整个解系呢？
【定理】 若有特解 \(x'\)，那么\(\begin{cases}x\equiv a_1\pmod{m_1}\\x\equiv a_2\pmod{m_2}\end{cases}\)的通解为：\(x'+k\times \text{lcm}(m_1,m_2)(k\in\mathbb{Z})\)
- 从线性代数的角度讲，这个通解的构造方式是十分平凡的。
对 \(\text{lcm}(m_1,m_2)\) 取模的结果，将整个整数集划分成了 \(\text{lcm}(m_1,m_2)\) 个等价类，哪个等价类里面有特解，那整个等价类肯定全都是解。
接下来唯一需要说明的事情就是：为什么任意一个完全剩余系里面，只会有一个解？这个问题等价于：
- 为什么在 \(0,1,2,\cdots ,\text{lcm}(m_1,m_2)\) 中，只有一个解？
设上述集合里面有 \(0\leq x,y\leq\text{lcm}(m_1,m_2)\) 满足：
\[\begin{cases}x\equiv a_1\pmod{m_1}\\x\equiv a_2\pmod{m_2}\end{cases},\begin{cases}y\equiv a_1\pmod{m_1}\\y\equiv a_2\pmod{m_2}\end{cases} \]- 不妨设 \(x\leq y\)，那么：
- 因此一个完全剩余系中，有且仅有一个解。
以上就是整个 \(\text{exCRT}\) 算法的全部数学基础。
- code：
#include<bits/stdc++.h>
using namespace std;
typedef __int128 ll;
ll x,y,d;
int n;
inline void Exgcd(ll &x,ll &y,ll a,ll b){
if(b)Exgcd(y,x,b,a%b),y-=a/b*x;
else d=a,x=1,y=0;
}
inline ll gcd(ll a,ll b){
ll temp;
while(b){
temp=b;
b=a%b;
a=temp;
}
return a;
}
inline ll lcm(ll a,ll b){
return a/gcd(a,b)*b;
}
ll a,b,A,B;
inline void merge(){
Exgcd(x,y,a,A);
ll c=B-b;
if(c%d)
puts("-1"),exit(0);
x=x*c/d%(A/d);
if(x<0)
x+=A/d;
ll mod=lcm(a,A);
b=(a*x+b)%mod;
if(b<0)
b+=mod;
a=mod;
}
int main(){
scanf("%d",&n);
for(long long i=1,_A,_B;i<=n;++i){
scanf("%lld %lld",&_A,&_B);
A=_A,B=_B;
if(i>1)
merge();
else
a=A,b=B;
}
printf("%lld\n",(long long)(b%a));
}
欧拉函数
- 定义：\(\varphi(n)\) 表示 \(1\sim n\) 中与 \(n\) 互质的数的个数。
它是一个积性函数，即当 \(m\bot n\) 时，\(\varphi(mn)=\varphi(m)\times\varphi(n)\)。
- 知道了欧拉函数时积性函数后，只需要考虑 \(\varphi(p^k)\) 的取值就能计算所有数的欧拉函数了，其中 \(p\) 为质数。
不难发现若一个数不和 \(p^k\) 互质，则这个数一定是 \(p\) 的倍数，那么 \(1\sim p^k\) 中有 \(\frac{p^k}{p}=p^{k-1}\) 个数和 \(p^k\) 不互质的。反过来和 \(p^k\) 互质的个数就是 \(\varphi(p^k)=p^k-p^{k-1}=p^k(1-\frac{1}{p})\)。
- 所以设 \(n=p_1^{k_1}p_2^{k_2}\cdots p_m^{k_m}\)，则：
- 可以 \(O(\sqrt{n})\) 计算 \(\varphi(n)\)
inline void solve(){
int phi=n;
for(int i=2;i*i<=n;++i)
if(n%i==0){ //如果i是n的一个质因子
phi=phi/i*(i-1); //乘上(1-1/p)
while(n%i==0) //把n中的i全部除掉
n/=i;
}
if(n!=1)
phi=phi/n*(n-1); //n有一个大于sqrt(n)的质因子的情况
}
- 也可以用线性筛求欧拉函数：
inline void solve(){
phi[1]=1;
for(int i=2;i<=n;++i){
if(!phi[i]) //phi(i)=0，说明i好没被筛去，是质数
p[++tot]=i,phi[i]=i-1; //phi(p)=p-1
for(int j=1;j<=tot and i*p[j]<=n;++j){
if(i%p[j]==0)[
phi[i*p[j]]=phi[i]*p[j]; //phi(p^k*x)=phi(p^(k-1)*x)*p
break;
}else
phi[i*p[j]]=phi[i]*(p[j]-1); //phi(p*x)=phi(x)*(p-1)
}
}
}
欧拉定理
- 设 \(a,n\in\mathbb{N^+}\)，若 \(a\bot n\)，则 \(a^{\varphi(n)}\equiv 1\pmod n\)。
当 \(n\) 是一个质数的时候，\(\varphi(n)=n-1\)。因此可以得到费马小定理：
- 若 \(p\) 为质数，且 \(p\nmid a\)，则\(a^{p-1}\equiv 1 \pmod p\)。
因此也可以用费马小定理求逆元：
- 若 \(p\) 是质数，就可以知道如果 \(a\) 不是 \(p\) 的倍数，就有 \(a\times a^{p-2}=a^{p-1}\equiv 1\pmod p\)，根据乘法逆元的定义，\(a\) 在模 \(p\) 意义下的逆元就是 \(a^{p-2}\)。否则逆元不存在。
inline int qPow(int a,int b,int p){//快速幂
int ans=1;
while(b>0){
if(b&1)ans=ans*a%p;
a=a*a%p;
b>>=1;
}
return ans;
}
inline int solve(int a,int p){
return qPow(a,p-2,p);
}
- 对于欧拉定理的证明，这里我会同步举例：
- 考虑 \(1\sim n\) 中与 \(n\) 互质的数组成的集合 \(X_n=\{x_1,x_2\cdots x_{\varphi(n)}\}\)
- 将 \(X_n\) 集合中的每一个数乘上 \(a\) 再\(\mod {n}\)，得到\(X_n'=\{ax_1\bmod{n},ax_2\bmod{n},\cdots\,ax_{\varphi(n)}\bmod{n}\}\)。
- 因为两个集合中的元素都是 \(1\sim n\) 中与 \(n\) 互质的数，并且原来不相同的两个数乘 \(a\) 后依然不可能相同，否则再乘上 \(a\) 后逆元就一样了。所以 \(X_n=X_n'\)。因此可以得出：
两边同乘\(\left(\prod_{x_i\in X_n}x_i\right)\)的逆元，于是有：
\[a^{\varphi(n)}\equiv 1\pmod n \]扩展欧拉定理
- 设 \(a,n\in\mathbb{N^+}\)，若 \(a\bot n\)，则 \(a^{2\varphi(n)}\equiv a^{\varphi(n)}\pmod n\)。
比如说，\(2\) 和 \(12\) 不互质，但：
\[2^{2\varphi(12)}=2^8=256\equiv4\pmod {12} \]\[2^{\varphi(12)}=2^4=16\equiv4\pmod {12} \]\[\therefore 2^{2\varphi(12)}\equiv2^{\varphi(12)}\pmod{12} \]- 两者依旧相同。
考虑扩展欧拉定理的证明。由唯一分解定理，设 \(n=\prod_{i=1}^m p_i^{k_i}\)。
若对于每个 \(i\)，我们都可以证明 \(a^{2\varphi(n)}\equiv a^{\varphi(n)}\pmod {p_i^{k_i}}\)，则 \(a^{2\varphi(n)}\) 和 \(a^{\varphi(n)}\) 都是下面这个同余方程的解：
\[\begin{cases}x\equiv a^{\varphi(n)}\pmod{ M_1^{k_1}}\\x\equiv a^{\varphi(n)}\pmod{ M_2^{k_2}}\\\cdots\\x\equiv a^{\varphi(n)}\pmod{ M_m^{k_m}}\end{cases} \]由中国剩余定理的唯一性，扩展欧拉定理也就成立了。下面来分情况讨论：
- 求证：\(a^{2\varphi(n)}\equiv a^{\varphi(n)}\pmod {p_i^{k_i}}\)
- 若 \(a\) 与 \(p_i^{k_i}\) 互质。
设 \(s=\dfrac{n}{p_i^{k_i}}\)，那么 \(s\) 和 \(p_i^{k_i}\) 显然也是互质的。由于 \(\varphi(x)\) 是积性函数，所以 \(\varphi(n)=\varphi(s)\varphi(p_i^{k_i})\)。
- 由欧拉定理：
两边同乘 \(a^{\varphi(n)}\)，即得到 \(a^{2\varphi(n)}\equiv a^{\varphi(n)}\pmod {p_i^{k_i}}\)。
- 若 \(a\) 与 \(p_i^{k_i}\) 不互质。
那么 \(a\) 一定是 \(p_i\) 的倍数，\(a^{\varphi(n)}\) 一定是 \(p_i^{\varphi(n)}\) 的倍数。同样设 \(s=\dfrac{n}{p_i^{k_i}}\)。
- 首先可以用数学归纳法证明 \(p_i^{k_i-1}\geq k_i\)。
当 \(k_i=1\) 时，\(p_i^{k_i-1}=p_i^0=1=k_i\).
当 \(k_i>1\) 时，\(p_i^{k_i-1}=p_i\times p_i^{k_i-2}\geq 2(k_i-1)\geq k_i\)。\(\left(p_i\geq 2,p_i^{k_i-2}\geq k_i-1\right)\)
- 于是有：
所以 \(p_i^{\varphi(n)}\) 是 \(p_i^{k_i}\) 的倍数，那么 \(a^{\varphi(n)}\) 也就是 \(p_i^{k_i}\) 的倍数。
- 因此得到：
综上所述，\(\forall i\in \mathbb{Z}\cap\left[1,m\right] ,a^{2\varphi(n)}\equiv a^{\varphi(n)}\pmod {p_i^{k_i}}\)，扩展欧拉定理成立。
- 将扩展欧拉定理的两边同乘 \(a^{\varphi(n)}\)，可以得到 \(a^{3\varphi(n)}\equiv a^{2\varphi(n)}\pmod n\)。
不断重复这个过程，可以得到：
\[a^{\varphi(n)}\equiv a^{2\varphi(n)}\equiv\cdots\equiv a^{k\varphi(n)}\equiv\cdots\pmod n \]- 所以对于任意的 \(a,k,n\)，就有：\(a^k\equiv a^{(k\text\ {mod}\ \varphi(n))+\varphi(n)}\pmod n\)。
BSGS
拔山盖世Baby Steps Giant Steps。
例题：P3846 [TJOI2007] 可爱的质数/【模板】BSGS
- \(\text{BSGS}\)适用于对于给定的两个互质的正整数 \(a,p\)，求一个非负整数 \(x\) 使得：
因为 \(a \bot p\)，所以根据欧拉定理：
\[a^{\varphi(p)}\equiv 1\pmod p \]又因为：
\[a^0 \equiv 1\pmod p \]所以 \(1\sim \varphi(p)\) 是一个循环节。
- 那么最多循环 \(\varphi(p)\) 次我们能找到答案，否则无解。暴力带么就很容易写出来了。
设 \(x=im-k\)，其中 \(0\leq k\leq m\)。
- 那么方程变为：
两边同乘 \(a^k\):
\[a^{im}\equiv a^kb\pmod p \]- 先计算右边 \(a^kb\mod p\) 的值，把它放到 \(\text{hash}\) 或者 \(\text{map}\) 里。再计算左边 \(a^{im}\mod p\) 的值，从小到大枚举所有可能的 \(i\) 值，再在 \(\text{hash}\) 或者 \(\text{map}\) 里查询。
这里的话倾向用 \(\text{hash}\)，因为 \(\text{map}\) 会多带一个 \(\log\)。
时间复杂度：\(O(\max\left(m,\dfrac{\varphi(p)}{m}\right))=O(\sqrt{p})\)
- code：
#define ll long long
inline ll qPow(ll a,ll b,ll p){
int ans=1;
while(b>0){
if(b&1)ans=ans*a%p;
a=a*a%p;
b>>=1;
}
return ans;
}
namespace work{
m=ceil(sqrt(p));
for(ll i=0,t=n;i<=m;i++,t=t*b%p)
vis[t]=i;
for(ll i=1,tt=qPow(b,m,p),t=tt;i<=m;i++,t=t*tt%p)
if(vis.count(t))
write(i*m-vis[t]),exit(0);
puts("no solution");
}
exBSGS
ex北上广深这玩意还有拓展，用于处理 \(\gcd(a,p)\neq 1\)的情况。
- 首先，上面那个式子可以等价写成：
- 设 \(\gcd(a,p)=d\)。由于裴蜀定理，这个式子有解的充要条件是：\(d\mid n\)。
那么上面的式子就是：
\[a^{x-1}\cdot\dfrac{a}{d}+k\cdot\dfrac{p}{d}=\dfrac{n}{d} \]此时 \(a^x\leftarrow a^{x-1},p\leftarrow\dfrac{p}{d},n\leftarrow\dfrac{n}{d}\)。那么一直递归，直到 \(a\bot p\)。
- 所以此时设一共递归了 \(\text{cnt}\) 次，这所有次递归的 \(d\) 的乘积是 \(D\)。原式就是：
那么我们就用 \(a'=a,n'=\dfrac{n}{D},p'=\dfrac{p}{D}\) 来跑一次 \(\text{BSGS}\) 就行了。
-
最后要再加上 \(\text{cnt}\)，并乘上 \(\dfrac{a^{cnt}}{D}\)的系数作为答案。
-
code：
inline int exBSGS(int a,int b,int p){
a%=p,b%=p;
if(b==1 or p==1)
return 0;
int d,ax=1,cnt=0,x,y;
while((d=exgcd(a,p,x,y))^1){
if(b%d)
return -1;
b/=d,p/=d,cnt++;
ax=1ll*ax*(a/d)%p;
if(ax==b)
return cnt;
}
exgcd(ax,p,x,y);
int inv=(x%p+p)%p;
b=1ll*b*inv%p;
_hash.clear();
int t=ceil(sqrt(p)),val=1;
for(int i=0;i<t;i++){
_hash[1ll*b*val%p]=i;
val=1ll*val*a%p;
}
a=val,val=1;
if(!a)
return !b?1+cnt:-1;
for(int i=0,j;i<=t;i++){
j=_hash.find(val)==_hash.end()?-1:_hash[val];
if(~j and i*t-j>=0)
return i*t-j+cnt;
val=1ll*val*a%p;
}
return -1;
}
Wilson
对于素数 \(p\) 有 \((p-1)!\equiv -1\pmod p\)。
- 证明：
当 \(p=2\) 时，命题显然成立。
以下设 \(p\geq 3\)，此时我们要证 \(\mathbf{Z}_p\) 中所有非零元素的积为 \(\overline{-1}\)。
我们知道 \(\mathbf{Z}_p\) 中所有非零元素 \(a\) 都有逆元 \(a^{-1}\)，于是 \(\mathbf{Z}_p\) 中彼此互逆的元素乘积为
\(\overline{1}\)。
但是要注意 \(a\) 和 \(a^{-1}\) 可能相等。若 \(a=a^{-1}\)，则 \(a^2\equiv 1\pmod p\)，
即 \(0\equiv a^2-1\equiv (a+1)(a-1)\pmod p\)，
从而 \(a\equiv 1\pmod p\) 或 \(a\equiv -1\pmod p\)。
这说明：\(\mathbf{Z}_p\setminus\{\overline{0},\overline{1},\overline{-1}\}\) 中所有元素的乘积为 \(\overline{1}\), 进而 \(\mathbf{Z}_p\) 中所有非零元素的积为 \(\overline{-1}\)。
对于整数 \(n\)，令 \((n!)_p\) 表示所有小于等于 \(n\) 但不能被 \(p\) 整除的正整数的乘积，即 \((n!)_p=n!/(\lfloor n/p\rfloor !p^{\lfloor n/p\rfloor})\)。
- Wilson 定理指出 \((p!)_p=(p-1)!\equiv -1\pmod p\) 且可被推广至模素数 \(p\) 的幂次。
应用
阶乘与素数
在某些情况下，有必要考虑以某个素数 \(p\) 为模的公式，包含分子和分母中的阶乘，就像在二项式系数公式中遇到的那样。
只有当阶乘同时出现在分数的分子和分母中时，这个问题才有意义。否则，后续项 \(p!\) 将减少为零。但在分数中，因子 \(p\) 可以抵消，结果将是非零。
因此，要计算 \(n! \bmod p\)，而不考虑阶乘中出现因数 \(p\)。写下素因子分解，去掉所有因子 \(p\)，并计算乘积模。
用 \((n!)_p\) 表示这个修改的因子。例如：
\[(7!)_p \equiv 1 \cdot 2 \cdot \underbrace{1}_{3} \cdot 4 \cdot 5 \underbrace{2}_{6} \cdot 7 \equiv 2 \bmod 3 \]- 这种修正的阶乘，可用于快速计算各种带取模和组合数的公式的值。
计算余数的算法
计算上述去掉因子 \(p\) 的阶乘模 \(p\) 的余数。
\[\begin{aligned}(n!)_p &= 1 \cdot 2 \cdot 3 \cdot \ldots \cdot (p-2) \cdot (p-1) \cdot \underbrace{1}_{p} \cdot (p+1) \cdot (p+2) \cdot \ldots \cdot (2p-1) \cdot \underbrace{2}_{2p} \\\\&\quad \cdot (2p+1) \cdot \ldots \cdot (p^2-1) \cdot \underbrace{1}_{p^2} \cdot (p^2 +1) \cdot \ldots \cdot n \pmod{p} \\\\&= 1 \cdot 2 \cdot 3 \cdot \ldots \cdot (p-2) \cdot (p-1) \cdot \underbrace{1}_{p} \cdot 1 \cdot 2 \cdot \ldots \cdot (p-1) \cdot \underbrace{2}_{2p} \\\\&\cdot 1 \cdot 2 \quad \cdot \ldots \cdot (p-1) \cdot \underbrace{1}_{p^2} \cdot 1 \cdot 2 \cdot \ldots \cdot (n \bmod p) \pmod{p} \end{aligned}\]- 可以清楚地看到，除了最后一个块外，阶乘被划分为几个长度相同的块。
- 块的主要部分 \((p-1)!\ \mathrm{mod}\ p\) 很容易计算，可以应用 Wilson 定理：
总共有 \(\lfloor \frac{n}{p} \rfloor\) 个块，因此需要将 \(\lfloor \frac{n}{p} \rfloor\) 写到 \(1\) 的指数上。可以注意到结果在 \(-1\) 和 \(1\) 之间切换，因此只需要查看指数的奇偶性，如果是奇数，则乘以 \(-1\)。除了使用乘法，也可以从结果中减去 \(p\)。
最后一个部分块的值可以在 \(O(p)\) 的时间复杂度单独计算。
这只留下每个块的最后一个元素。如果隐藏已处理的元素，可以看到以下模式：
\[(n!)_p = \underbrace{ \ldots \cdot 1 } \cdot \underbrace{ \ldots \cdot 2} \cdot \ldots \cdot \underbrace{ \ldots \cdot (p-1)} \cdot \underbrace{ \ldots \cdot 1 } \cdot \underbrace{ \ldots \cdot 1} \cdot \underbrace{ \ldots \cdot 2} \cdots \]这也是一个修正的阶乘，只是维数小得多。它是：
\[\left(\left\lfloor \frac{n}{p} \right\rfloor !\right)_p \]因此，在计算修改的阶乘 \((n!)_p\) 中，执行了 \(O(p)\) 个操作，剩下的是计算 \((\lfloor \frac{n}{p} \rfloor !)_p\)，于是有了递归，递归深度为 \(O(\log_p n)\)，因此算法的总时间复杂度为 \(O(p \log_p n)\)。
- 可以预先计算阶乘 \(0!, 1!, 2!, \dots, (p-1)!\pmod p\)，那么时间复杂度为 \(O(\log_p n)\)。
计算余数算法的实现
具体实现不需要递归，因为这是尾部递归的情况，因此可以使用迭代轻松实现。在下面的实现中预计算了阶乘 \(0!, 1!, \dots, (p-1)!\)。
因此时间复杂度为 \(O(p + \log_p n)\)。如果需要多次调用函数，则可以在函数外部进行预计算，于是计算 \((n!)_p\) 的时间复杂度为 \(O(\log_p n)\)。
- code
inline int factmod(int n, int p){
vector<int>f(p);
f[0]=1;
for(int i=1;i<p;i++)
f[i]=f[i-1]*i%p;
int res=1;
while(n>1){
if((n/p)%2)
res=p-res;
res=res*f[n%p]%p;
n/=p;
}
return res;
}
如果空间有限，无法存储所有阶乘，也可以只存储需要的阶乘，对它们进行排序，然后计算阶乘 \(0!,~ 1!,~ 2!,~ \dots,~ (p-1)!\) 而不显式存储它们。
Lucas
Lucas 定理内容如下：对于质数 \(p\)，有：
\[\dbinom{n}{m}\bmod p=\dbinom{\left\lfloor n/p\right\rfloor}{\left\lfloor m/p\right\rfloor}\cdot\dbinom{\left\lfloor n\bmod p\right\rfloor}{\left\lfloor m\bmod p\right\rfloor}\bmod p \]观察上述表达式，可知 \(n\bmod p\) 和 \(m\bmod p\) 一定是小于 \(p\) 的数，可以直接求解。\(\dbinom{\left\lfloor n/p\right\rfloor}{\left\lfloor m/p\right\rfloor}\) 可以继续用 Lucas 定理求解。这也就要求 \(p\) 的范围不能够太大，一般在 \(10^5\) 的级别。
边界条件：当 \(m=0\) 的时候，返回 \(1\)。
时间复杂度为 \(O(f(p) + g(n)\log n)\)，其中 \(f(n)\) 为预处理组合数的复杂度，\(g(n)\) 为单次求组合数的复杂度。
Lucas证明
- 考虑
\(\displaystyle\binom{p}{n}\bmod p\)的取值，注意到 \(\displaystyle\binom{p}{n}=\frac{p!}{n!(p-n)!}\)，分子的质因子分解中 \(p\) 的次数恰好为 \(1\)。
因此只有当 \(n = 0\) 或 \(n = p\) 的时候 \(n!(p-n)!\) 的质因子分解中含有 \(p\)，即 \(\displaystyle\binom{p}{n} \bmod p = [n = 0 \vee n = p]\)。进而我们可以得出：
\[\begin{aligned} (a+b)^p &= \sum_{n=0}^p \binom pn a^n b^{p-n}\\ &\equiv \sum_{n=0}^p [n=0\vee n=p] a^n b^{p-n}\\ &\equiv a^p + b^p \pmod p \end{aligned}\]注意过程中没有用到费马小定理，因此这一推导不仅适用于整数，亦适用于多项式。因此我们可以考虑二项式 \(f^p(x)=(ax^n + bx^m)^p \bmod p\) 的结果：
\[\begin{aligned} (ax^n + bx^m)^p &\equiv a^p x^{pn} + b^p x^{pm} \\ &\equiv ax^{pn} + bx^{pm}\\ &\equiv f(x^p) \end{aligned}\]考虑二项式 \((1+x)^n \bmod p\)，那么
\(\displaystyle\binom n m\) 就是求其在 \(x^m\) 次项的取值。使用上述引理，我们可以得到：
注意前者只有在 \(p\) 的倍数位置才有取值，而后者最高次项为 \(n\bmod p \le p-1\)，因此这两部分的卷积在任何一个位置只有最多一种方式贡献取值，即在前者部分取 \(p\) 的倍数次项，后者部分取剩余项，即：
\[\dbinom{n}{m}\bmod p=\dbinom{\left\lfloor n/p\right\rfloor}{\left\lfloor m/p\right\rfloor}\cdot\dbinom{\left\lfloor n\bmod p\right\rfloor}{\left\lfloor m\bmod p\right\rfloor}\bmod p \]Code
#include<bits/stdc++.h>
using namespace std;
typedef long long LL;
inline LL read(){
LL s=0,f=1;char ch=getchar();
while(ch<'0' or ch>'9'){if(ch=='-')f=-1;ch=getchar();}
while(ch>='0' and ch<='9'){s=(s<<1)+(s<<3)+(ch^48);ch=getchar();}
return f*s;
}
inline void write(LL num){
if(num<0)putchar('-'),num=-num;
if(num>9)write(num/10);
putchar(num%10+48);
}
const int N=1e5+5;
LL fac[N];
int T,a,b,m;
inline LL Pow(LL a,LL n,LL ma){
LL ans=1;a%=m;
while(n){
if(n&1)ans=(ans*a)%m;
a=(a*a)%m;
n>>=1;
}
return ans;
}
inline LL inverse(LL a,LL m){
return Pow(fac[a],m-2,m);
}
inline LL C(LL n,LL r,LL m){
if(r>n)return 0;
return (fac[n]*inverse(r,m))%m*inverse(n-r,m)%m;
}
inline LL lucas(LL n,LL r,LL m){
if(r==0)return 1;
return C(n%m,r%m,m)*C(n/m,r/m,m)%m;
}
int main(){
T=read();
while(T--){
a=read(),b=read(),m=read();
fac[0]=1;
for(int i=1;i<=m;i++)
fac[i]=(fac[i-1]*i)%m;
write(lucas(a+b,a,m));
puts("");
}
return 0;
}