学习笔记——数学
死也学不会系列。
数学题单 1
1、质数相关
1.筛法求质数
埃筛
当我们找出一个质数时，将它在范围内的所有倍数全部标记，因为小于它的数已经处理过，所以无需标记，复杂度 \(O(n \log\log n)\)。
inline void get_pr1(int n){
memset(vis,0,sizeof(vis));
memset(pr,0,sizeof(pr));
for(int i=2;i<=n;i++){
if(vis[i]) continue;
pr[++pr[0]]=i;
for(int j=1;i*j<=n;j++){
vis[i*j]=1;
}
}
}
线性筛
埃筛的不足在于，一个质因数较多的数会被多次标记到，因此我们把标记由 \(0/1\) 改为最小的质因数，就可以避免此情况的出现。
inline void get_pr2(int n){
memset(vis,0,sizeof(vis));
memset(pr,0,sizeof(pr));
for(int i=2;i<=n;i++){
if(!vis[i]){
vis[i]=i;
pr[++pr[0]]=i;
}
for(int j=1;j<=pr[0];j++){
if(i*pr[j]>n||pr[j]>vis[i]) break;
vis[i*pr[j]]=pr[j];
}
}
}
2.唯一分解定理
对于一个正整数 \(N\)，一定有：
\[N=p_1^{c_1}p_2^{c_2}\dots p_m^{c_m} \]其中 \(p_i\) 是质数，\(c_i\) 是非负整数。
由唯一分解亦可得：
\(N\) 正约数的个数：
\[(c_1+1)\times (c_2+1)\times \dots \times(c_m+1)=\prod_{i=1}^m(c_i+1) \]\(N\) 正约数的和：
\[(1+p_1+p_1^2+\dots+p_1^{c_1})\times \dots\times(1+p_m+p_m^2+\dots +p_m^{c_m})=\prod_{i=1}^m \sum_{j=0}^{c_i} p_i^j \]2、欧几里得算法相关
1.欧几里得算法
用来求最大公约数和最小公倍数。
更相减损术
\(\gcd(a,b)=\gcd(a,a-b)=\gcd(b,a-b)\)
在处理大整数的最大公约问题时应用较多。
辗转相除法
\(\gcd(a,b)=\gcd(b,a\bmod b)\)
证明：
若 \(a
若 \(a\ge b\)，设 \(a=p\times b +q\)，此时 \(q=a\bmod b\)，对于任意一个公约数 \(d\)，都有 \(d\mid a\) 且 \(d\mid b\)，所以就有 \(d\mid (p\times b)\)，即 \(d\mid (a-p\times b)\)，也就是 \(d\mid q\)，此情况成立。
inline int gcd(int a,int b){
if(!b) return a;
else return gcd(b,a%b);
}
最小公倍数
设两正整数 \(a=p_1^{k_1}p_2^{k_2}\dots p_m^{l_m},b=p_1^{l_1}p_2^{l_2}\dots p_m^{l_m}\)
由唯一分解可得：
\[\gcd(a,b)=p_1^{\min(k_1,l_1)}p_2^{\min(k_2,l_2)}\dots p_m^{\min(k_m,l_m)} \]\[\operatorname{lcm}(a,b)=p_1^{\max(k_1,l_1)}p_2^{\max(k_2,l_2)}\dots p_m^{\max(k_m,l_m)} \]\[a\times b=p_1^{k_1+l_1}p_2^{k_2+l_2}\dots p_m^{k_m+l_m} \]因此有：
\[\operatorname{lcm}(a,b)=\dfrac{a\times b}{\gcd(a,b)} \]多个数求gcd以及lcm
\[\gcd(a,b,c)=\gcd(\gcd(a,b),c) \]\[\operatorname{lcm}(a,b,c)=\operatorname{lcm}(\operatorname{lcm}(a,b),c) \]求 \(\operatorname{lcm}\) 时，先求出前两个数的 \(\operatorname{lcm}\)，再与下一个数一起求 \(\operatorname{lcm}\)，以此类推。
2.扩展欧几里得算法
用于求 \(ax+by=\gcd(a,b)\) 解。
求一组可行解
由欧几里得算法得，\(\gcd(a,b)=\gcd(b,a\bmod b)\)
所以列含参方程得：
\[\begin{cases} ax_1+by_1=\gcd(a,b)\\bx_2+(a\bmod b)y_2=\gcd(b,a \bmod b)\end{cases} \]于是有:
\[ax_1+by_1=bx_2+(a\bmod b)y_2 \]\[ax_1+by_1=bx_2+(a-\lfloor \dfrac{a}{b} \rfloor\times b)y_2 \]\[ax_1+by_1=ay_2+b(x_2-\lfloor \dfrac{a}{b} \rfloor\ y_2) \]所以：\(x_1=y_2,y_1=(x_2- \lfloor \frac{a}{b}\rfloor y_2)\)，不断递归求解可以得到。
inline int exgcd(int a,int b,int &x,int &y){
if(!b){
x=1,y=0;
return a;
}
int d=exgcd(b,a%b,x,y);
int t=x;
x=y,y=t-a/b*y;
return d;
}
对于上述的方程是一定有解的，并且存在一组 \((x,y)\) 满足 \(|x|\le b,|y| \le a\)，我们称其为裴蜀定理。
求任意一组通解
设已经求出的可行解为 \((x_0,y_0)\)，带入到方程中有：
\[\begin{cases} ax_0+by_0=\gcd(a,b)\\ax+by=\gcd(a,b)\end{cases} \]联立并整理得：
\[a(x_0-x)=b(y-y_0) \]\[\dfrac{a}{\gcd(a,b)}(x_0-x)=\dfrac{b}{\gcd(a,b)}(y-y_0) \]容易发现这个等式左右两边的前一个因式是互质的，因此后一个因式与前一个因式是对应整除关系（根据唯一分解定理可以得到），设这个商为 \(k\)，就有：
\[\begin{cases}x_0-x =\dfrac{b}{\gcd(a,b)}\times k\\\\y-y_0=\dfrac{a}{\gcd(a,b)}\times k\end{cases} \]整理得到：\(x=x_0-\frac{b}{\gcd(a,b)}\times k,y=y_0+\frac{a}{\gcd(a,b)}\times k\)，其中 \(k\in \mathbb{Z}\)。
推广到不定方程
对于不定方程 \(ax+by=c\)，解决方式如下：
求解 \(ax+by=\gcd(a,b)\) 的一组解 \((x_0,y_0)\)，同乘 \(\frac{c}{\gcd(a,b)}\)，就能得到不定方程的一组可行解 \((x,y)=(x_0\times\frac{c}{\gcd(a,b)},y_0\times \frac{c}{\gcd(a,b)})\)，而我们把它代入到通解的式子中，得到最终的结果：
\[\begin{cases} x=\dfrac{c}{\gcd(a,b)}\times x_0- \dfrac{b}{\gcd(a,b)}\times k \\\\ y=\dfrac{c}{\gcd(a,b)}\times y_0+ \dfrac{a}{\gcd(a,b)}\times k\end{cases} \]特别地，\(x\) 的最小整数解为 \((x \bmod \frac{b}{\gcd(a,b)}+\frac{b}{\gcd(a,b)})\bmod\frac{b}{\gcd(a,b)}\)
例题
1.P4549 裴蜀定理
其实就是要求一个 \(a_1x_1+a_2x_2+\dots+a_nx_n\) 的最小值，发现前两项和的最小值为 \(\gcd(a_1,a_2)\)，因此一边读入一边求 \(\gcd\) 即可。
inline ll gcd(ll a,ll b){
if(!b) return a;
else return gcd(b,a%b);
}
int main(){
n=read();
for(int i=1;i<=n;i++){
int a=read();
ans=gcd(ans,llabs(a));
}
printf("%lld\n",ans);
return 0;
}
2.P5656 二元一次方程（exgcd）
需要推式子……
设 \(a'=\frac{a}{\gcd(a,b)},b'=\frac{b}{\gcd(a,b)}\)，为保证 \(x \ge 1\)，需要有 \(x_0-b'k\ge1\)，得到：
\[k\ge \lceil\frac{1-x_0}{b'}\rceil,x_{min}=x_0+\lceil\frac{1-x_0}{b'}\rceil\times b' \]其对应的 \(k\) 可以得到 \(y_{max}\)，而\(y_{min}=y_{max}\bmod a'\)，对应可以得到 \(x_{max}\)，总之就是很“模板”。
3.求解同余方程
方程 \(ax+by=c\) 与同余方程 \(ax \equiv c\pmod b\) 是等价的，因此求解只需要按照上面的方法就可以了。
例题
1.P1082 同余方程
板子……
inline int exgcd(int a,int b,int &x,int &y){
if(!b){
x=1,y=0;
return a;
}
int d=exgcd(b,a%b,x,y);
int t=x;
x=y,y=t-a/b*y;
return d;
}
int main(){
int a=read(),b=read();
int x=0,y=0;
int d=exgcd(a,b,x,y);
x=(x+b)%b;
printf("%d\n",x);
return 0;
}
2.P1516 青蛙的约会
因为是首尾相连的路径所以题目是在求一个：
\[pos_1+sp_1\times k \equiv pos_2+sp_2\times k\pmod l \]我们整理一下，可以得到：
\[(sp_1-sp_2)\times k-l\times q=(pos_2-pos_1) \]其中 \(q\in \mathbb{Z}\)，设 \(a=sp_1-sp_2,b=l,c=pos_2-pos_1\)，就得到了方程：
\[ak+bq=c \]求 \(k\) 的最小整数解即可。
inline ll exgcd(ll a,ll b,ll &x,ll &y){
if(!b){
x=1,y=0;
return a;
}
ll d=exgcd(b,a%b,x,y);
ll t=x;
x=y,y=t-a/b*y;
return d;
}
ll pos1,pos2,sp1,sp2,l;
int main(){
pos1=read(),pos2=read(),sp1=read(),sp2=read(),l=read();
ll a=sp1-sp2,b=l,c=pos2-pos1,x=0,y=0;
if(a<0) a=-a,c=-c;
ll d=exgcd(a,b,x,y);
if(c%d){
printf("Impossible\n");
}
else{
printf("%lld\n",((x*c/d)%(b/d)+(b/d))%(b/d));
}
return 0;
}
3.P1292 倒酒
4.乘法逆元
如果有 \(ax\equiv 1\pmod b\)，就称 \(x\) 是 \(a\bmod \ b\) 的逆元，记作 \(a^{-1}\)。
扩展欧几里得求逆元
因为 \(ax\equiv 1\pmod b\) 等价于 \(ax+by=1\)，因此可以解这个不定方程，来得到逆元的值。
费马小定理
若 \(p\) 为质数，且 \(\gcd(a,p)=1\)，则有：\(a^{p-1}\equiv 1\pmod p\)，也可写作：\(a^p\equiv a\pmod p\)。
快速幂求逆元
由费马小定理可得：
\[ax\equiv 1\equiv a^{b-1} \pmod b \]此时，\(b\) 是质数，所以可以得到：
\[x \equiv a^{b-2}\ \pmod b \]用快速幂可以 \(O(logn)\) 求出一个数的逆元。
线性求 \([1,n]\) 的逆元
因为 \(1\) 在任何模定义下的逆元都为 \(1\)，所以可以根据这个性质递推得解。
已知同余方程 \(i\times x\equiv 1\pmod p\)，其中 \(x\) 的值为 \(\bmod \ p\) 下 \(i^{-1}\) 的值。
设 \(k=\lfloor\frac{p}{i}\rfloor,j=p \bmod i,p=ki+j\)，所以就有同余方程 \(ki+j\equiv0 \pmod p\)。
方程两边同乘 \(i^{-1}\times j^{-1}\) 并移项可得：\(i^{-1}\equiv-kj^{-1} \pmod p\)
将 \(k=\lfloor\frac{p}{i}\rfloor,j=p \bmod i\) 代入得到：\(i^{-1}\equiv-\lfloor\frac{p}{i}\rfloor\times (p\bmod i)^{-1} \pmod p\)
因为 \(p\bmod i\) 的逆元是已知的，直接代入求解即可。
求任意 \(n\) 个的逆元
我们先用求出 \(pord(i)\)，表示前 \(i\) 个数的前缀积，用快速幂求出 \(prod(n)\) 的逆元，设为 \(pinv(n)\)，也就是 \(pinv(n)= prod(n)^{b-2}\)，接下来对于每一个 \(i\ (i\in[1,n-1])\)，\(pinv(i)=pinv(i+1)\times a(i+1)\)，而每个数的逆元就是 \(inv(i)=pinv(i)\times prod(i-1)\)。
ll n,b;
ll a[maxn],prod[maxn],prod_inv[maxn],inv[maxn];
inline ll q_pow(ll x,ll p){
ll ans=1;
x=(x%b+b)%b;
while(p){
if(p&1){
ans=ans*x%b;
}
x=x*x%b;
p>>=1;
}
return ans;
}
int main(){
n=read(),b=read();
prod[0]=1;
for(int i=1;i<=n;i++){
a[i]=read();
prod[i]=prod[i-1]*a[i]%b;
}
prod_inv[n]=q_pow(prod[n],b-2);
for(int i=n-1;i>=1;i--){
prod_inv[i]=prod_inv[i+1]*a[i+1]%b;
}
for(int i=1;i<=n;i++){
inv[i]=prod[i-1]*prod_inv[i]%b;
printf("%lld\n",inv[i]);
}
return 0;
}
例题
1.P3811 乘法逆元
模板题，代码同上。
2.P5431 乘法逆元2
同上上，需要注意一下卡常。
3.P2054 AHOI2005 洗牌
我们通过找规律可以发现形如 \({1,2,3,4,5,6,7,8,9,10}\) 的数列在经过一次洗牌后变作 \({6,1,7,2,8,3,9,4,10,5}\)，也就是对于洗牌前的位置 \(i\)，洗牌后的 \(i'\) 为：
\[i'=\begin{cases} i\times 2 \ (i\le\dfrac{n}{2})\\i\times 2-(n+1) \ (i>\dfrac{n}{2})\end{cases} \]也就有：
\[i'=(i\times 2)\bmod (n+1) \]那么根据题意求 \(n\) 纸牌在进行 \(m\) 次变换后第 \(l\) 张的初始位置（初始位置和其数值是等价的），设这个初始位置为 \(x\)，根据式子可以得到：
\[x\times 2^m \equiv\ l\pmod {n+1} \]为了得到 \(x\)，我们需要在方程两边乘上 \(2\) 在 \(\bmod \ n+1\) 定义下的逆元 \(inv\)，来约去 \(2^m\)，因此就是：
\[x\times 2^m\times inv^m\equiv x\equiv l\times inv^m\pmod {n+1} \]求这个式子即可，注意这里不保证模数 \(n+1\) 是质数，因此不能用快速幂求逆元，同时需要使用快速幂与龟速乘结合。
inline ll exgcd(ll a,ll b,ll &x,ll &y){
if(!b){
x=1,y=0;
return a;
}
ll d=exgcd(b,a%b,x,y);
ll t=x;
x=y,y=t-a/b*y;
return d;
}
inline ll s_mul(ll x,ll p){
ll ans=0;
while(p){
if(p&1){
ans=(ans+x)%(n+1);
}
x=(x+x)%(n+1);
p>>=1;
}
return ans%(n+1);
}
inline ll q_pow(ll x,ll p){
ll ans=1;
while(p){
if(p&1){
ans=s_mul(ans,x);
}
x=s_mul(x,x);
p>>=1;
}
return ans%(n+1);
}
ll inv,ans,tmp;
int main(){
n=read(),m=read(),l=read();
exgcd(2,n+1,inv,tmp);
inv=(inv%(n+1)+n+1)%(n+1);
inv=q_pow(inv,m);
ans=s_mul(l,inv);
printf("%lld\n",ans);
return 0;
}
4.P2613 有理数取余
求 \(\dfrac{a}{b} \bmod 19260817\)，并不需要大整数，因为是逆元的问题，在快读时一边读入一边取模就可以。
inline ll mread(){
ll x=0,w=1;char c=getchar();
while(c<'0'||c>'9'){if(c=='-')w=-1;c=getchar();}
while(c<='9'&&c>='0'){x=x*10%mod+(c-'0')%mod;c=getchar();}
return x%mod*w;
}
3、欧拉定理相关
1.欧拉函数
\(\varphi(N)\) 表示小于等于 \(N\) 且与 \(N\) 互质的正整数个数，即
\[\varphi(N)=\sum_{i=1}^{N}\ [\gcd(i,N)=1] \]特别地，\(\varphi(p)=p-1\)，\(p\) 为质数。
由容斥原理可以得到：
\[\varphi(N)=N\times \dfrac{p_1-1}{p_1}\times \dfrac{p_2-1}{p_2} \times \dots \times \dfrac{p_m-1}{p_m}=N\times \prod_{p\mid N}(1-\dfrac{1}{p}) \]inline int euler1(int n){
int ans=n;
for(int i=2;i<=sqrt(n);i++){
if(n%i==0){
ans=ans/i*(i-1);
while(n%i==0){
n/=i;
}
}
}
if(n>1){
ans=ans/n*(n-1);
}
return ans;
}
欧拉函数具有积性函数的性质。
若 \(\gcd(a,b)=1\)，则 \(f(a)\times f(b)=f(a\times b)\)，则 \(f\) 是积性函数，而对于 \(\varphi(a)\times \varphi(b)\) 而言，因为 \(a\) 与 \(b\) 互质，带入到求欧拉函数的公式当中后，每一项 \(\frac{p-1}{p}\) 都不重复，所以性质成立。
求 \(n\) 个数的欧拉函数
int phi[maxn];
inline void get_euler1(int n){
for(int i=1;i<=n;i++){
phi[i]=i;
}
for(int i=2;i<=n;i++){
if(phi[i]==i){
for(int j=1;i*j<=n;j++){
phi[i*j]=phi[i*j]/i*(i-1);
}
}
}
}
例题
1.P2158 SDOI2008 仪仗队
不难发现，当且仅当 \(\gcd(x,y)=1\) 时，\((x,y)\) 能够被看到，因此答案应为 \(1+2\times \sum_{i=1}^{n-1} \varphi(i)\)，筛出欧拉函数值即可。
int phi[maxn];
inline void get_euler(int n){
for(int i=1;i<=n;i++){
phi[i]=i;
}
for(int i=2;i<=n;i++){
if(phi[i]==i){
for(int j=1;i*j<=n;j++){
phi[i*j]=phi[i*j]/i*(i-1);
}
}
}
}
int ans;
int main(){
n=read()-1;
if(!n){
printf("0\n");
return 0;
}
get_euler(n);
for(int i=1;i<=n;i++){
ans+=phi[i];
}
ans=ans*2+1;
printf("%d\n",ans);
return 0;
}
2.P2303 SDOI2012 Longge的问题
\[\begin{aligned} \sum_{i=1}^n\gcd(i,n)&= \sum_{d\mid n}d\cdot \sum_{i=1}^n[\gcd(i,n)=d]\\&= \sum_{d\mid n}d\cdot \sum_{i=1}^{n/d}\left[\gcd\left(i,\dfrac nd\right)=1\right]\\&= \sum_{d\mid n}d\cdot\varphi\left(\dfrac nd\right) \end{aligned}\]3.P2155 SDOI2008 萨拉公主的困惑
要求 \([1,n!]\) 中整数与 \(m!\ (m\le n)\) 互质的个数。
先考虑 \([1,m!]\) 的部分：
\[\begin{aligned} \sum_{i=1}^{m!}[\gcd(i,m!)=1]=\varphi(m!) \end{aligned}\]由欧几里得算法得：
\[\gcd(a,b)=\gcd(a\bmod b,b) \]所以有：
\[\gcd(m!\times k+i,m!)=\gcd(i,m!) \]因此 \([m!\times (k-1)+1,m!\times k]\ (1\le k \le \frac{n!}{m!})\) 与 \([1,m!]\) 的部分都是等价的，所以有:
\[\begin{aligned} \sum_{i=1}^{n!} [\gcd(i,m!)=1]&=\dfrac{n!}{m!}\times \sum_{i=1}^{m!} [\gcd(i,m!)=1]\\&=\dfrac{n!}{m!}\times \varphi(m!) \end{aligned}\]现在考虑 \(\varphi(m!)\) 怎么求，已知 \(\varphi(i)\) 的递推式：
\[\varphi(i)=i\times \prod_{p\mid i}\dfrac{p-1}{p} \]所以可以得到：
\[\varphi(i!)=\begin{cases} i\times (i-1)!\times \prod_{p\mid i!} \left(\frac{p-1}{p}\right)\times \frac{i-1}{i} &(i\in \mathbb{P})\\\\ i\times (i-1)!\times \prod_{p\mid i!}\left(\frac{p-1}{p}\right) &(i\not\in\mathbb{P}) \end{cases}\]也就是说 \(i\) 是否为素数只会影响到求积中是否会添加一项，化简就有：
\[\varphi(i!)=\begin{cases} \varphi((i-1)!)\times (i-1) &(i\in \mathbb{P})\\\\ \varphi((i-1)!)\times i &(i\not\in \mathbb{P}) \end{cases}\]需要注意的是题目中的范围没有明确给出 \(R>n\)，所以可能存在 \(R\mid n!\)，使逆元不存在，因此需要处理时模去 \(R\)，并且特判是否答案为 \(0\)。
inline void get_pr(){
vis[0]=vis[1]=1;
for(int i=2;i<=base;i++){
if(!vis[i]){
pr[++pr[0]]=i;
}
for(int j=1;j<=pr[0];j++){
if(i*1ll*pr[j]>base) break;
vis[i*pr[j]]=1;
if(i%pr[j]==0) break;
}
}
}
inline ll q_pow(int x,int p){
ll ans=1;
while(p){
if(p&1){
ans=ans*1ll*x%r;
}
x=x*1ll*x%r;
p>>=1;
}
return ans;
}
int main(){
t=read(),r=read();
get_pr();
fac[0]=1,fphi[0]=1;
for(int i=1;i<=base;i++){
int tmp=i;
while(tmp%r==0){
tmp/=r;
}
fac[i]=fac[i-1]*1ll*tmp%r;
tmp=vis[i]?i:i-1;
while(tmp%r==0){
tmp/=r;
}
fphi[i]=fphi[i-1]*1ll*tmp%r;
}
while(t--){
n=read(),m=read();
if(n/r>m/r){
printf("0\n");
continue;
}
ll ans=fac[n]*1ll*fphi[m]%r*q_pow(fac[m],r-2)%r;
printf("%lld\n",ans);
}
return 0;
}
4、组合数学
\[\mathrm{C}_n^m=\dfrac{n!}{m!(n-m)!}\quad \mathrm{A}_n^m=\dfrac{n!}{(n-m)!} \]卢卡斯定理
用于解决 \(n,m\) 值较大的排列数，要求模数 \(p\in\mathbb{P}\)。
\[\dbinom{n}{m}\bmod p=\dbinom{\lfloor n/p\rfloor}{\lfloor m/p\rfloor}\cdot\dbinom{n\bmod p}{m\bmod p}\bmod p \]inline ll C(int n,int m){
if(n==0||m==0) return 1;
return fac[n]*finv[m]%mod*finv[n-m]%mod;
}
inline ll Lucas(int n,int m){
if(!m) return 1;
return (C(n%mod,m%mod)*Lucas(n/mod,m/mod))%mod;
}
二项式定理
二项式系数公式：
\[(a+b)^n=\sum_{i=0}^n \dbinom{n}{i} a^{n-i}b^i \]错排列
已知 \(n\) 个数和 \(n\) 个位置，且每个位置有一个对应的数不能放置，求排列情况数。
设情况数 \(D_n\)，模拟可以得到 \(D_1=0,D_2=1,D_3=2\)，考虑递推，要求一步得到符合条件的序列，若前 \(n-1\) 个是符合条件的序列，则与任意一个调换都可以，因此有 \((n-1)\times D_{n-1}\)；若有且只有一个不符合条件，则需要和它调换位置，因此有 \((n-1)\times D_{n-2}\)，综上递推式为 \(D_n=(n-1)\times (D_{n-1}+D_{n-2})\)。
圆排列
在排列的基础上连链成环，因此要除去重复情况（显然 \({1,3,2,4}\) 与\({2,3,4,1}\) 在首尾相连成环以后是等价的），因此圆排列数为 \(A_n^n/n=(n-1)!\)
例题
1.P3233 HNOI2012 排队
分情况讨论，是由男生隔开老师还是由女生隔开老师。
若男生隔开老师，则情况数为男生人数的排列 \(n!\) 乘上老师可选择的位置与个数的排列 \(A_{n+1}^2\)，以及女生可选择位置与个数的排列 \(A_{n+3}^m\)。
若女生隔开老师，则情况数为男生人数的排列 \(n!\) 乘上隔开老师的女生情况数 \(m\) 乘上老师个数的排列 \(2!\)，将老师与这名女生视为一个整体，可以放的位置有 \((n+1)\) 个，剩余女生可选择的位置与个数的排列为 \(A_{n+2}^{m-1}\)。
所以最终结果：
\[\operatorname{ans}=n!\times \mathrm{A}_{n+1}^2\times \mathrm{A}_{n+3}^m + n!\times m\times 2\times (n+1)\times \mathrm{A}_{n+2}^{m-1} \]化简用大整数求解。
2.P4701 SDOI2016 排列计数
错排列的板子，答案为 \(\mathrm{C}_n^m\times D_{n-m}\)
int main(){
t=read();
fac[0]=finv[0]=1;
for(int i=1;i<=maxn;i++){
fac[i]=fac[i-1]*i%mod;
}
finv[maxn]=q_pow(fac[maxn],mod-2);
for(int i=maxn-1;i>=1;i--){
finv[i]=finv[i+1]*(i+1)%mod;
}
d[1]=0,d[2]=1,d[3]=2;
for(int i=4;i<=maxn;i++){
d[i]=(i-1)*(d[i-2]+d[i-1])%mod;
}
while(t--){
n=read(),m=read();
if(m+1==n){
printf("0\n");
continue;
}
if(m==n){
printf("1\n");
continue;
}
if(m==0){
printf("%lld\n",d[n]);
continue;
}
ll ans=fac[n]%mod*finv[m]%mod*finv[n-m]%mod*d[n-m]%mod;
printf("%lld\n",ans);
}
return 0;
}
3.P2606 ZJOI2010 排列计数
发现满足条件的数列满足一个小根堆的性质，我们将数列扩大一倍，在原有的叶子节点下加上一层。
设以 \(i\) 为根的子堆有 \(dp(i)\) 种符合条件的情况，设 \(siz(i)\) 表示以 \(i\) 为根子堆的大小，那么就有：
\[dp(i)=dp(i\times 2)\times dp(i\times 2+1)\times \dbinom{siz(i)-1}{siz(i\times 2)} \]因为小根堆与二叉搜索树有明显的区别，因此只要是子树中的节点无论怎样排列都满足性质，我们不需要考虑它具体的分布。
inline ll q_pow(ll x,ll p){
ll ans=1;
while(p){
if(p&1){
ans=ans*x%mod;
}
x=x*x%mod;
p>>=1;
}
return ans;
}
inline ll C(ll n,ll m){
return fac[n]*q_pow(fac[m]*fac[n-m]%mod,mod-2)%mod;
}
inline ll Lucas(ll n,ll m){
if(!m) return 1;
return (C(n%mod,m%mod)*Lucas(n/mod,m/mod))%mod;
}
int main(){
n=read(),mod=read();
fac[0]=1;
for(int i=1;i<=n;i++){
fac[i]=fac[i-1]*i%mod;
siz[i]=1;
}
for(int i=n;i>=1;i--){
siz[i>>1]+=siz[i];
}
for(int i=n+1;i<=2*n+1;i++){
dp[i]=1;
}
for(int i=n;i>=1;i--){
dp[i]=dp[i<<1]*dp[i<<1|1]%mod*Lucas(siz[i]-1,siz[i<<1])%mod;
}
printf("%lld\n",dp[1]);
return 0;
}
4.集合选数
题意：设一个集合 \(S\) 有 \(n\) 个元素，可以得到它有 \(2^n\) 个子集（包括 \(\varnothing\)），现在从中选出任意多个子集（至少一个），使他们交集的元素为 \(k\) ,求方案数，答案对 \(1e9+7\) 取模。
首先选 \(k\) 个数的方案数为 \(\mathrm{C}_n^k\)，而答案就应该为选数的方案数乘 \(n-k\) 个元素没有交集的方案数，显然会用到到容斥原理。
设 \(f(i)\) 表示交集大于等于 \(i\) 的元素个数，根据上面可以得到方案数也应为 \(\mathrm{C}_n^i\) 乘上一个选择方案数，那么现在考虑这个选择方案数。
可以得到剩下的 \(n-i\) 个数组成的集合 \(S'\) 的子集个数为 \(2^{n-i}\)，那么我们选择这些子集的方案数为 \(2^{2^{n-i}}\)，注意到我们是不能选空集的，方案数应当为 \(2^{2^{n-i}}-1\)，综上 \(f(i)=\mathrm{C}_n^i\times 2^{2^{n-1}}\)。
接下来考虑容斥，应当为“一加一减”，注意 \(f(i)\) 中 \(i\in[0,n-k]\)，所以是加上偶数项减去奇数项。
实现的话，费马小定理线性求逆元以及欧拉定理降幂即可，注意 \(f(i)\) 中的 \(n\) 实际是排除已经选择的 \(k\) 个元素后的 \(n-k\)。
ll n,k;
ll fac[maxn],inv[maxn];
ll f[maxn];
inline ll q_pow(ll x,ll p,ll md){
ll ans=1;
x=x%md;
while(p){
if(p&1){
ans=ans*x%md;
}
x=x*x%md;
p>>=1;
}
return ans;
}
inline ll C(ll n,ll m){
return fac[n]*inv[m]%mod*inv[n-m]%mod;
}
inline ll Lucas(ll n,ll m){
if(!m) return 1;
return C(n/mod,m/mod)*Lucas(n%mod,m%mod)%mod;
}
ll ans,sum;
int main(){
n=read(),k=read();
fac[0]=inv[0]=1;
for(int i=1;i<=n;i++){
fac[i]=fac[i-1]*i%mod;
}
inv[n]=q_pow(fac[n],mod-2,mod);
for(int i=n-1;i>=1;i--){
inv[i]=inv[i+1]*(i+1)%mod;
}
for(int i=0;i<=n-k;i++){
ll power=q_pow(2,n-k-i,mod-1);
f[i]=(C(n-k,i)*(q_pow(2,power,mod)-1))%mod;
if(i&1) f[i]=-f[i];
sum=(sum+f[i]+mod)%mod;
}
ans=(C(n,k)*sum)%mod;
ans=(ans+mod)%mod;
printf("%lld\n",ans);
return 0;
}
5.P6191 Bulls And Cows S
求有 \(i\) 头公牛时，除去母牛后剩下位置与公牛数的组合数，i的边界需要分类讨论。
假设 \(1\) 头公牛搭配 \(k\) 头奶牛，每组就是 \((k+1)\) 头牛，若 \((k+1)\mid n\)，公牛最多有 \(k/(n+1)\) 头，而如果 \((k+1)\nmid n\)，剩余的位置还可以放 \(1\) 头公牛，因此公牛最多有 \(k/(n+1)+1\) 头。
这里为了简便完成分类讨论，我们可以设置循环范围 \([0,(n-1)/(k+1)+1]\)，可以发现因为设置一定会 \(+1\)，所以当且仅当出现整除情况时与实际情况不符，因此把数量减小，不整除的情况数量不会改变，整除的情况就会与实际情况相符。（分块常用求块的个数）
\[\operatorname{ans}=\sum_{i=0}^{(n-1)/(k+1)+1}\ \dbinom{n-(i-1)\times k}{i} \]PS：考虑到即便整除，在不排满的情况下，也可以选择将 \((i-1)\times k\) 头奶牛夹在 \(i\) 头公牛中间，所以组合数应选择 \(n-(i-1)\times k\)。
int main(){
n=read(),k=read();
fac[0]=inv[0]=1;
for(int i=1;i<=n;i++){
fac[i]=fac[i-1]*i%mod;
}
inv[n]=q_pow(fac[n],mod-2);
for(int i=n-1;i>=1;i--){
inv[i]=inv[i+1]*(i+1)%mod;
}
for(int i=0;i<=(n-1)/(k+1)+1;i++){
ans=(ans+Lucas(n-(i-1)*k,i))%mod;
}
printf("%lld\n",ans);
return 0;
}
6.序列统计
题意：给定三个整数 \(n\),\(\ l\),\(\ r\) 求出长度从 \(1\) 到 \(n\) 的，元素均在 \([l,r]\) 内的单调不下降序列数，答案对 \(1e6+3\) 取模。
借本题详细的说一说此类问题。
首先来看一个化简的题面：求长度为 \(n\)，元素均在 \([1,m]\) 内的单调不下降序列个数。其实所为“单调不下降”是无需要考虑的，我们假设把这 \(m\) 个数字对应箱子从小到大依次排列，将 \(n\) 个球放入其中，所得到的的序列一定是“单调不下降”的。
那设第 \(i\) 个箱子放的球的个数为 \(x_i\)，可以得到一个线性不定方程：
\[n=\sum_{i=1}^{m} x_i \]这个方程的正整数解有多少种呢？用隔板法得到答案为 \(\binom{n-1}{m-1}\)，而因为有的数可以不选取，所以真正所求是方程的非负整数解，不妨设 \(y_i=x_i+1\)，方程可化为：
\[m+n=\sum_{i=1}^{m} y_i \]这里的 \(x_i\) 就可以作为 \(0\) 存在了，同理方案数为 \(\binom{m+n-1}{m-1}\)
这样一来，我们的最终答案为 \(\sum_{i=1}^{n} \binom{m+i-1}{m-1}\)
根据杨辉三角可以得到 \(\binom{x}{y}=\binom{x-1}{y}+\binom{x-1}{y-1}\)。
在原式基础上添加一个 \(\binom{m}{m}\)，就能得到：
\[\begin{aligned} \sum_{i=1}^n \dbinom{m+i-1}{m-1}+\dbinom{m}{m}&= \dbinom{m}{m}+\dbinom{m}{m-1}+\dbinom{m+1}{m-1}+\cdots+\dbinom{m+n-2}{m-1}+\dbinom{m+n-1}{m-1}\\&= \dbinom{m+1}{m}+\dbinom{m+1}{m-1}+\cdots+\dbinom{m+n-2}{m-1}+\dbinom{m+n-1}{m-1}\\&= \dbinom{m+n-1}{m}+\dbinom{m+n-1}{m-1}\\&= \dbinom{m+n}{m} \end{aligned}\]因此最终答案：
\[\operatorname{ans}=\dbinom{m+n}{m}-1 \]7.P2154 SDOI2009 虔诚的墓主人
部分分
暴力前缀和+组合数
\[\operatorname{ans}=\sum_{i=1}^n\sum_{j=1}^m \dbinom{u}{k}\cdot\dbinom{d}{k}\cdot\dbinom{l}{k}\cdot\dbinom{r}{k}\ (u,d,l,r\ge k) \]正解
首先因为地图很大而数量很少，所以考虑离散化，同时我们发现在同一横坐标的两棵树之间，其上下树的数量是不变的，只需要考虑左右树的数量，这个乘积需要支持查询和修改，考虑树状数组。
我们只需要预处理出组合数的值再判断并运算即可。
int n,m,w,k;
int x[maxn],y[maxn];
int fac[maxn][15],cnt[maxn][2];
pair p[maxn];
struct BIT{
int num[maxn];
inline int lowbit(int x){
return x&(-x);
}
inline void update(int x,int k){
k%=mod;
while(x<=w){
num[x]=(num[x]+k)%mod;
x+=lowbit(x);
}
}
inline int query(int x){
int res=0;
while(x){
res=(res+num[x])%mod;
x-=lowbit(x);
}
return res;
}
}tree;
int ans;
int cntw[maxn];
int main(){
n=read(),m=read(),w=read();
for(int i=1;i<=w;i++){
p[i].first=x[i]=read(),p[i].second=y[i]=read();
}
k=read();
fac[0][0]=1;
for(int i=1;i<=w;i++){
fac[i][0]=1;
for(int j=1;j<=k;j++){
fac[i][j]=(fac[i-1][j-1]+fac[i-1][j])%mod;//预处理出(0,0)到(w,k)的组合数
}
}
sort(x+1,x+w+1);
sort(y+1,y+w+1);
for(int i=1;i<=w;i++){
p[i].first=lower_bound(x+1,x+w+1,p[i].first)-x;
p[i].second=lower_bound(y+1,y+w+1,p[i].second)-y;
cnt[p[i].first][0]++,cnt[p[i].second][1]++;//统计前缀和并离散化
}
//横坐标为第一关键字，纵坐标为第二关键字排序
sort(p+1,p+w+1);
int cntn=0,cnts=1;
for(int i=1;i 8.P2467 SDOI2010 地精部落
话说这是孟周师兄暑假讲的一道题。他称其为“排列dp”，状态与排名有关。
首先可以发现，上升与下降实际上是等价的，因为将每个数用 \(n+1\) 减去后既不重复又恰好是与原数列相反的。
考虑dp，\(dp(i,j)\)表示前 \(i\) 个数中第 \(i\) 个数在该数列的排名（从小到大）为 \(j\)，我们只考虑开头为山峰的数列，那么分两种情况：
- 当 \(i\) 为奇数时，\(h_i>h_{i-1}\)，所以 \(dp(i,j)\) 可以由在前 \(i-1\) 个数中排名为 \([1,j-1]\) 的转移而来，注意到这里取 \(j\) 是不可行的，因为放入一个大于 \(h_{i-1}\) 的数，\(h_{i-1}\) 的排名不会改变。
- 当 \(i\) 为偶数时，\(h_i
初始状态为 \(dp(i,j)=1\)，接着得到递推式：
\[dp(i,j)=\begin{cases}\sum_{k=1}^{j-1} dp(i-1,k)&j=2\times z+1\\\\ \sum_{k=j}^{i-1} dp(i-1,k)&j=2\times z +2 \end{cases}\]其中 \(z\in \mathbb{N}\)（只是为了表示奇偶数）
发现单次转移是 \(O(n)\) 的，状态数是 \(O(n^2)\) 的，显然不过，考虑用前后缀和来优化（根据转移方程）
int main(){
n=read(),p=read();
dp[1][1]=1;
for(int i=2;i<=n;i++){
sum=0;
if(i&1){
dp[i][1]=0;
for(int j=2;j<=i;j++){
sum=(sum+dp[i-1][j-1])%p;
dp[i][j]=sum;
}
}
else{
dp[i][i]=0;
for(int j=i-1;j>=1;j--){
sum=(sum+dp[i-1][j])%p;
dp[i][j]=sum;
}
}
}
for(int i=1;i<=n;i++){
ans=(ans+dp[n][i])%p;
}
ans=(ans*2)%p;
printf("%d\n",ans);
return 0;
}
9.P3330 ZJOI2011 看电影
10.P3214 HNOI2011 卡农
先不考虑“排列”，“重复”的问题，设 \(dp(i)\) 表示 \(i\) 个数的合法方案数。
由题可得对一个合法的集合有如下三个要求：
- 集合中不含表示空集的元素
- 集合中不含重复的元素
- 集合中所有元素内的元素之和应都为偶数
若不考虑三者，可以发现，我们第 \(i\) 个集合可视为前 \(i-1\) 个的补充，也就是哪些元素目前还是奇数个，就在其中添加，因此有 \(f(i)=\mathrm{A}_{2^n-1}^{i-1}\)。
接着考虑满足条件与否，因为不能包含空集，也就是说前 \(i-1\) 个集合中必定包含奇数个的元素，否则是不合法的，而这个不合法的集合显然对于 \(i-1\) 个集合是满足要求的，所以这里的不合法情况数恰恰为 \(f(i-1)\)。接着若有重复元素又该如何去掉，考虑乘法原理，不合法情况数=与第 \(i\) 个集合重复集合的位置情况数 \(\times\) 重复集合的情况数 \(\times\) 去重后的情况数。重复位置的情况数自然是 \(i-1\)，因为目前只有一组集合重复，所以重复的集合种类要排除其余 \(i-2\) 个，也就是 \(2^n-1-(i-2)\)，最后剩余的合法情况不难得出 \(f(i-2)\)。
综上，\(f(i)=\mathrm{A}_{2^n-1}^{i-1}-f(i-1)+(i-1)\times (2^n-1-(i-2)\times f(i-2)\)，除此之外考虑排列后的去重，直接乘上 \(m!\) 在 \(\pmod{1e8+7}\) 定义下的逆元即可。
int main(){
n=read(),m=read();
power=(q_pow(2,n)-1+mod)%mod;
for(int i=1;i<=m;i++){
fac=fac*i%mod;
}
c[0]=1,dp[0]=1,dp[1]=0;
for(int i=1;i<=m;i++){
c[i]=c[i-1]*(power-i+1+mod)%mod;
}
for(int i=2;i<=m;i++){
dp[i]=c[i-1];
dp[i]=(dp[i]-dp[i-1]+mod-dp[i-2]*(i-1)%mod*(power-i+2)%mod+mod)%mod;
}
printf("%lld\n",dp[m]*q_pow(fac,mod-2)%mod);
return 0;
}
5、中国剩余定理相关
1.中国剩余定理
求解如下的线性同余方程：
\[\begin{cases} x\equiv a_1\pmod{n_1}\\ x\equiv a_2\pmod{n_2}\\ \vdots\\ x\equiv a_k\pmod{n_k} \end{cases}\](\(n_1,n_2\dots n_k\) 两两互质)
求解步骤如下：
- 求出所有模数的乘积 \(n=\prod_{i=1}^k n_i\)；
- 计算每一个方程的 \(m_i=\frac{n}{n_i}\)，并求出 \(m_i\) 在模 \(n_i\) 定义下的逆元 \(m_i^{-1}\)，并求出 \(c_i=m_im_i^{-1}\)（不取模）；
- 方程的唯一解 \(x=\sum_{i=1}^k a_ic_i \pmod{n}\)
证明如下：
把答案写成与方程系数相关的数应当为：
\[x=\sum_{i=1}^k a_i\times \frac{n}{n_i}\times inv(\frac{n}{n_i},n_i) \]当这个式子代入 第 \(i\) 个方程时，考虑到 \(\frac{n}{n_i}\times inv(\frac{n}{n_i},n_i)\) 在 \(\pmod{n_i}\) 下的值为 \(1\)，因此答案为 \(a_i\times 1=a_i\)；
而当带入到第 \(j \ (i\neq j)\) 个方程中，因为 \(\frac{n}{n_i} \bmod n_j\) 的值为 \(0\)，所以对答案是没有贡献的。
例题
1.P1495 曹冲养猪
套板子……
inline ll exgcd(ll a,ll b,ll &x,ll &y){
if(!b){
x=1,y=0;
return a;
}
ll d=exgcd(b,a%b,x,y);
ll tmp=x;
x=y,y=tmp-a/b*y;
return d;
}
ll p[15],a[15],c[15];
ll px=1;
ll ans;
int main(){
n=read();
for(int i=1;i<=n;i++){
p[i]=read();
a[i]=read();
px*=p[i];
}
for(int i=1;i<=n;i++){
ll x=0,y=0,tmp=px/p[i];
exgcd(tmp,p[i],x,y);
x=(x%p[i]+p[i])%p[i];
c[i]=tmp*x;
}
for(int i=1;i<=n;i++){
ans=(ans+a[i]*c[i])%px;
}
printf("%lld\n",ans);
return 0;
}
2.P3868 TJOI2009 猜数字
单拿出一个式子 \((n-a_i)\mid b_i\)，这个式子等价于 \(n-a_i \equiv0\pmod{b_i}\)，也就是 \(n\equiv a_i\pmod{b_i}\)，即整除式都等价于线性同余方程组，套板子即可。（需要龟速乘）
3.P2480 SDOI2010 古代猪文
求 \(g^{\sum_{d\mid n} \mathrm{C}_n^d} \bmod 9999116589\)。
这个式子可以化成：\(g^{\sum_{d\mid n}\mathrm{C}_n^d \bmod 999911658} \bmod 999911659\)
\(999911658\) 这个数很奇怪对吧，发现 \(999911658=2\times 3\times 4679\times 35617\)，所以我们列一个线性同余方程组：
\[\begin{cases} x\equiv \sum_{d\mid n} \mathrm{C}_n^d \pmod{2}\\ x\equiv \sum_{d\mid n} \mathrm{C}_n^d \pmod{3}\\ x\equiv \sum_{d\mid n} \mathrm{C}_n^d \pmod{4679}\\ x\equiv \sum_{d\mid n} \mathrm{C}_n^d \pmod{35617} \end{cases}\]然后求出 \(x\) 并快速幂求解。
//2 3 4679 35617
ll n,g;
ll pr[5]={0,2,3,4679,35617};
ll modx=999911658;
ll fac[40005];
ll a[5];
inline ll q_pow(ll x,ll p,ll md){
ll ans=1;
while(p){
if(p&1){
ans=ans*x%md;
}
x=x*x%md;
p>>=1;
}
return ans;
}
inline ll C(ll n,ll m,ll md){
if(n2.扩展中国剩余定理
求解如下的线性同余方程：
\[\begin{cases} x\equiv a_1\pmod{n_1}\\ x\equiv a_2\pmod{n_2}\\ \vdots\\ x\equiv a_k\pmod{n_k} \end{cases}\](\(n_1,n_2\dots n_k\) 不保证两两互质)
我们单单来看前两个式子，将 \(x\) 表示为 \(x=pn_1+a_1=qn_2+a_2\)，整理一下就有 \(pn_1-qn_2=a_2-a_1\)，首先可以得到当 \(\gcd(n_1,n_2)\mid (a_2-a_1)\) 时有解，而后解这个不定方程得出 \(x=pn_1+a_1\) 的值，随后代入到式子中，得出一个新的方程 \(x\equiv pn_1+a_1 \pmod{\operatorname{lcm}(n_1,n_2)}\)，以此类推不断合并。
例题
1.P4777 扩展中国剩余定理（EXCRT）
ll n;
ll a[maxn],p[maxn];
inline ll exgcd(ll a,ll b,ll &x,ll &y){
if(!b){
x=1,y=0;
return a;
}
ll ans=exgcd(b,a%b,x,y);
ll tmp=x;
x=y;
y=tmp-a/b*x;
return ans;
}
inline ll s_mul(ll x,ll p,ll md){
ll ans=0;
while(p){
if(p&1){
ans=(ans+x)%md;
}
x=(x+x)%md;
p>>=1;
}
return ans;
}
inline ll excrt(){
ll ans=a[1],px=p[1],x,y;
for(int i=2;i<=n;i++){
ll c=((a[i]-ans)%p[i]+p[i])%p[i];
ll g=exgcd(px,p[i],x,y),md=p[i]/g;
if(c%g) return -1;
x=s_mul(x,c/g,md);
ans+=x*px;
px*=md;
ans=(ans%px+px)%px;
}
return (ans%px+px)%px;
}
int main(){
n=read();
for(int i=1;i<=n;i++){
p[i]=read(),a[i]=read();
}
printf("%lld\n",excrt());
return 0;
}
3.扩展卢卡斯定理
放在这里的原因是，这东西需要的前置知识太多了……
首先，目的是要求一个 \(\binom{n}{m} \mod P\)，而 \(P\) 不一定为质数，简略的实现原理是把 \(P\) 唯一分解然后构造线性同余方程组，再用中国剩余定理合并并求出答案，下面是具体过程。
将 \(P\) 质因数分解使其等于 \(p_1^{\alpha_{1}}p_2^{\alpha_{2}}\cdots p_{k}^{\alpha_{k}}\)，之后可以求出：
\[\begin{cases} \mathrm{C}_{n}^{m} \bmod p_1^{\alpha_{1}}\\ \mathrm{C}_{n}^{m} \bmod p_2^{\alpha_{2}}\\ \vdots\\ \mathrm{C}_{n}^{m} \bmod p_k^{\alpha_{k}}\\ \end{cases}\]不难发现，这可以构造成线性同余方程组：
\[\begin{cases} \mathrm{C}_{n}^{m} \equiv a_1 \pmod{p_1^{\alpha_{1}}}\\ \mathrm{C}_{n}^{m} \equiv a_2 \pmod{p_2^{\alpha_{2}}}\\ \vdots\\ \mathrm{C}_{n}^{m} \equiv a_k \pmod{p_k^{\alpha_{k}}}\\ \end{cases}\]只需要求出 \(a_1,a_2\dots a_k\) ，就能通过中国剩余定理合并得到结果。（明显这 \(k\) 个模数是两两互质的）
之后我们单单拿出一个式子来看 \(\mathrm{C}_{n}^m \bmod p^\alpha\)，也就是求：
\[\dfrac{n!}{m!(n-m)!} \bmod p^\alpha \]这里由于 \(p^\alpha\) 并非对所有数都存在逆元，所以我们把式子写成：
\[\dfrac{\dfrac{n!}{p^{g(n)}}}{\dfrac{m!}{p^{g(m)}}\dfrac{(n-m)!}{p^{g(n-m)}}}\times p^{g(n)-g(m)-g(n-m)} \bmod p^\alpha \]其中 \(g(n)\) 表示 \(n\) 中所含 \(p\) 质因子的个数。
不难发现，分子和分母的三个分式都是等价的，只需研究：
\[\dfrac{n!}{p^{g(n)}} \bmod p^\alpha \]而后因为 \(n!\) 中必然含 \(\lfloor \frac{n}{p}\rfloor\) 个 \(p\) 质因子，所以可以变形成如下形式：
\[\begin{aligned} n!&=1\times 2\times \cdots \times n\\ &=(p\times 2p\times \cdots \times \lfloor \frac{n}{p}\rfloor\times p)\times (1\times 2\times \cdots \times n)\\ &=p^{\lfloor \frac{n}{p}\rfloor}\times (1\times 2\times \cdots\times \lfloor \frac{n}{p}\rfloor)\times \left(\prod_{i=1,i\not\equiv0\pmod{p}}^n i\right)\\ &=p^{\lfloor \frac{n}{p}\rfloor}\times (\lfloor \tfrac{n}{p}\rfloor)!\times \left(\prod_{i=1,i\not\equiv0\pmod{p}}^n i\right) \end{aligned}\]显然后面的求积是可以每几个构成一组，实际是相同的 \(\lfloor \frac{n}{p^\alpha}\rfloor\) 组以及一组剩余数字，就有：
\[\begin{aligned} n!&=p^{\lfloor \frac{n}{p}\rfloor}\times(\lfloor \tfrac{n}{p}\rfloor)!\times \left(\prod_{i=1,i\not\equiv0\pmod{p}}^n i\right)\\ &=p^{\lfloor \frac{n}{p}\rfloor}\times (\lfloor \tfrac{n}{p}\rfloor)!\times \left(\prod_{i=1,i\not\equiv0\pmod{p}}^{p^\alpha} i\right)^{\tfrac{n}{p^\alpha}}\times \left(\prod_{i=p^\alpha\cdot\lfloor\tfrac{n}{p^\alpha}\rfloor,i\not\equiv0\pmod{p}}^n i\right) \end{aligned}\]我们令 $f(n)=\dfrac{n!}{p^{g(n)}} $，由于式子中的 \(p^{\lfloor \frac{n}{p}\rfloor}\)，是可以消掉的，因此把式子化简成：
\[f(n)=f(\lfloor\tfrac{n}{p}\rfloor)\times \left(\prod_{i=1,i\not\equiv0\pmod{p}}^{p^\alpha} i\right)^{\tfrac{n}{p^\alpha}}\times \left(\prod_{i=p^\alpha\cdot\lfloor\tfrac{n}{p^\alpha}\rfloor,i\not\equiv0\pmod{p}}^n i\right) \]接下来看前面的 \(g(n)\)，会影响到的也就是 \(p^{\lfloor \frac{n}{p}\rfloor}\times(\lfloor \tfrac{n}{p}\rfloor)!\) 的结果，因此需要递归求解，即：
\[g(n)=\lfloor\tfrac{n}{p}\rfloor +g(\lfloor\tfrac{n}{p}\rfloor) \]回带到最后一个式子中：
\[\operatorname{ans}=\dfrac{f(n)}{f(m)f(n-m)} \times p^{g(n)-g(m)-g(n-m)} \bmod p^\alpha \]求解就可以了。（大喘气）
例题
1.P4720 扩展卢卡斯定理/exLucas
inline ll q_pow(ll x,ll p,ll md){
ll ans=1;
x%=md;
while(p){
if(p&1){
ans=ans*x%md;
}
x=x*x%md;
p>>=1;
}
return ans;
}
inline ll exgcd(ll a,ll b,ll &x,ll &y){
if(!b){
x=1,y=0;
return a;
}
ll d=exgcd(b,a%b,x,y);
ll tmp=x;
x=y,y=tmp-a/b*y;
return d;
}
inline ll get_inv(ll a,ll md){
ll x,y;
exgcd(a,md,x,y);
return (x%md+md)%md;
}
inline ll f(ll n,ll md,ll px){
if (!n) return 1;
ll num1=1,num2=1;
for(ll i=1;i<=px;i++){
if(i%md){
num1=num1*i%px;
}
}
num1=q_pow(num1,n/px,px);
for(ll i=px*(n/px);i<=n;i++){
if(i%md){
num2=num2*(i%px)%px;
}
}
return f(n/md,md,px)*num1%px*num2%px;
}
inline ll g(ll n,ll md){
if(n1){
a[++tot]=C(n,m,tmp,tmp);
b[tot]=tmp;
}
ll ans=0;
for(ll i=1;i<=tot;i++){
ll tmpx=md/b[i],inv=get_inv(tmpx,b[i]);
ans=(ans+a[i]*tmpx%md*inv%md)%md;
}
return ans;
}
ll n,m,p;
int main(){
n=read(),m=read(),p=read();
printf("%lld\n",exLucas(n,m,p));
return 0;
}
2.P2183 礼物
式子是非常好推的，我们记 \(sum=\sum_{i=1}^m w(i)\)，那么分两种情况讨论
- 若 \(sum>n\) 则礼物一定是不够分的，输出“\(\text{Impossible}\)”；
- 否则我们看任意一组个数为 \(sum\) 的组，方案数为 \(\mathrm{C}_{sum}^{w(1)}\cdot \mathrm{C}_{sum-w(1)}^{w(2)}\cdots \mathrm{C}_{w(m)}^{w(m)}\)，而这样的组一共有 \(\mathrm{C}_{n}^{sum}\)个，答案乘上即可。
考虑到模数不一定未知数，需要扩展卢卡斯。
int main(){
md=read();
n=read(),m=read();
for(int i=1;i<=m;i++){
w[i]=read();
sum+=w[i];
}
if(sum>n){
printf("Impossible\n");
return 0;
}
ans=exLucas(n,sum,md);
for(int i=1;i<=m;i++){
ans=ans*exLucas(sum,w[i],md)%md;
sum-=w[i];
}
printf("%lld\n",ans);
return 0;
}
6、矩阵乘法相关
对于 \(n\times m\) 的矩阵 \(A\)，\(m\times p\) 的矩阵 \(B\)，相乘后就有：
\[C_{i,j}=\sum_{k=1}^m A_{i,k}\times B_{k,j} \]矩阵乘法常用于优化递推式，尤其是重复的累乘问题。
例题
1.P1962 斐波那契数列
如果直接递推 \(fib(i)=fib(i-1)+fib(i-2)\) 是绝对会超时的，我们将矩阵 \(\begin{bmatrix}fib(n-1)&fib(n-2)\end{bmatrix}\)，通过乘法转移到 \(\begin{bmatrix}fib(n)&fib(n-1)\end{bmatrix}\)，因为 \(f(n)=1\times f(n-1)+1\times f(n-2),f(n-1)=1\times f(n-1)+0\times f(n-2)\)，因此可以写为：
\[\begin{bmatrix}f(n-1)&f(n-2)\end{bmatrix}\times \begin{bmatrix}1&1\\1&0\end{bmatrix}=\begin{bmatrix}f(n)&f(n-1)\end{bmatrix} \]不难发现这个可以用快速幂优化，设这个矩阵为 \(f\)，最终结果就是 \(\begin{bmatrix}f(n-1)&f(n-2)\end{bmatrix}\times f^{n-2}\)。
struct mat{
ll num[3][3];
mat operator *(const mat &tmp)const{
mat res;
memset(res.num,0,sizeof(res.num));
for(int i=1;i<=2;i++){
for(int j=1;j<=2;j++){
for(int k=1;k<=2;k++){
res.num[i][j]=(res.num[i][j]+num[i][k]*tmp.num[k][j]%m)%m;
}
}
}
return res;
}
}f,ans;
inline mat q_pow(mat x,ll p){
mat ans=x;
p--;
while(p){
if(p&1){
ans=ans*x;
}
x=x*x;
p>>=1;
}
return ans;
}
int main(){
f.num[1][1]=1,f.num[1][2]=1,f.num[2][1]=1,f.num[2][2]=0;
ans.num[1][1]=ans.num[1][2]=1;
n=read();
if(n<=2){
printf("1\n");
}
else{
ans=ans*q_pow(f,n-2);
printf("%lld\n",ans.num[1][1]);
}
return 0;
}
2.P4159 SCOI2009 迷路
考虑如果矩阵是 \(0-1\) 的那么只需要求矩阵 \(f\) 的 \(t\) 次幂，又发现矩阵只含 \([0,9]\) 的元素，考虑拆点，快速幂即可。
3.佳佳的 Fibonacci
首先令 \(S(i)=\sum_{j=1}^i F(i)\)，\(P(i)=\sum_{j=1}^i S(i)\)，不难发现，\(P(i)=\sum_{j=1}^i (n-i+1)\times F(i)= (i+1)\times S(i)-T(i)\)，考虑用矩阵乘法去优化。
可以发现的是：
于是得到：
\(\begin{bmatrix}1&1&0&0\\1&0&0&0\\1&1&1&0\\1&1&1&1\end{bmatrix}\) 的矩阵，快速幂即可。
最终答案 \(T(n)=(n+1)\times S(n)-P(n)\)
ll n,m;
struct mat{
ll num[5][5];
mat operator *(const mat &tmp)const{
mat res;
memset(res.num,0,sizeof(res.num));
for(int i=1;i<=4;i++){
for(int j=1;j<=4;j++){
for(int k=1;k<=4;k++){
res.num[i][j]=(res.num[i][j]+num[i][k]*tmp.num[k][j]%m)%m;
}
}
}
return res;
}
}f,ans;
inline mat q_pow(mat x,ll p){
mat ans=x;
p--;
while(p){
if(p&1){
ans=ans*x;
}
x=x*x;
p>>=1;
}
return ans;
}
int main(){
n=read(),m=read();
if(n==1){
ll res=1%m;
printf("%lld\n",res);
}
else if(n==2){
ll res=3%m;
printf("%lld\n",res);
}
else{
f.num[1][1]=1,f.num[1][2]=1,f.num[1][3]=0,f.num[1][4]=0;
f.num[2][1]=1,f.num[2][2]=0,f.num[2][3]=0,f.num[2][4]=0;
f.num[3][1]=1,f.num[3][2]=1,f.num[3][3]=1,f.num[4][3]=0;
f.num[4][1]=1,f.num[4][2]=1,f.num[4][3]=1,f.num[4][4]=1;
ans.num[1][1]=1,ans.num[2][1]=0,ans.num[3][1]=1,ans.num[4][1]=1;
ans=q_pow(f,n-1)*ans;
printf("%lld\n",(((n+1)*ans.num[3][1]%m)-ans.num[4][1]+m)%m);
}
return 0;
}