同余&逆元

1. 同余

1. 同余的基本概念及性质

若$x$%$m=a$即m是 x-a 的一个因子, 则称x与a关于m同余,记作:$$x \equiv a(mod ;m)$$
同余基本性质:

○1. 自反性:$a \equiv a(mod\;m)$

○2. 对称性:$a \equiv b(mod\;m) \rightarrow b \equiv a(mod\;m)$

○3. 传递性:$a \equiv b(mod\;m),b \equiv c(mod\;m) \rightarrow a \equiv c(mod\;m)$

○4. 同加性:$a \equiv b(mod\;m) c \equiv d(mod\;m) \rightarrow a+c \equiv
b+d(mod\;m)$

○5. 同乘性:$a \equiv b(mod\;m) c \equiv d(mod\;m) \rightarrow ac \equiv bd(mod\;m)$

○6. 同幂性:$a \equiv b (mod\;m) \rightarrow a^n \equiv b^n(mod m)$n是自然数

○7. 若$a \equiv b(mod\;m),n|m$ 则 $a \equiv b(mod\;n)$

○8. 若$ac \equiv bc(mod\;m),(c,m)=d$则$a \equiv b(mod\;\dfrac{m}{d})$

○9. 若$(m,n)=1,a \equiv b(mod\;m),a \equiv b(mod\;n) \Leftrightarrow a \equiv
b(mod\;mn)$

2. 求解线性同余方程 $ax≡c(mod\;b)$

可转化为求解方程: $ax+by=c$

(同余方程和线性方程的关系很重要,经常用到!!)

1.预处理:

	if(a<0) a=-a,c=-c;

	while(c<0) c+=b;//保证 a,c 为正

2.第一步: 检验是否有解

	int gcd=Gcd(a,b);

	if((c%gcd)!=0) return -1;//若c不是gcd(a,b) 的倍数，则无解

	//可以转化为直线上的整点来理解

3.第二步:求解同余方程:$ax≡gcd(a,b)(mod b)$即$ax+by=gcd(a,b)$

扩展欧几里得算法

inline int ex_gcd(int a,int b)

{

	if(b==0) {x=1;y=0;return a;}

	int gcd=ex_gcd(b,a%b);

	int tmp=x;x=y;

	y=tmp-a/b*y;//扩欧

	return gcd;

}

//最后得到的x即为原同余方程的一个可行解;

扩欧的证明:

当最后 b'0 时a'gcd(a,b) 则此时 $a'x+b'y=gcd(a,b)*x$

令$x=1,y=0$即得最后的一组解。

考虑从后往前推出$ax+by=c$的解:

设我们前一步求出的解为$x_1,y_1,此时a,b的值就表示为a,b,当前的解表示为x,y$

那么因为$gcd(a,b)=gcd(b,a$%$b),a$%$b=a-(a/b)*b$有:$$bx_1+(a-(a/b)b)y_1=gcd(a,b)$$

则:

\[b*x_1+(a-(a/b)*b)y_1=ax+by
\]

整理得:

\[ay_1+b(x_1-(a/b)y_1)=ax+by
\]

容易看出:

\[x=y_1,y=x_1-(a/b)y_1
\]

即证。

4.第四步:根据题意得出答案

若要求出最小正整数解：

    while(x<0) x+=b;x%=b;

	b/=gcd;//mod 要变成 mod/gcd ;(mod 即为 b)

	x=x*c/gcd;//同余的同乘性质

    while(x<0) x+=b;x%=b;//最小正整数解

若要求出解的个数(或所有解)

    int tot=gcd(a,b);// 只有gcd(a,b) 个解

    //要求出每个解，只需对其不断加 b/gcd 即可(同时y-=a/gcd)

3.求解单变量模线性方程组(中国剩余定理)

有如下方程:

\[\begin{cases}
x \equiv a_1 (mod\;m_1)\\
x \equiv a_2 (mod\;m_2)\\
x \equiv a_3 (mod\;m_3)\\
\dots\dots \dots\\
x \equiv a_n (mod\;m_n)\\
\end{cases}
\]

其中$(m_1\;m_2\;m_3\dots m_n$两两互质$)$

为了方便表示,将x设为S

(1)设$M=\Pi^n_{i=1}m_i$，设$M_i=M/m_i$

(2)可知对于每一个$i有:(M_i,m_i)=1$

即:

$\qquad\qquad\qquad\qquad\qquad M_ix+m_iy=1$

那么$x为M_1$的逆元,用$t_i$表示

两边同时扩大$a_i倍$

$\qquad\qquad\qquad\qquad\qquad M_ia_it_i+m_ia_iy=a_i$

而$y$的取值与求解无关,可将$a_iy$视为$y$,则:

$\qquad\qquad\qquad\qquad\qquad M_ia_it_i+m_iy=a_i$

那么易知 $S=M_ia_it_i$

则原同余方程组通解为:

\[x=a_1t_1M_1+a_2t_2M_2+....+a_nt_nM_n+kM,k∈Z
\]

为什么把每一个加起来就行了呢?

因为每一个$M_i$都含有因子$m_j(j\ne i)$，对于其他的同余方程不产生影响。

若要求最小正整数解，则对$M$取模即可。

代码如下:

//中国剩余定理求解单变量模线性同余方程组

int CRT(int a[],int m[],int h)

{

	int ans=0;int M=1;

	for(int i=1;i<=h;i++) M*=m[i];//求M

	for(int i=1;i<=h;i++) {

		ll Mi,ti;

		Mi=M/m[i];//求Mi

		ti=Rev(Mi,m[i]);//求Mi的逆元

		ans+=((a[i]*Mi%M)*ti)%M;//累加答案

		if(ans>=M) ans-=M;//取模

	}

	return ans;

}

4.扩展中国剩余定理

这同样是用来解决单变量模线性方程组的,但是能够应用于模数不互质的情况

其实这个和中国剩余定理没有什么关系,CRT是用构造法,而EXCRT则基于扩展欧基里德算法

做法:

还是如下方程:

\[\begin{cases}
x≡a_1 (mod\;m_1)\\
x≡a_2 (mod\;m_2)\\
x≡a_3 (mod\;m_3)\\
\dots\dots \dots\\
x≡a_n (mod\;m_n)\\
\end{cases}
\]

其中$m_1\ m_2\ m_3\ m_4 \dots m_n$不一定互质

我们可以发现左边的式子都是相同的,于是有了同余方程合并这种操作

既然是合并,我们只要讨论两个式子的时候的情况

对于:

\[\begin{cases}
x \equiv a_1\ (mod\ m_1)\\
x \equiv a_2\ (mod\ m_2)\\
\end{cases}
\]

可以看做是两个方程:

\[\begin{cases}
x =a_1+x_1m_1\\
x =a_2+x_2m_2\\
\end{cases}
\]

合并一下得到:$$a_1+x_1m_1=a_2+x_2m_2$$

移项:$$x_1m_1=a_2-a_1+x_2m_2$$

假定$a_2 \geq a_1$,设为$c$,再化为同余方程:$$m_1x_1\equiv c\ (mod\ m_2)$$

令$gcd(m_1,m_2)=d$,该同于方程有解当且仅当$d|c$,所以如果$d$不整除$c$则整个同余方程组无解

反之,由同余的性质得:$$\frac{m_1}{d}x_1\equiv \frac{c}{d}\ (mod \frac{m_2}{d})$$

记$d_1=\dfrac{m_1}{d},d_2=\dfrac{m_2}{d},c_2=\dfrac{c}{d}$

由于此时$d_1,d_2$ 一定互质,所以$d_1$在模$d_2$的意义下一定有逆元,记为$d_1^{-1}$,那么可以解出$x_1$(其实就是扩欧)

\[x_1=d_1^{-1}c_2+d_2x_2
\]

回代进一开始的方程:

\[x=c+(d_1^{-1}c_2+d_2x_2)m_1
\]

展开化简得:

\[x=d_1^{-1}c_2m_1+c+\frac{m_1m_2}{d}x_2
\]

于是我们可以得到一个新的同余方程:

\[x\equiv d_1^{-1}c_2m_1+c \ (mod\ \frac{m_1m_2}{d})
\]

于是就这样一直合并下去,最后的解就直接出来了(注意最后的模数会变成$lcm(m_1,m_2,...,m_n)$)

中间结果注意防溢出

函数代码:

inline void EXCRT()

{

	ll p1,b1,p2,b2;scanf("%lld %lld",&p1,&b1);

	for(register int i=2;i<=n;++i) {

		scanf("%lld %lld",&p2,&b2);

		ll d=gcd(p1,p2);ll c=b2-b1;

		if(c%d) return void(puts("no solution"));

		ll d1=p1/d,d2=p2/d,lcm=p1/d*p2;// 这里根据的是负数也能取模 , 可以简化代码

		c/=d;ll inv,y;exgcd(d1,d2,inv,y);

		ll x1=Mul(inv,c,d2);

		b1=(b1+Mul(x1,p1,lcm))%lcm;p1=lcm;

	}

	printf("%lld\n",b1%p1);

	return;

}

5.卢卡斯定理(大组合数取模)

对于组合数取模,如$$C^m_n\ mod\ p$$

其中p是一个质数,有如下定理:

卢卡斯定理:组合数$C^m_n$在模意义下等价于把n和m看成一个p进制数,对每一位分别求出组合数后乘起来

比如说假设:

$n=a_1*p^0+a_2*p^1+a_3*p^3+\dots a_k*p^k,m=b_1*p^0+b_2*p^1+b_3*p^3+\dots b_k*p^k$

那么:$$C^m_n; mod; p=\prod_{i=0}^k C^{b_i}_{a_i} ; mod;p$$

显然如果p很大的话没有什么鸟用

但是当p不是特别大的话,我们可以发现通过这个定理我们要求的组合数的n,m都不会超过p,可以使用阶乘来解决,并且这时阶乘一定和p是互质的,一定存在逆元,通过阶乘逆元我们可以直接算出组合数

复杂度是$O(p\ log_pn)$,预处理阶乘逆元就直接是$O(p+log_pn)$了,看上去还是很有用的

主要代码:

inline ll C(ll n,ll m,ll p){

	if(m>n) return 0;

	return fac[n]*fpow(fac[m],p-2,p)%p*fpow(fac[n-m],p-2,p)%p;

}

ll Lucas(ll n,ll m,ll p)

{

	if(!m) return 1;

	return C(n%p,m%p,p)*Lucas(n/p,m/p,p)%p;

}

6.扩展卢卡斯定理

还是这个东西:$$C^m_n\ mod\ p$$

但是p不一定是质数

这个其实和卢卡斯定理也没有什么很大的关系

我们先把p给质因数分解:$$p=p_1^{k_1}p_2{k_2}p_3^{k_3}\dots p_n^{k_n}$$

可以看出,如果所有的k都是1的话,我们设$C_n^m=x$,对每一个$p_i$分别用卢卡斯定理求出组合数$C_i$,那么就变成了求解一系列的同余方程组:

\[\begin{cases}
x \equiv C_1 \ mod\ p_1 \\
x \equiv C_2 \ mod\ p_2 \\
x \equiv C_3 \ mod\ p_3 \\
x \equiv C_4 \ mod\ p_4 \\
x \equiv C_5 \ mod\ p_5 \\
\end{cases}
\]

于是我们可以用中国剩余定理进行合并,求出最后的x,显然在模p意义下最后只会有唯一解

问题在于我们现在的p可能不是一次,而是有k次,要想用CRT来进行合并只能靠求出$C_n^m\ mod\ p_i^{k_i}$

因为要把同余式$x\equiv a\ (mod\ p_1p_2)$拆开必须要满足$p_1$和$p_2$互质,显然两个相同的数不会互质

于是问题转化为快速求出$C_n^m\ mod\ p_i^{k_i}$

为方便,我们设现在考虑的模数是$p^k$,$p$是一个质数

还是考虑用阶乘来解决:$$C^m_n=\dfrac{n!}{m!(n-m)!}$$

我们仔细观察可以发现$n!$中含有的$p$这个因子的个数一定会不少于$m!(n-m)!$中p的个数,我们可以很容易得到一个阶乘中含有的p的个数的递推公式:$$f(n)=f(\lfloor { \frac{n}{p} }\rfloor)*\lfloor { \frac{n}{p} }\rfloor$$

例如我们要求:$9!$中$2$的个数:

\[9!=1\times2\times3\times4\times5\times6\times7\times8\times9 \\
=1\times3\times5\times7\times9\times[2\times(1\times2\times3\times4)]
\]

这就比较直观了,有了这个的话,我们假设求出阶乘中不含$p$的项的积,这样就可以通过逆元来进行组合数计算了,只需要最后把因该有的$p$给再乘上去就行了

于是问题再次变为如何快速求出阶乘

其实方法在上面$9!$的变换中就可以发现了,由于我们不管$p$有多少,发现提出来一个$p$之后,右边那部分的阶乘可以递归进行计算,就是$\lfloor { \frac{n}{p} }\rfloor!$,于是关键在于计算左边

由于是模p意义下,在上面的式子中,可以发现左边其实全部都是1,手玩一下其他的发现显然这个东西是以p个一循环的,并且可能会剩下不超过p个数

所以循环部分算出一个然后快速幂$n/p$次,最后还会剩下$n\%p$个,直接暴力算这些

所以对于一次的阶乘要算的次数也不会超过$O(p)$次,总共有$log_pn$层

所以计算一个阶乘的复杂度为$O(p\ log_pn)$

总复杂度也就是把所有模数的复杂度加起来,最高也不超过最大质因子的复杂度

所以我们就解决了这个问题

剩下的就是算出逆元,求出组合数,处理多余的质因子p,然后CRT合并

代码:

#include<bits/stdc++.h>

using namespace std;

typedef long long ll;

ll n,m;int p;

inline void exgcd(int a,int b,int&x,int &y){

	if(!b) {x=1;y=0;return;}exgcd(b,a%b,x,y);

	int tmp=x;x=y;y=tmp-a/b*y;return;

}

template<class T>inline int fpow(int x,T k,int mod){int ret=1;for(;k;k>>=1,x=(ll)x*x%mod) if(k&1) ret=(ll)ret*x%mod;return ret;}

template<class T>inline void Inc(T&x,int y,int mod){x+=y;if(x>=mod) x-=mod;return;}

inline int Fac(ll n,int pi,int pk){

	if(n<=1) return 1;

	int ret=1,rsub=Fac(n/pi,pi,pk),res=n%pk;

	if(n>=pk) {for(int i=2;i<=pk;++i) if(i%pi) ret=(ll)ret*i%pk;ret=fpow(ret,n/pk,pk)%pk;} // 统计长度为模数且不含该模数质因子的阶乘

	for(int i=2;i<=res;++i) if(i%pi) ret=(ll)ret*i%pk;//模意义下,直接枚举模了之后的未统计数也可以

	return (ll)ret*rsub%pk;

}

inline int Inv(int a,int b) {int x,y;exgcd(a,b,x,y);x=(x+b)%b;return x;}

inline ll Count(ll n,int d){return n<d? 0:(Count(n/d,d)+n/d);}

inline int C(ll n,ll m,int pi,int pk)

{

	if(m>n) return 0;

	int facn=Fac(n,pi,pk),facm=Fac(m,pi,pk),facmn=Fac(n-m,pi,pk);//求解三个阶乘

	ll num=Count(n,pi)-Count(m,pi)-Count(n-m,pi);//把p这个质因子都提出来单独算(其实算阶乘的时候没有处理这些数)

	return (ll)facn*Inv(facm,pk)%pk*Inv(facmn,pk)%pk*fpow(pi,num,pk)%pk;

}

inline int EXLucas(ll n,ll m,int p){

	if(n<m) return 0;if(n==m||!m) return 1;

	int ans=0;int x=p;

	for(int i=2;i<=p;++i)

		if(!(x%i)){

			int pk=1;while(!(x%i)) pk*=i,x/=i;

			int res=C(n,m,i,pk);

			Inc(ans,(ll)res*(p/pk)%p*Inv(p/pk,pk)%p,p);

		}return ans;

}

int main(){scanf("%lld %lld %d",&n,&m,&p);printf("%d\n",EXLucas(n,m,p));}

Upd: 稍微改了一下上面的模板 , 下面的模板是预处理阶乘后的 , 这样 $p$ 这部分的复杂度就不用带 $log$ 了。

int FC[4][N],Pr[4]={0,3,11,100003},mo[4]={0,81,121,100003};//这里用于预处理到模数(不含其质因子)的阶乘,存放质因子和分解后模数

inline int gcd(int a,int b){return b? gcd(b,a%b):a;}

inline int Count(int n,int d){return n<d? 0:(n/d+Count(n/d,d));}//计算阶乘中的该因子个数

inline void Exgcd(int a,int b,int &x,int &y){

	if(!b) {x=1;y=0;return;}

	Exgcd(b,a%b,x,y);int tmp=x;x=y;y=tmp-a/b*y;return;

}

inline int Inv(int a,int mod){a%=mod;int x,y;Exgcd(a,mod,x,y);return (x+mod)%mod;}

inline int Fac(int n,int id){//阶乘

	return n<=Pr[id]? FC[id][n]:((ll)fpow(FC[id][mo[id]],n/mo[id],mo[id])*FC[id][n%mo[id]]%mo[id]*Fac(n/Pr[id],id))%mo[id];

}

inline int Comb(int n,int m,int id){//组合数

	int facn=Fac(n,id),facm=Fac(m,id),facnm=Fac(n-m,id);

	int Pop=fpow(Pr[id],Count(n,Pr[id])-Count(m,Pr[id])-Count(n-m,Pr[id]),mo[id]);

	if(!Pop) return 0;

	return (ll)facn*Inv(facm,mo[id])%mo[id]*Inv(facnm,mo[id])%mo[id]*Pop%mo[id];

}

inline int ExLucas(int n,int m){

	if(n<m) return 0;if(!n&&!m) return 1;int ret=0;

	for(int i=1;i<=3;++i) {int C=Comb(n,m,i);Inc(ret,(ll)C*Inv(mod/mo[i],mo[i])%mod*(mod/mo[i])%mod);}

	return ret%mod;

}

//压行真的好看多了

7.二次剩余

求解 $x$ 使得:

\[x^2\equiv n \ (mod\ p)
\]

$p$ 如果是 $2$ 那么 $x=n$ 一定是 $n$的一个二次剩余 , 所以这里讨论的 $p$都是奇质数。

首先二次剩余并不一定存在 , 其存在的充要条件是: $n^{\frac{p-1}{2}} \equiv 1\ (mod\ p)$

证明(以下运算均在模意义下):

首先我们有 $x^{p-1}=1$

把要求的东西左右各 $\frac{p-1}{2}$ 次方那么就是 $n^{\frac{p-1}{2}}=x^{p-1}=1,$证毕。

求解方法:

法1:

我们可以利用质数 $p$ 的原根 $g$ 来解决这个问题

我们找到一个 $d$ 使得，$g^d=n$，因为 $g$ 是 $p$ 的原根那么一定是有解的，可以使用 $BSGS$ 求解。

那么问题变成 $x^2=g^d$ ，可以证明 $d$ 一定是一个偶数，所以 $x=g^{\frac{d}{2}}$ ，就做完了。

复杂度为 $O(\sqrt p)$ ，因为要使用 $BSGS$ 算法求解离散对数。

法2:

这就是一种很数学的方法了。

我们随机一个数 $a$ 满足 $a^2-n$ 没有二次剩余，这样的期望次数为 2 。

然后令 $\delta=\sqrt{a^2-n}$，定义一个新的数域，所有数可以表示为 $a+b\delta$

结论是 $(a+\delta)^{\frac{p+1}{2}}$ 就是 $n$ 的二次剩余。

证明:

因为 $a^2-n$ 没有二次剩余，所以 $(a^2-n)^{\frac{p-1}{2}}\neq 1$

由费马小定理得到 $(a^2-n)^{p-1}=1$ , 那么由代数基本定理这个形如 $x^2=1$ 的方程只有 $+1,-1$ 两个根，所以 $(a^2-n)^{\frac{p-1}{2}}=-1$

也就是 $\delta^{p-1}=-1$

然后考虑 $(a+\delta)^{p}$

\[(a+\delta)^p=\sum_{i=0}^p a^i\delta^{p-i}{p\choose i}
\]

因为 $p$ 是一个奇质数，所以组合数在$\mod p$ 意义下时当且仅当 $i=0$ 或 $i=p$ 时才不为0

所以:$(a+\delta)^p=a^p+\delta^p$

因为 $a^p=a,\delta^p=-\delta$

$(a+\delta)^p=a-\delta$

$(a+\delta)^{p+1}=a-\delta$

代码：

inline int Get_sqrt(int n){//二次剩余

	if(n<=1) return n;if((int)sqrt(n)*(int)sqrt(n)==n) return (int)sqrt(n);

	srand(time(NULL));int a,delta,D=(mod-1)>>1;

	while(233) {a=rand()%mod;delta=Dif((ll)a*a%mod,n);if(fpow(delta,D)!=1) break;}

	D=(mod+1)>>1;int b=1;int c=1,d=0;

	while(D){

		if(D&1) {int _c=c;c=((ll)a*c+(ll)b*d%mod*delta)%mod,d=((ll)a*d+(ll)b*_c)%mod;}

		int _a=a;a=((ll)a*a+(ll)b*b%mod*delta)%mod;b=(2ll*b*_a)%mod;D>>=1;

	}c=min(c,mod-c);

	return c;

}

2.逆元

1.定义:在 mod m 的意义下,设b是a的逆元,则:$a*b≡1 (mod\;m)$

2.求解方法:

1.线性递推法

逆元的递推公式:$$inv[i]=(P-P/i)*inv[P%i]%P$$

阶乘逆元的递推公式:$$inv[i]=inv[i+1]*(i+1)%P$$

上面的$inv[i]$表示 $i!$ 的逆元

2.扩展欧基里德法:

在上面我们用扩欧求了如下同余方程的解:$$ax \equiv b\ (mod\ p)$$

而逆元是求的这个同余方程:$$ax \equiv 1\ (mod\ p)$$

所以把b直接设为1然后扩欧解出x的最小正整数解就是a的逆元了

可以发现a,p一定要互质,不然同余方程无解,即a在a和p不互质的情况下是没有逆元的

3.费马小定理

费马小定理:

当a和p互质且p是质数时,有

\[a^{p-1} \equiv 1 \ (mod\ p)
\]

相当于是:$$a^{p-2}*a \equiv 1\ (mod\ p)$$

所以$a^{p-2}$就是a的逆元了

4.欧拉定理

欧拉定理:

当a和p互质时,有

\[a^{\varphi(p)} \equiv 1\ (mod\ p)
\]

所以$a$的逆元是$a^{\varphi(p)-1}$,其实费马小定理是欧拉定理在p是质数时的一个特殊情况

补充:

扩展欧拉定理:

\[a^b \equiv\begin{cases}
a^{b \%\varphi(p)+\varphi(p)}\ (mod\ p) & ,b> \varphi(p)\\
a^b & ,b\leq \varphi(p)\\
\end{cases}
\]

证明不会

小结:逆元在a和p互质的情况下才有

巴特西

同余&逆元简单总结