简单数论
-1 写在前面
阅读本文您可能需要一些基本的数学常识，如质数、合数、约数等的定义等。
如果没有特殊声明，本文所讨论的范围都在整数以内。
本文同步发表在我的 Cnblogs。
-1.0 长文警告，本文目前大约 2.2w，请谨慎阅读
-1.1 更新日志
- 2022.03.19 完成数论基础和组合计数的更新
- 2022.03.20 完成 crt 和 BSGS 的更新
- 2022.03.21 完成简单容斥和卡特兰数的更新
- 2022.03.22 完成欧拉函数的更新
- 2022.03.24 完成莫比乌斯函数的更新
- 2022.03.26 完成二项式反演的更新并改正了一些错误
-1.2 可能涵盖的内容
- 目前已有的
0 数论基础（欧拉筛、乘法逆元和 exgcd、crt 和 excrt、BSGS 和 exBSGS）
1 组合计数
2 简单容斥
3 卡特兰数
4 欧拉函数
5 莫比乌斯反演
6 二项式反演
- 未来可能有的
斯特林数
贝尔数
多项式计数
-1.3 参考文献
- 基础数论
各个百度百科，在没出地方都有放链接，这里就不放了。
Nickel_Angel's nest 大佬的 [数学]求解一类方程的方法 和 Luogu 日报
- 组合计数和容斥
单凭回忆和之前代码瞎口胡
- 卡特兰数
- 欧拉函数
handsomecui 大佬的 博客
- 莫比乌斯反演
神仙 command_block 的博客 的 博客
- 二项式反演
Rabbit_lb's Blog 大佬的 博客
如有侵权问题，请联系我，我会马上声明或修改。
如文章中有问题，请大佬们私信指出，蒟蒻不胜感激！
0 数论基础
我们先学习一些简单的数论基础为后面的学习做铺垫.
0.1 欧拉筛
问题 0.1.1：线性时间内求出 \(n\) 以内所有质数。
首先根据质数的定义，我们可以很容易打出一个 \(O(n\log n)\) 的埃氏筛。
bool flg[N];
int pri[N],tot;
void get_pri(int n){
flg[1]=1;
for(int i=2;i<=n;i++){
for(int j=i;i*j<=n;j++)flg[i*j]=1;
}
for(int i=1;i<=n;i++)if(!flg[i])pri[++tot]=i;
}
然后我们考虑如何优化这个过程。
我们发现中间 flg 的赋值有很多重复，比如 \(12=2\times6=3\times4\)，就被计算了两次。
这样子浪费了很多时间，我们能否构造一种算法使其每个合数只被赋值一次呢？
我们可以把一个合数分解为它的最小质因子乘以另外一个数，
对于每一个 \(i\)，从小到大枚举质数，直到枚举的质数能够整除 \(i\)，然后就停止枚举。
这样的话我们发现每一个数最多只会被乘一次，就可以 \(O(n)\) 筛质数了，这就是欧拉筛。
同时，我们可以再欧拉筛的过程中预处理出一些积性函数的值，如后面的欧拉函数和莫比乌斯函数。
bool flg[N];
int pri[N],tot;
void get_pri(int n){
for(int i=2;i<=n;i++){
if(!flg[i])pri[++tot]=i;
for(int j=1;j<=tot&&i*pri[j]<=n;j++){
flg[i*pri[j]]=1;
if(i%pri[j]==0)break;//这里是break，否则复杂度还是一个log
}
}
}
然后你就可以顺利通过 P3383 【模板】线性筛素数。
我们也可以用欧拉筛优化一些过程，比如质因数分解。
我们可以在欧拉筛的过程记录下来每一个数的最小质因子，然后质因数分解的时候时间复杂度就只有 \(O(质因数个数)\) 了。
bool flg[N];
int pri[N],tot,num[N];
void get_pri(int n){
for(int i=2;i<=n;i++){
if(!flg[i])pri[++tot]=i,num[i]=i;
for(int j=1;j<=tot&&i*pri[j]<=n;j++){
flg[i*pri[j]]=1,num[i*pri[j]]=pri[j];
if(i%pri[j]==0)break;//这里是break，否则复杂度还是一个log
}
}
}
int p[N],e[N],cnt;
void fen(int x){
cnt=0;
while(x!=1){
p[++cnt]=num[x];
e[cnt]=0;
while(x%p[cnt]==0)x/=p[cnt],++e[cnt];
}
}
0.2 乘法逆元和 exgcd
我们都知道很多题目为了防止出现精度问题或爆 long long，都会要求我们对某个数 \(p\) 取模。
如果在运算中有除法操作，一般来说 \(p\) 都是质数，那么我们如何执行这个操作呢？
乘法逆元就是用来解决这个问题的。
要学习乘法逆元，我们需要初步了解群论，这里就不介绍基础群论知识了。
乘法逆元的定义（摘自百度百科）：
定义：乘法逆元，是指数学领域群 \(G\) 中任意一个元素 \(a\)，都在 \(G\) 中有唯一的逆元 \(a^{-1}\)，具有性质 \(a\times a^{-1}=a^{-1}\times a=e\)，其中 \(e\) 为该群的单位元。
我们讨论的域是在模素数 \(p\) 域下的。
因此，我们有 \(\frac ab\equiv\frac ab\times b\times b^{-1}\equiv ab^{-1}\)，然后，我们就可以用乘法逆元代替除法运算了。
我们一般采用费马小定理来求一个数的逆元。
引理 0.2.1（费马小定理）：如果 \(a\) 不是 \(p\) 的倍数，那么 \(a^{p-1}\equiv 1(\bmod p)\)。
证明可以采用剩余系的方法，这里不展开了，如果有兴趣可以去百度百科。
这其实是欧拉定理的一种特殊情况，我们在后面介绍欧拉函数性质的时候会介绍。
既然如此，我们就得到了 \(a\times a^{p-2}\equiv 1\)，所以我们可以通过快速幂计算 \(a^{p-2}\% p\) 来得到 \(a\) 的逆元。
接下来我们来学习扩展欧几里得定理。
问题 0.2.2：给定 \(a,b,c\)，求出一组整数解 \(ax+by=c\) 或判无解。
首先，我们先判一组方程是否有解。
引理 0.2.2（裴蜀定理）：若不定方程 \(ax+by=c\) 有整数解，当且仅当 \(\gcd(a,b)|c\)。
证明可以用数学归纳法，这里不过多赘述，若对此感兴趣，可移步百度百科。
然后，我们就可以用递归的方法求解。
令 \(d=\gcd(a,b)\)，则我们先求出 \(ax+by=d\) 的解，然后 \(x,y\) 同时乘以 \(\frac c{d}\) 就可以了。
如果 \(b=0\)，那么 \(a=d\)，其中一组特解是 \(x=1,y=0\)。
否则我们先求出 \(bx+(a\%b)y=d\) 的一组解 \(x',y'\)，然后我们推一下式子：
所以，我们让 \(x=y',y=x'-\lfloor\frac ab\rfloor y'\) 即可。
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 tmp=x;
x=y,y=tmp-a/b*y;
return d;
}
然后，根据这组特解，我们可以解决很多问题，如正整数解个数等。
我们就可以顺利通过 P5656 【模板】二元一次不定方程 (exgcd)。
exgcd 还可以求解线性同余方程的问题：
问题 0.2.3 给定 \(a,b,m\)，求解 \(ax\equiv b(\bmod m)\) 或判无解。
我们转化一下问题的式子：\(ax+my=b\)，其中 \(y\) 可以是任意整数，这样就满足了同余的条件了。
我们发现，当 \(b=1\)时，我们就可以用这种方法求出 \(a\) 的逆元。
int inv(int a,int p){
int x,y;
int d=exgcd(a,p,x,y);
if(d!=1)return -1;
return (x%p+p)%p;
}
这种求逆方法在 \(p\) 不为质数时也适用。
根据这个，我们可以推出一种 \(p\) 为质数是的线性递推求逆的方法。
首先，我们有 \(1^{-1}=1\)，
然后，对于一个数 \(i\)，我们假设 \(p=ki+r\)，那么就有：
然后我们就可以愉快递推了。
int inv[N];
inv[1]=1;
for(int i=2;i<=n;i++)inv[i]=inv[i%p]*(p-p/i)%p;
因为这个笔者比较懒问题较简单，所以就不放例题了，直接放练习。
0.3 crt 和 excrt
crt，又称中国剩余定理，解决如下的问题：
问题 0.3.1：\(\begin{cases}x\equiv a_1(\bmod m_1)\\x\equiv a_2(\bmod m_2)\\\dots\\x\equiv a_n(\bmod m_n)\end{cases}\)
求出以上同余方程的最小非负整数解（保证 \(m_i\) 两两互质）。
我们假设 \(M=\prod m_i,M_i=\frac{M}{m_i}，t_i\equiv M_i^{-1}(\bmod m_i)\)，
\(M\) 是所有模数的乘积，\(M_i\) 时所有模数除了当前模数的乘积，\(t_i\) 是在模 \(m_i\) 的意义下 \(M_i\) 的逆元。
我们可以发现能够构造出一组特解 \(x'=\sum a_iM_it_i\)，我们发现每方程的 \(a_iM_it_i\) 在模其他模数的意义下都是 \(0\)，所以全部加起来这个式子仍然成立。
所以，我们能够得到方程通解：\(x=x'+kM\)，其中 \(k\) 可以是任意整数，所以我们就能很轻松的求出最小非负整数解。
然后，结合上面 0.2 的求逆方法，我们就能轻松过掉 P1495 【模板】中国剩余定理(CRT)/曹冲养猪。
signed main(){
scanf("%lld",&n);
M=1;
for(int i=1;i<=n;i++){
scanf("%lld%lld",&m[i],&a[i]);
M*=m[i];
}
for(int i=1;i<=n;i++){
Mi[i]=M/m[i];
int x=0,y=0;
exgcd(Mi[i],m[i],x,y);
ans+=a[i]*Mi[i]*x;
}
}
接下来，我们考虑扩展到 \(m_i\) 不一定互质的情况，这个时候，我们发现 \(M_i\) 不一定存在逆元，所以中国剩余定理不再适用。
对于一组方程：
问题 0.3.2：\(\begin{cases}x\equiv a_1(\bmod m_1)\\x\equiv a_2(\bmod m_2)\\\dots\\x\equiv a_n(\bmod m_n)\end{cases}\)
我们可以考虑合并两个线性同余方程，直到最后 \(x\equiv A(\bmod M)\)，就能够得出答案。
对于两个同余方程：\(\begin{cases}x\equiv a_1(\bmod m_1)\\x\equiv a_2(\bmod m_2)\end{cases}\)，我们发现可以这样表示：
我们发现这个东西有点眼熟，乍一看，不就是前面的不定方程吗，直接上 exgcd，当然，根据裴蜀定理，也有可能无解！
然后假设我们得出了一个特解 \(x'\)，通解是什么，又如何转化成同余方程的形式继续和别的方程合并呢？
我们发现，对 \(x\) 增减若干个 \(t=\text{lcm}(m_1,m_2)\)，是不会影响两个同余方程的结果的，
所以，我们可以把 \(x\) 满足的要求写成这样： \(x\equiv x'(\bmod \text{lcm}(m_1,m_2))\)，然后我们就能继续合并了。
这就是扩展中国剩余定理（excrt）的大致思路了，至此，我们能够顺利通过P4777 【模板】扩展中国剩余定理（EXCRT）。
signed main(){
scanf("%lld",&n);
for(int i=1;i<=n;i++)scanf("%lld%lld",&s[i],&r[i]);
M=s[1],X=r[1];
for(int i=2;i<=n;i++){
int d=exgcd(M,s[i],x,y),c=(r[i]-X%s[i]+s[i])%s[i];
if(c%d!=0)return puts("No solution"),0;//并没卵用的特判
x=mul(x,c/d,s[i]);//注意可能爆long long，要写龟速乘
X+=x*M;
M=M/d*s[i];
X%=M;
}
printf("%lld\n",(X%M+M)%M);
}
几道练习：
0.4 BSGS 和 exBSGS
BSGS（Baby Step Giant Step，中文名叫大步小步法），用于解决高次同余方程的问题：
问题 0.4.1（P3846 [TJOI2007] 可爱的质数/【模板】BSGS）：给出 \(a,b,p\)，求最小的 \(x\) 满足 \(a^x\equiv b(\bmod p)\)，或判无解（满足 \(p\in prime\)）。
本题 \(p<2^{31}\)，可以用 BSGS。
整体思路大概就是一个分块和哈希表。要满足的条件不一定是 \(p\in prime\)，只要满足 \(\gcd(a,p)=1\) 即可。
设 \(t=\lfloor\sqrt p\rfloor,x=it-j,(0\le j<t)\)，那么，我们有
然后我们把所有 \(0\le j<t\) 的 \(a^jb\) 放到一个 hash 里，然后对于所有 \((a^t)^i\)，看看有没有符合条件的 \(a^jb\) 即可。
然后你就可以顺利通过这道题了。
因为这题我是用 exBSGS 写的，就不放代码了。
问题 0.4.2（P4195 【模板】扩展 BSGS/exBSGS）：在问题 0.4.1 的基础上，扩展到 \(\gcd(a,p)\not=1\) 的情况。
我们发现可能到某一步 \(a^k\) 就变成了 \(0\)，无法除过去。
我们设 \(d=\gcd(a,p),a=a'd,p=p'd\)，
那么原来的方程等价于 \(a^{x-1}a'\equiv \frac bd(\bmod p')\)，
容易发现，当 \(d∤b\) 的时候，方程无解。
这个时候，\(a'\) 是有逆元的，我们把它移项，然后方程变成：\(a^{x-1}\equiv \frac bda'^{-1}(\bmod p')\)，然后继续循环，直到满足 \(\gcd(a,p')=1\)位置，然后就可以用 BSGS了。
因为 exBSGS 代码复杂度不高，但范围比 BSGS 大，所以笔者一般都写 exBSGS。
struct HashTable{
const int mod=1e6+7;
int mp[1000700],hs[1000700];
int find(int x){
int t=x%mod;
for(;mp[t]!=x&&mp[t]!=-1;)t=(t+107)%mod;
return t;
}
void insert(int x,int y){
int f=find(x);
mp[f]=x,hs[f]=y;
}
bool check(int x){
int f=find(x);
return mp[f]==x;
}
int query(int x){
int f=find(x);
return hs[f];
}
void clear(){
memset(hs,-1,sizeof hs);
memset(mp,-1,sizeof mp);
}
}ha;
int exBSGS(int a,int b,int m){
if(b%m==1)return 0;
int tmp=1,d=gcd(a,m),k=0;
while(d!=1){
if(b%d!=0)return -1;
b/=d,m/=d,++k;
tmp=tmp*(a/d)%m;
if(tmp==b)return k;
d=gcd(a,m);
}
ha.clear();
int t=ceil(sqrt(m));
ha.insert(b,0);
for(int i=1;i<=t;i++)b=b*a%m,ha.insert(b,i);
a=qpow(a,t,m);
for(int i=1;i<=t;i++){
tmp=tmp*a%m;
int j=ha.query(tmp);
if(j!=-1)return i*t+k-j;
}
return -1;
}
BSGS 和 exBSGS 的用途还有很多，可以求模意义下的二次剩余等，这里不展开吸收。
几道练习：
1 组合计数
1.1 加法原理和乘法原理
我们先通过一个例子来简单介绍一下加法原理和乘法原理。
例 1.1.1：奶牛从广州到浙江有 \(3\) 种高铁路线，\(3\) 种动车路线，\(2\) 种飞机路线，求总共有多少种方式从广州到浙江。
很显然，总共有 \(3+3+2=8\) 种方式。
像这样，选项之间互相不干预，或者说是互相独立的，并列事件的选择（也就是我们常说的分类计算），我们应该在它们之间使用加号统计，我们把这称为加法原理。
例 1.1.2：奶牛从广州到浙江有 \(3\) 条路径，从浙江到北京有 \(2\) 条路径，求总共有多少种方式从广州到北京。
很显然，总共由 \(3\times 2=6\) 种方式。
像这样，选项之间互相独立的，不同时进行的选择（也就是我们常说的分步计算），我们应该用乘号统计，因为可以理解为每一个前一种选择就对应后一种选择，我们把这称为乘法原理。
现在，如果你已经掌握了使用加法原理和乘法原理来计算方案数的方法，那就尝试做做下面的简单练习把！
练 1.1.1：奶牛从广东到浙江有 \(3\) 条路径，从浙江到北京有 \(4\) 条路径；不经过浙江，直接从广州到达北京有 \(3\) 条路径。求总共有多少条路径从广州到北京。
答案：\(3\times 4+3=15\)，加乘原理的混合运用。
接下来，我们开始将问题升级。
1.2 排列组合
我们考虑经典照相问题：
例 1.2.1：有 \(n\) 个人，选 \(m\) 个人出来排成一排，从左到右依次站好，然后对着 \(m\) 个人拍照，求最后的照片有多少种不同的结果。
注意：每个人互不相同。
我们这样考虑，首先我们给 \(n\) 个人标号 \(1,2,\dots,n\)，\(m\) 个位置标号 \(1,2,\dots,m\)。
先选择一个人去 \(1\) 号位置，有 \(n\) 种选择，此时这个被选中的人已经选定了位置，不能参与接下来的选择，剩下 \(n-1\) 个人。
再选择一个人去 \(2\) 号位置，有 \(n-1\) 种选择，此时这个被选中的人已经选定了位置，不能参与接下来的选择，剩下 \(n-2\) 个人。
以此类推，最后 \(m\) 号位置，有 \(n-m+1\) 种选择。
我们发现这个过程是分步计算的，适用乘法原理，所以最后有 \(n\times (n-1)\times \dots \times (n-m+1)\) 种不同结果。
我们发现这个式子可以写简单一点，如果对于一个非负整数 \(n\)，记 \(n!=n\times(n-1)\times(n-2)\times\dots\times 1\)，（特殊地，\(0!=1\)）（阶乘，则这个答案可以写作 \(\frac{n!}{(n-m)!}\)。
我们习惯把这样的问题叫做排列，我们可以这样描述 1.2.1 中的问题：
例 1.2.1'：形式化地，从 \(n\) 个互不相同的元素中，有序地选出 \(m\) 个的方案数 \(A_n^m=\frac{n!}{(n-m)!}\)。
特殊地，如果 \(m=n\)，那么我们把这个问题称为全排列，答案就是 \(n!\)。
例 1.2.2：有 \(n\) 个人，要选出 \(m\) 个人参加比赛，请问有多少中不同的方案数。
这一题的区别在于选出来的人是无序的，也就是说，对于 1 2 3 和 1 3 2，这两种选法我们认为它们是相同的。
我们先假设它们不同，那么方案数就是 \(A_n^m\)，但是每一种原来的选法会对应有 \(A_m^m\) 种原本是相同但是被认为是不同的方案，所以我们计算出的答案是真实的答案的 \(A_m^m\) 倍，除掉就可以了。
我们把这样的问题乘做组合，我们可以这样描述 1.2.2 中的问题。
例 1.2.2'：形式化地，从 \(n\) 个互不相同的元素中，无序地选出 \(m\) 个的方案数 \(C_n^m=\frac{n!}{m!(n-m)!}\)。
当然，组合数也可以写成二项式系数的形式：
1.3 组合数的性质
我们已经初步认识了排列数与组合数，接下来我们来探究一下与组合数有关的性质。
性质 1.3.1：若 \(n<m\)，则 \(\dbinom{n}{m}=0\)。
很显然，这个可以从组合数的定义入手，如果 \(n<m\)，我们无法选出一个合法方案，方案数为 \(0\)。
性质 1.3.2：\(\forall n\ge 0,\dbinom n0=\dbinom nn=1\)。
显然，两者一个是全选，一个是全部选，用 1.2 的式子也可以得出同样的结果。
性质 1.3.3：\(\dbinom nm=\dbinom {n-1}m+\dbinom{n-1}{m-1}\)。
这个性质可以直接暴力拆开式子通分得证，但我们还可以从组合意义入手。
解答 1.3.3：我们从 \(n\) 个学生选出 \(m\) 个学生，相当于不选班长，从其它 \(n-1\) 个学生种选 \(m\) 个，也可以选班长，从剩下 \(n-1\) 个学生种选择 \(m-1\) 个。式子得证。
我们发现这个式子非常眼熟，乍一看，不就是杨辉三角嘛。
我们可以根据这个式子做很多事情，这些我们以后再谈，先继续发觉组合数的性质。
性质 1.3.4：\(\sum_{i=0}^n \dbinom ni=2^{n}\)。
当然暴力拆式子也是可证的，但这不就太无趣了嘛。
解答 1.3.4：左边的式子相当于从 \(n\) 个元素中任意分成两组的方案数，相当于我们对于每一个元素可以选择第一组，也可以选择第二组，根据乘法原理，方案数是 \(2^n\)。
性质 1.3.5：\(\dbinom nm=\dbinom n{n-m}\)。
这个也很好理解，就是我们选择 \(m\) 个取出，相当于选择 \(n-m\) 个保留。
性质 1.3.6：\(\sum_{i=0}^n (-1)^i\dbinom ni=[n=0]\)。
这个的组合意义就没这么好考虑了，所以我们从杨辉三角入手。
\(n=0\) 的时候显然等于 \(1\)，考虑证明 \(n\not=0\) 的情况该式的值为 \(0\)。
解答 1.3.6：我们发现 \(i\) 是奇数的时候式子刚好覆盖了第 \(n-1\) 行的杨辉三角，是偶数的时候也一样，所以两者恰好相等，可以自行画图理解，
也可以暴力拆式子。
性质 1.3.7（二项式定理）：\((x+1)^n=\sum_{i=0}^n \dbinom ni x^i\)。
对于每一项考虑，对于第 \(i\) 项，我们相当于从 \(i\) 个 \((x+1)\) 项里面选 \(x\)，其它选 \(1\) 全部加起来。
当然如果你会生成函数，可以用性质 1.3.4 入手。
设第 \(i\) 行的组合数对应的普通型生成函数 \(C_i(x)=\sum_{i=0}^n\dbinom ni x^i\)，则根据性质 1.3.4 \(C_i(x)=C_{i-1}(x)\times(x+1)\)，因此性质 1.3.7 得证。
看不懂没关系，以后就懂了。
引理 1.3.7'（广义二项式定理）：当 \(n\) 是正整数，\(x\in[-1,1]\) 时，\((1-x)^{-n}=\sum_{r=0}^{\infty}\dbinom nr x^{r}\)。
这个式子后面生成函数会常用，这里不细讲。
性质 1.3.8（Lucas 定理）：\(\dbinom nm \equiv \dbinom {n\%p}{m\%p}\dbinom {\lfloor\frac np\rfloor}{\lfloor\frac mp\rfloor}(\bmod p)\)，其中 \(p\) 是质数。
这个式子非常重要，如果不会证明背结论就行了。
再证明这个式子之前我们需要做一点准备工作，在此之中我们默认 \(\equiv\) 操作是在 \(\bmod p\) 意义下进行的。
性质 1.3.8.1：\(\forall i\in[1,p-1],\dbinom pi\equiv 0\)。
证明：\(\dbinom pi=\frac {p!}{i!(p-i)!}\)，因为 \(p\) 是质数，所以分母中肯定不含因子 \(p\)，所以这个数肯定是 \(p\) 的倍数，所以 \(\bmod p=0\)。
性质 1.3.8.2：\((x+1)^p\equiv 1+x^p\)。
我们把左侧用二项式定理展开，发现说了首尾两项其它全部系数都是 \(p\) 的倍数。
回归性质 1.3.8，
我们令 \(a=\lfloor\frac np\rfloor\)，\(b=\lfloor\frac mp \rfloor\)，\(l=n\%p\)，\(r=m\%p\)，
那么 \(n=ap+l,m=bp+r\)，也就是说，我们要证明 \(\dbinom nm\equiv \dbinom ab \dbinom lr\)。
我们继续考虑二项式定理：
然后我们观察 \(x^m\) 项系数。
左边：\(\dbinom nmx^m\)，
我们发现右边的左边（也就是 \((1+x^p)^a\)）只能得到指数是 \(p\) 的倍数的，而 \((1+x)^l\) 只能得到指数 \(\le p-1\) 的。
所以右边：\(\dbinom ab\dbinom lrx^{bp+r}=\dbinom ab\dbinom lrx^m\)。
定理得证。
1.4 多种方法求组合数
如何求组合数的方法有很多，因为组合数数字是呈指数级增长的，所以一半会模个什么数，不然就老老实实打高精度吧。
对于不同的数据范围，我们可以应用不同的方法。
为了方便，在这里我们讨论的问题都是求 \(\dbinom nm\%p\)。
范围 1.4.1：\(n,m\le 1000\)。
根据性质 1.3.4，打个杨辉三角的表即可做到 \(O(n^2)\) 预处理，\(O(1)\) 查询。
范围 1.4.2：\(n,m\le 10^6,p\in prime\)。
我们直接从组合数最初在 1.2 里的公式。
我们预处理出 \(n\) 的阶乘，阶乘的逆元，然后就可以做到 \(O(n)\) 预处理，\(O(1)\) 查询。
范围 1.4.3：\(n,m\le10^{18},p\le10^5,p\in prime\)。
可以考虑使用 Lucas 定理，对于 \(n<p\)，用 1.4.2 的方法或者直接暴力计算，否则用 Lucas 定理递归处理。
单次时间复杂度 \(O(\log_p n)\)。
至此，你可以用这种方法同过 P3807 【模板】卢卡斯定理/Lucas 定理。
int C(int n,int m,int p){
int x=1,y=1;
for(int i=1;i<=m;i++){//因为懒得打表，直接暴力计算
x=x*(n-i+1)%p;
y=y*i%p;
}
// printf("%lld %lld\n",x,y);
return x*qpow(y,p-2,p)%p;
}
int lucas(int n,int m,int p){
if(!m)return 1;
return C(n%p,m%p,p)*lucas(n/p,m/p,p)%p;
}
范围 1.4.4：\(n,m\le10^6,p\le10^5\)。
我们发现 Lucas 定理不适用了，所以可以上扩展 Lucas 定理。
吐槽：感觉两者算法上并没有多大的关系。
因为 \(p\) 可能不是质数，所以上述方法全部废了。（1.4.1 时间复杂度显然不能接受,1.4.2 有可能出现没有逆元的情况,1.4.3 Lucas 定理不适用）
我们考虑转化问题：
设 \(p=p_1^{e_1}p_2^{e_2}\dots p_k^{e_k},(p_i\in prime)\)
那么我们只需要分别求出 \(C_n^m\) 对 \(p_1^{e_1},p_2^{e_2},\dots,p_k^{e_k}\)，然后 excrt 就好了。
那么我们如何求 \(C_n^m\%p^e\) 呢？
也就是 \(\frac{n!}{m!(n-m)!}\% p^e\)。
我们发现还是无法直接求，所以我们可以考虑把 \(p\) 的因子提出来单独算，其它的就可以求逆元了。
也就是说，问题转化成我们需要知道 \(n!\) 种含 \(p\) 的因子数量，以及其它因子的乘积。
我们可以这样算：
对于某个阶乘，假设 \(\ \ kp\le n<(k+1)p\)
我们发现前面的部分已经剔除了 \(p\) 的因子，可以求逆元了。而后面的可以每一项提一个 \(p\) 出来，这样我们获得了 \(k\) 个因子 \(p\)，剩下的就是 \(k!\)，可以递归下去。
那么这个问题就被解决了。
至此，你就可以拿下 P4720 【模板】扩展卢卡斯定理/exLucas。
int fac(int n,int p,int k){
if(!n)return 1;
int ans=1;
for(int i=2;i<=k;i++)
if(i%p)ans=ans*i%k;
ans=qpow(ans,n/k,k);
for(int i=2;i<=n%k;i++)
if(i%p)ans=ans*i%k;
return ans*fac(n/p,p,k)%k;
}
int C(int n,int m,int p,int k){
if(n<m)return 0;
int a=fac(n,p,k),b=fac(m,p,k),c=fac(n-m,p,k);
int cnt=0;
for(int i=p;i<=n;i*=p)cnt+=n/i;
for(int i=p;i<=m;i*=p)cnt-=m/i;
for(int i=p;i<=n-m;i*=p)cnt-=(n-m)/i;
return a*inv(b,k)%k*inv(c,k)%k*qpow(p,cnt,k)%k;
}
int crt(int n,int m){
return n*(p/m)%p*inv(p/m,m)%p;
}
int exlucas(int n,int m,int p){
int t=p,ans=0;
for(int i=2;i*i<=t;i++){
if(t%i)continue;
int tmp=1;
while(t%i== 0){
tmp*=i;
t/=i;
}
ans=(ans+crt(C(n,m,i,tmp),tmp))%p;
}
if(t>1)ans=(ans+crt(C(n,m,t,t),t))%p;
return ans%p;
}
如果模数 \(p\) 固定的话，可以把求阶乘的部分提前预处理一下，这样会更快一些。
1.5 简单组合计数问题
组合计数是信息学数论中笔者认为比较困难的一个板块了。
我们先从一些经典问题入手吧。
问题 1.5.1（圆排列问题）：\(n\) 个人选 \(m\) 个人站成一圈的方案数。
我们发现这个一个环我们可以选择 \(n\) 个地方断开，也就是会在总排列中重复计算 \(n\) 次，所以答案为 \(\frac{A_n^m}n\)。
问题 1.5.2（插板法）：有 \(n\) 个完全相同的球，放进 \(m\) 个互不相同的盒子里，每个盒子至少放一个的方案数。
我们可以把 \(n\) 个球排成一排，然后会形成 \(n-1\) 个缝隙。
我们发现我们可以把 \(m-1\) 个隔板放进去，然后 \(m-1\) 个隔板就对应 \(m\) 个有序的箱子，相邻两个隔板中间的球放进一个箱子，每个箱子至少有一个球。（因为不能在一个缝隙放两个隔板）
发现两个事情是等价的，而放隔板的方案是 \(C_{n-1}^{m-1}\)，所以问题 1.5.2 的答案也是这个。
问题 1.5.2'（插板法）：有 \(n\) 个完全相同的球，放进 \(m\) 个互不相同的盒子里的方案数。
我们发现我们没有办法计算盒子不加限制的方案数，所以我们可以给每个盒子强制放一个球，分配完之后再拿走，这样就可以很好的处理这个限制了，答案就是 \(C_{n+m-1}^{m-1}\).
问题 1.5.2''（插板法）：有 \(n\) 个完全相同的球，放进 \(m\) 个互不相同的盒子里，每个盒子至少放 \(k\) 个的方案数。
同理，如果至少放 \(k\) 个，我们先每个箱子都拿出 \(k-1\) 个，分配完再放回去即可。（如果 \(k=0\)，就是问题 1.5.2'，每个箱子先放一个，分配完再拿走）所以答案为 \(C_{n-(k-1)m-1}^{m-1}\)。
问题 1.5.3（捆绑法）：\(n\) 个人拍照，排成一排坐着，但是有 \(m\) 对情侣 \((2m\le n)\)，每一对情侣都要求坐在一起拍照，求方案数。
我们把每一对情侣捆绑在一起，排完之后一对情侣 AB 坐在了一起，但是他们两个可以换座位，所以每一对情侣会让答案乘以一个 \(2\)（乘法原理），所以答案就是 \((n-m)!2^m\)。
问题 1.5.4 （容斥）：\(n\) 个人拍照，排成一排坐着，有一对死对头，他们不能坐在一起，求方案数。
我们考虑容斥，答案等于总方案数减去死对头坐在一块的方案数，即 \(n!-2(n-1)!=(n-2)\times(n-1)!\)。
还有很多与组合数有关的容斥问题，不过比较复杂，我们在后面会见到的。
问题 1.5.5（Oiclass 大包子的算式）：有一个长度为 \(n\) 的数字串，里面放 \(k\) 个加号，算式两遍不允许有加号，求所有合法表达式的值的和对 \(1000000009\) 取模的结果，\(0\le k<n\le 10^5\)。
组合数入手好题，多种做法，这里介绍笔者想到的做法。
如果 \(k=0\) 直接判掉，接下来考虑 \(k\not=0\) 的情况。
我们把问题转化一下，求每一段数字的值乘上它出现次数。
对于一个长度为 \(i\) 的子串，如果它靠边，那么它的出现次数为 \(C_{n-i-1}^{k-1}\)，因为这一段中间不能填加号，边上要填一个，其他地方随意；
同理，如果不靠边，那它两遍都要放加号，出现次数为 \(C_{n-i-2}^{k-2}\)。
对于靠边的，靠左边或右边可以 \(O(n)\) 算出来，不靠边的，要对一样长度的子串一起考虑。
也就是我们要快速算出 \(\sum_{l=1}^{n-i+1}num[l,l+i-1]\)，我们可以考虑增量，假设我们已经求出了 \(i-1\) 对应的 \(s'\)，那么 \(s\) 相当于先扔掉最后面的一个，然后每一个剩下的再后面加一位。
预处理出前缀和，然后随便统计一下就好了。
signed main(){
scanf("%lld%lld%s",&n,&k,s+1);
if(k==0){
int sum=0;
for(int i=1;i<=n;i++)sum=(sum*10+s[i]-'0')%p;
printf("%lld\n",sum);
return 0;
}
fac[0]=ifac[0]=t[0]=1;
for(int i=1;i<=n;i++)fac[i]=fac[i-1]*i%p;
ifac[n]=qpow(fac[n]);
for(int i=n-1;i>=1;i--)ifac[i]=ifac[i+1]*(i+1)%p;
for(int i=1;i<=n;i++)t[i]=(t[i-1]*10)%p;
for(int i=1;i<=n;i++)sum[i]=(sum[i-1]*10+s[i]-'0')%p;
for(int i=1;i<=n;i++)sum2[i]=(sum2[i-1]+s[i]-'0')%p;
int ans=0,S=0;
for(int i=1;i<n;i++){
S=(S-(sum[n-1]-sum[n-i]*t[i-1]%p+p)%p+p)%p*10%p+(sum2[n-1]-sum2[i]+p)%p;
ans=(ans+S*C(n-i-2,k-2))%p;
}
// printf("%lld\n",ans);
for(int i=1;i<n;i++){
ans=(ans+sum[i]*C(n-i-1,k-1))%p;
}
// printf("%lld\n",ans);
for(int i=n;i>1;i--){
ans=(ans+(sum[n]-sum[i-1]*t[n-i+1]%p+p)%p*C(i-2,k-1))%p;
}
printf("%lld\n",ans);
}
问题1.5.6（数位dp，P6669 [清华集训2016] 组合数问题），给出 \(n,m,k\)，求有多少 \(0\le i\le n,0\le j\le\min(i,m)\) 满足 \(C_i^j\) 是 \(k\) 的倍数，答案对 \(1000000009\) 取模。共 \(t\) 组数据，\(k\) 不变，题目保证 \(k\) 是质数，\(1\le n,m\le10^{18},1\le t,k\le 100\)。
其实就是 \(C_i^j\equiv 0(\bmod k)\) 的个数。
然后用 Lucas 定理拆开：\(C_n^m\equiv \prod C_{a_i}^{b_i}\)。
我们发现满足条件必须要这些组合数中必须要有某个 \(C_{a_i}^{b_i}=0\)，因为每个 \(a_i,b_i<k\)。
所以要有 \(a_i<b_i\)。
我们发现问题可以转化成两个 \(p\) 进制数 \(i,j\)，两个有上限，要求 \(i\ge j\)，某一位 \(a_l<b_l\)，然后跑数位 dp 就行了。
后面有机会讲 dp 的时候可能也会有这道题。
int n,m,t,k,f[N][2][2][2][2];
int an[N],am[N],t1,t2;
void add(int &x,int y){x=(x+y)%p;}
int dfs(int len,bool isok,bool tpn,bool tpm,bool tpij){
if(!len)return isok;
if(f[len][isok][tpn][tpm][tpij]!=-1)return f[len][isok][tpn][tpm][tpij];
int res=0;
int upn=tpn?an[len]:k-1,upm=tpm?am[len]:k-1;
for(int i=0;i<=upn;i++)
for(int j=min(upm,(tpij?i:upm));j>=0;j--)
add(res,dfs(len-1,isok|(i<j),tpn&(i==upn),tpm&(j==upm),tpij&(i==j)));
return f[len][isok][tpn][tpm][tpij]=res;
}
练习：
2 简单容斥
学习完基础数论和组合计数之后，我们就可以学习一些简单容斥了。
容斥不是一个固定的算法，而是一种思想，总量减去多余量，全部减去重复计算的，这就是容斥的主要思想。
相信很多同学在小学的时候都学过类似问题，这里举一个例子。
问题 2.1：学校有社团活动，参加
A社的同学有 \(a\) 人，参加B社的有 \(b\) 人，参加C社的有 \(c\) 人，同时参加AB两社的同学有 \(d\) 人，同时参加AC两社的有 \(e\) 人，同时参加BC两社的有 \(f\) 人，同时参加ABC三社的有\(g\) 人，求学校有多少人参加了此次社团活动。
我们可以用三个集合 \(A,B,C\) 来表示分别参加了三个社团的学生，给出来的条件相当于每个集合的大小，两两交集大小，三个交集大小，要求三个集合并集大小，如果对集合有一定的了解，就能够轻松解决这个问题。
如果不是很了解，那么也没关系，可以画图来看会比较直观。
我们发现，\(a\) 代表的是图中红色的圆，也就是 \(1,4,6,7\) 部分，\(b\) 代表的是 \(2,4,5,7\)，\(c\) 代表的是 \(3,5,6,7\)，\(d\) 代表的是 \(4,7\)，\(e\) 代表的是 \(5,7\)，\(f\) 代表的是 \(6,7\)，\(g\) 代表的是 \(7\)，而我们要求的就是 \(1,2,3,4,5,6,7\) 每个部分各算一次的和。
我们考虑怎么求和，我们发现，\(1,2,3\) 部分只有 \(a,b,c\) 包含了，所以这些事肯定要加起来的，剩下的部分中，我们发现 \(4,5,6\) 都被统计了 \(2\) 次，但是我们只要 \(1\) 次，除了 \(a,b,c\) 之外，只有 \(d,e,f\) 包含，所以只能通过减掉 \(d,e,f\) 使 \(4,5,6\) 部分的系数保持为 \(1\) ，这个时候我们发现 \(7\) 的系数为 \(0\)，那就再加上 \(g\) 即可算出答案。
从特殊到一般，如果有多个集合，我们发现被偶数个集合同时包含的部分要乘以系数 \(-1\)，被奇数个集合同时保安的部分要乘以系数 \(1\)，可以统一成 \((-1)^{i-1}\)，我们把这个数成为容斥系数。不同的条件下容斥系数可能不同，需要通过推理得出。
问题 2.2：在问题 2.1 的条件下，增加一个条件：总共有 \(n\) 个学生。求有多少个学生没参加此次社团活动。
这个问题比较简单，通过上面问题 2.1 的解答，我们已经算出了参加至少一项社团活动的学生数（记为 \(m\)），那么一个都没有得学生总数就是 \(n-m\)。这种用全部减去非法状态的思想也属于容斥的范畴，很多比较困难的题目也会用到这种正难则反的思路。
问题 2.3（P6076 [JSOI2015]染色问题）：有一个 \(n\times m\) 的棋盘，要把它染成 \(c\) 种颜色，每个格子可以选择染成其中一种，或者不染，要求每一行至少有一个格子染色，每一列至少有一个格子染色，每一种颜色至少被染一次。求合法方案数对 \(1000000007\) 取模的值。
这里我们发现有关至少的条件比较难考虑，我们考虑至多，如果至多有 \(n-i\) 行， \(m-j\) 列有格子被染色，至多有 \(c-k\) 种颜色被使用的方案数，相当于有 \(i\) 行，\(j\) 列，\(c\) 种颜色强制不用，我们先选出 \(n-i\) 行，\(m-j\) 列，\(c-k\) 种颜色，方案数用组合数统计，然后我们每个格子有 \(c-k+1\) 种方案（除了染某种颜色外，还可以选择不染），这样就统计出来了。但是，我们发现会有重复统计。比如空两行的，在空一行的时候会重复两次，所以也要乘上某个容斥系数。
经过推算（套路），我们会发现容斥系数是 \((-1)^{i+j+k}\) 次方，然后全部加起来就好了。
for(int i=0;i<N;i++){
C[i][0]=1;
for(int j=1;j<=i;j++)C[i][j]=(C[i-1][j]+C[i-1][j-1])%p;
}
for(int i=0;i<=n;i++){
for(int j=0;j<=m;j++){
for(int k=0;k<=c;k++){
ans=(ans+((i+j+k)%2?-1:1)*C[n][i]*C[m][j]%p*C[c][k]%p*qpow(c-k+1,(n-i)*(m-j))%p)%p;
ans=(ans+p)%p;
}
}
}
几道类似的练习：
ABC 242 F Black and White Rooks
3 卡特兰数
卡特兰数（英语：Catalan number），又称卡塔兰数、明安图数，是组合数学中一种常出现于各种计数问题中的数列。（摘自百度百科）
3.1 卡特兰数的通项公式
引理 3.1.1：\(\begin{cases}C_n=1(n=0)\\C_{n+1}=\sum_{i=0}^{n}C_iC_{n-i}(n\ge 1)\end{cases}\)
这个就是卡特兰数的递推公式。
根据这个递推公式我们可以做很多东西，我们先求一下它的通项公式。
我们考虑它的意义，它可以表示成这样一个问题：
问题 3.1.2：在一个二维平面中，你初始在 \((0,0)\)，每次可以使横坐标或纵坐标 \(+1\)，要求保证每一个经过的坐标 \((x,y)\) 都有 \(x\ge y\)，求走到 \((n,n)\) 的方案数。
我们发现这个的答案其实就是卡特兰数。
我们考虑一个经典容斥，我们先不考虑 \(x\ge y\) 的限制，那么方案数就是 \(\dbinom{2n}{n}\)，在 \(2n\) 步里选 \(n\) 步向上走，其他向右走。
我们发现对于一种方案，它合法当且仅当每一个数都在图中的红线右下方（包括线上），那么我们考虑对于一种不合法方案，第一个越界的点为 \((x,x+1)\)，然后把剩下的操作全部取反（走上变成走右，走右变成走上），也就是剩下部分关于绿色的线对称，然后我们发现最后就会走到 \((n-1,n+1)\)。我们发现一条不合法路径都能对应一条从 \((0,0)\) 走到 \((n-1,n+1)\) 的路径，而且每一条路径也能对应回一个不合法路径（因为我们能够找到第一个越界点）。所以这两个事情是双射的。我们可以通过统计从 \((0,0)\) 走到 \((n-1,n+1)\) 的方案数代替不合法方案数，我们发现这个方案数非常好统计，就是 \(\dbinom{2n}{n-1}\)，方法和刚刚一样。
所以，我们就得出了卡特兰数的通项公式：
然后，我们还可以对这个式子进行化简来加速运算。
3.2 卡特兰数的运用
我们发现卡特兰数有很多应用：
- \(n\) 对括号匹配的方案数。
- \(n\) 个节点构成的不同二叉树个数。
- \(n+1\) 个叶子（\(n\) 个非叶子节点）的满二叉树形态数。
- \(n\) 个数入栈出栈排列的总数。
- 凸 \(n+2\) 变形的不同三角形分割的方案数。
- \(n\) 层的阶梯切割为 \(n\) 个矩形的方案数。
类似的应用还有很多，如何证明的方法很简单，基本都可以转化成我们在 3.1 讨论的问题。
据此，我们能够过掉很多模板题。
基本都是模板了
3.3 卡特兰数的扩展
我们考虑对上面问题 3.1.2 的扩展：
问题 3.3.1：在一个二维平面中，你初始在 \((0,0)\)，每次可以使横坐标或纵坐标 \(+1\)，要求保证每一个经过的坐标 \((x,y)\) 都有 \(x\ge y\)，求走到 \((n,m)\) 的方案数。
因为 \(n<m\) 显然答案是 \(0\)，所以我们只讨论 \(n\ge m\) 的情况。
我们仍然可以用上述方法经典容斥，对于不合法的方案，找到第一个越界的点，剩下的部分翻转，然后会到达 \((m-1,n+1)\)，然后我们对这个部分统计方案数，也就是 \(\dbinom{n+m}{n+1}\)，而所有方案的和是 \(\dbinom{n+1}n\)，两个相减。
我们能够据此通过一堆模板题：
问题 3.3.2（SOJ 22 游走）：从 \((0,0)\) 开始随机游走，每次随机走上下左右四个方向，求路径中经过的点距离 \((0,0)\) 的曼哈顿距离的最大值的期望。
首先曼哈顿距离 \(|x|+|y|=\max(|x+y|,|x-y|)\)，每次走一步可以看做是两个分别等概率 \(+1\) 或者 \(-1\)，
我们可以计算出所有方案最大值得和，最后再除以总方案数 \(4^n\)，即可算出答案。
考虑容斥，对于每一个 \(d\)，统计 \(\le d\) 的答案，这样，它的答案就是对于一个随机变量 \(a\)，求出它 \(\le d\) 的概率。由于两者可以看做是独立的，根据乘法原理，直接想成就得到原问题答案。问题转化成求 \(a\in[-d,d]\) 的概率。
记 \(L=\lceil\frac{n-d}2\rceil,R=\lfloor\frac{n+d}2\rfloor\)，如果只考虑最后一个点落在 \([-d,d]\) 的话，我们可以通过枚举 \(+1\) 的操作数，利用组合数轻松求出答案，也就是 \(\sum_{i=L}^R\dbinom ni\)。
对于计算最后一个点在范围之内，但是过程中有超出范围的不合法方案，我们可以考虑卡特兰数的经典翻折容斥，我们把第一次经过 \(d+1\) 或 \(-d-1\) 的方案后的操作取反，那么最后的 \(|a|\) 会在 \([d+2,3d+2]\) 中取值，正负两边一样，所以 \(a\le d\) 的方案数的式子为：
预处理除 \(\dbinom ni\) 的前缀和，能在 \(O(\frac nd)\) 的时间内计算出单个答案，总时间复杂度 \(O(n\log n)\)。
#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=1e6+10;
int n,p,fac[N],ifac[N],s[N],ans;
int qpow(int x,int y=p-2){
int m=1;
for(;y;y>>=1,x=x*x%p)if(y&1)m=m*x%p;
return m;
}
int C(int n,int m){
if(m>n)return 0;
return fac[n]*ifac[m]%p*ifac[n-m]%p;
}
signed main(){
scanf("%lld%lld",&n,&p);
fac[0]=ifac[0]=1;
for(int i=1;i<=n;i++)fac[i]=fac[i-1]*i%p;
ifac[n]=qpow(fac[n]);
for(int i=n-1;i>=1;i--)ifac[i]=ifac[i+1]*(i+1)%p;
s[0]=1;
for(int i=1;i<=n;i++)s[i]=(s[i-1]+C(n,i))%p;
for(int d=1,lst=0;d<=n;d++){
int L=(n-d+1)/2,R=(n+d)/2;
int res=(s[R]-(L<=0?0ll:s[L-1])+p)%p;
for(int j=1;R>=d+1;j++){
L-=d+1,R-=d+1;
res=(res+(s[R]-(L<=0?0ll:s[L-1])+p)%p*((j&1)?p-2:2ll))%p;
}
ans=(ans+(res*res%p-lst*lst%p+p)%p*d%p)%p;
lst=res;
}
printf("%lld\n",ans*qpow(qpow(4,n))%p);
}
我们发现有很多题并不直接考卡特兰数，而是考一些关于卡特兰数翻折方法的运用，所以要对卡特兰数深入理解，并不能一言不合就默写。
另：卡特兰数的通项公式也可以用生成函数来求，不过内容会在下文有提到。
4 欧拉函数
4.1 欧拉函数与它的性质
欧拉函数的定义（摘自 百度百科）：
显然，\(\varphi(1)=1\)，\(\varphi(p)=p-1(p\in prime)\)。
若 \(m=m_1m_2\)，我们考虑如何通过 \(\varphi(m_1),\varphi(m_2)\) 表示出 \(\varphi(m)\) 的值。
引理 4.1.2：\(\varphi(m)=\varphi(m_1)\cdot \varphi(m_2)(\gcd(m_1,m_2)=1)\)。
我们可以利用容斥来证，左边表示的是 \(m\) 以内的正整数有多少个与 \(m_1\)，\(m_2\) 均互质，反过来就表示 \(m-\) 与 \(m_1,m_2\) 其中一个的 \(\gcd>1\) 的个数。右边就是 \((m_1-\gcd(m_1,x)>1的个数)(m_2-\gcd(m_2,x)>1的个数)\) 。
拆括号，发现等于 \(m-m_2(\gcd(m_1,x)>1)的个数-m_1(\gcd(m_2,x)>1)的个数+(\gcd(m_1,x)>1)的个数\cdot(\gcd(m_2,x)>1)的个数\)。
我们可以这样理解，\(m_2(\gcd(m_1,x)>1)的个数\)相当于 \(m\) 以内不与 \(m_1\) 互质的数的个数， 另一个同理，但是有可能·有的数同时是两者的倍数，被减去了两次，要加回去，所以佳慧去后面的一项。
所以最终结果就是 \(m-\) 与 \(m_1,m_2\) 其中一个的 \(\gcd>1\) 的个数。所以引理 4.1.2 得证。
引理 4.1.3：如果 \(m_1,m_2\) 有相同质因数，则 \(\varphi(m)=m_1\varphi(m_2)\)。
这个结论证明方法类似，留给读者思考。
具体可以看这里：知乎。
引理 4.1.4：\(\varphi(m)=m\prod_{p|m,p\in prime}(1-\frac 1p)\)。
这个结论就是我们准备欧拉筛求欧拉函数的原理。
我们可以用引理 4.1.2 和 4.1.3 来证明。
我们先考虑对 \(m\) 进行分解质因数，设 \(m=\prod_{i=1}^k p_i^{e_i}(p_i \in prime),m'=\prod_{i=1}^k p_i\)，
那么根据引理 4.1.2 可知，\(\varphi(m')=\prod_{i=1}^k \varphi(p_i)=\prod_{i=1}^k (p_i-1)=m'\prod_{i=1}^k (1-\frac 1{p_i})\)。
然后再把重复因子乘进来，根据引理 4.1.3, 可知，\(\varphi(m)=\varphi(m')\prod_{i=1}^k p_i^{e_i-1}=m\prod_{p|m,p\in prime}(1-\frac 1p)\)。
引理得证。
引理 4.1.5（欧拉定理）：如果有 \(\gcd(a,m)=1\)，则 \(a^{\varphi(m)}\equiv 1(\bmod m)\)。
当 \(p=m\) 是质数时，它就是我们之前说的费马小定理：\(a^{p-1}\equiv 1(\bmod p)\)。
证明的话要引入既约剩余系，也叫简化剩余系或缩系。
我们设模 \(m\) 意义下的一个既约剩余系 \(R\{r_1,r_2,\dots,r_{\varphi(m)}\}\)，如果 \(\gcd(a,m)=1\)，那么 \(a\) 也必然在 \(R\) 里。
所以 \(R'{ar_1,ar_2,\dots,ar_{\varphi(m)}}\) 也是 \(m\) 的既约剩余系。
所以欧拉定理得证。
根据欧拉定理，我们也能得到：\(a\cdot a^{\varphi(m)-1}\equiv 1(\bmod m)\)，所以 \(a^{-1}\equiv a^{\varphi(m)-1}\)，所以我们就可以利用欧拉定理快速求一个数的逆元。
引理 4.1.6：\(m=\sum_{d|m}\varphi(d)\)。
我们考虑这个式子的实际意义，对于一个数对 \((x,y),(x\le y,y|m)\)，如果 \(\gcd(x,y)=1\)，那么它会被统计在 \(\sum_{d|m} \varphi(d)\) 里，而且我们假设 \(ky=m\)，则 \(\frac xy=\frac {kx}m\)，恰好能写成一个 \(m\) 为分母的真分数。我们发现样写刚好互不重复，并且恰好有 \(m\) 个，所以引理得证。
4.2 求欧拉函数
问题 4.2.1：给出一个数，求它的欧拉函数。
如果只求一个数的欧拉函数，我们可以直接分解质因数，然后用引理 4.1.4 直接求即可。
int phi(int x){
int ph=x;
for(int i=2;i*i<=x;i++){
if(x%i)continue;
ph=ph/i*(i-1);
while(x%i==0)x/=i;
}
if(x!=1)ph=ph/x*(x-1);
return ph;
}
问题 4.2.2：在线性时间内求出 \(1\) 到 \(n\) 所有正整数的欧拉函数。
根据引理 4.1.2 和引理 4.1.3 我们可以知道欧拉函数是一个积性函数。
我们可以利用欧拉筛来求欧拉函数。
对于欧拉筛，如果是 \(1\) 或者质数我们可以直接求出。
对于筛的时候遇到的一个数，如果，我们知道 \(i\cdot pri_j\) 其中一个质因子是 \(pri_j\)，并且此时我们已经知道了 \(\varphi(i)\) 的值。
我们可以知道，如果 \(i\) 中含质因子 \(pri_j\)，那么 \(\varphi(i\cdot pri_j)=pri_j\cdot \varphi(i)\)，否则 \(\varphi(i\cdot pri_j)=(pri_j-1)\cdot \varphi(i)\)，
然后我们就可以愉快地筛了。
bool flg[N];
int pri[N],phi[N],tot;
void get_pri(int n){
phi[1]=1;
for(int i=2;i<=n;i++){
if(!flg[i])pri[++tot]=i,phi[i]=i-1;
for(int j=1;j<=tot&&i*pri[j]<=n;j++){
flg[i*pri[j]]=1;
if(i%pri[j]==0){
phi[i*pri[j]]=phi[i]*pri[j];
break;
}
phi[i*pri[j]]=phi[i]*(pri[j]-1);
}
}
}
4.3 欧兰函数的运用
问题 4.3.1（UVA10820 交表 Send a Table）（P2158 [SDOI2008] 仪仗队）。
两道题其实都是一样的，发现都是求 \(\sum_{i=1^n}\sum_{j=1}^i [\gcd(i,j)=1]\)，然后再针对题目做一些处理。
可以用莫比乌斯反演去做（后面会讲），但是我们发现对欧拉函数做个前缀和即可。
问题 4.3.2（UVA11426 拿行李（极限版） GCD - Extreme (II)）：给定 \(n\)，求 \(\sum_{i=1}^{n-1}\sum_{j=i+1}^{n}\gcd(i,j)\)。
本题有好几倍经验，题解里面有，这里不列举了。
同样可以莫反，但是 \(i,j\) 上限都是 \(n\) 的性质让我们也可以用欧拉函数。
我们发现我们要求的答案 \(Ans=\frac{\sum_{i=1}^n\sum_{j=1}^n\gcd(i,j)-\sum_{i=1}^n\gcd(i,i)}2\)。
然后我们然后分子两边，右边求个和即可，左边枚举 \(\gcd\)，然后就能得到 \(\sum_{d=1}^nd\sum_{i=1}^n\sum_{j=1}^n[\gcd(i,j)=d]\)，
然后把 \(i,j\) 同时提出一个 \(d\) 出来，就能得到 \(\sum_{d=1}^nd\sum_{i=1}^{\lfloor\frac nd\rfloor}\sum_{j=1}^{\lfloor\frac nd\rfloor}[\gcd(i,j)=1]\)。
发现和问题 4.3.1 一样了，可以做个前缀和之后 \(O(n)\)，也可以再做一个前缀和整除分块达到 \(O(\sqrt n)\)。
4.4 扩展欧拉定理及其运用
引理 4.4（扩展欧拉定理）：如果 \(a,m\) 不一定互质，有：
\[a^c\equiv\begin{cases} a^{c\ \bmod \phi(m)}&\gcd(a,m)=1\\ a^c&\gcd(a,m)\not=1\and c<\phi(m)\\ a^{c\ \bmod \phi(m)+\phi(m)}&\gcd(a,m)\not=1\and c\ge \phi(m) \end{cases} \]
当 \(m=1\) 时显然成立，下面讨论 \(m\not=1\) 的情况。
- 若 \(\gcd(a,m)=1\)，根据欧拉定理，有 \(a^{\varphi(m)}\equiv 1\)，所以我们可以除掉若干个 \(a^{\varphi(m)}\)，所以 \(a^c\equiv a^{c\ \bmod \varphi(m)}\)。
- 若 \(\gcd(a,m)\not=1\and c\ge\varphi(m)\)，我们继续讨论，我们先证明 \(a\in prime\) 的情况。
引理 4.4.1：\(\forall p\in prime,r\in Z^+,\varphi(p^r)=(p-1)\cdot p^{r-1}\)。
这个用引理 4.1.4 带进去就能得到。
引理 4.4.2：\(\forall k\in Z,\exist a,b,x,y\in Z^+\)，对于 \(x^a\cdot y^b=k\)，都有 \(a,b\le\varphi(k)\)。
- 若 \(k=p^r(p\in prime)\)，那么根据引理 \(1\)，\(\varphi(k)=(p-1)\cdot p^{r-1}\)，问题等同于证明 \((p-1)\cdot p^{r-1}\ge r\)。
当 \(p=2\) 时，手算 \(r=1,2,3\) 时成立，\(r>3\) 时用归纳法，显然也成立。
当 \(p>2\) 时，我们发现不等式右边不变的情况下，左边变大，所以肯定也成立。
- 若 \(k\) 时多个质数幂的成绩，我们发现不等式左边相当与全部乘起来，右边取 \(\max\)，显然也成立。
引理 4.4.3：\(\exist r\le c,a^{\varphi(m)+r}\equiv a^r(\bmod m)\)。
我们对于一个 \(m\)，令 \(m=t\times a^r\)，并且使 \(\gcd(a,t)=1\)，那么根据引理 4.1.4，显然有 \(\varphi(t)|\varphi(m)\)，
根据欧拉定理，有 \(a^{\varphi(t)}\equiv1(\bmod t)\)，同余式同时乘以 \(a^r\)，我们就得到 \(a^{\varphi(m)+r}\equiv a^r(\bmod m)\)。
因此，我们就有：\(a^c\equiv a^{c-r+r}\equiv a^{c-r+\varphi(m)+r}\equiv a^{c+\varphi(m)}(\bmod m)\)。
用归纳法，我们可知：\(a^c\equiv a^{c+k\varphi(m)}\)，要求 \(k\) 是整数，而且指数要是正数。
我们发现最小合法的指数：\(c\% \varphi(m)+\varphi(m)\)，所以 \(a\) 是质数时得证。
这时候工作已经差不多了。当 \(a\) 是一个质数的幂次时（假设是 \(p^k\)），我们把上述过程提一个 \(k\) 出来，完了之后再放回去就好了：\(a^c\equiv p^{ck}\equiv p^{ck+k\varphi(m)}\equiv(p^k)^{c+\varphi(m)}(\bmod m)\)。
当 \(a\) 是多个质数的幂次时，假设 \(a=\prod p_i^{e_i}\)，那么 \(a^c\equiv\prod (p_i^{e_i})^c\equiv\prod (p_i^{e_i})^{c+\varphi(m)}\equiv(\prod p_i^{e_i})^{c+\varphi(m)}\equiv a^{c+\varphi(m)}(\bmod m)\)。
至此，扩展欧拉定理得证。
扩展欧拉定理很好用，直接上板子即可，这里列举两道：
问题 4.4.4（P4139 上帝与集合的正确用法）：求 \(2^{2^{2^{\dots}}}(\bmod p)\)。
列个方程 \(a\equiv 2^a\)，然后发现这个指数很大，肯定大于 \(\varphi(p)\)，直接上公式即可。
问题 4.4.5（CF17D Notepad）：求 \((a-1)a^{n-1}(\bmod p)\)。
\(n\) 特别大，边读入边取模就行，注意特判 \(n-1<\varphi(p)\) 的情况。
5 莫比乌斯反演
5.1 狄利克雷卷积与数论函数
ps：这一部分看不懂并不会影响莫比乌斯反演的使用，但是还是建议尽量看一下。
定义 5.1.1：数论函数，就是值域为整数的函数。
因此，我们知道刚刚研究的欧拉函数 \(\varphi\) 就是一个数论函数。
定义 5.1.2：定义两个数论函数 \(f(n),g(n)\) 的狄利克雷卷积：
\((f\times g)(n)=\sum_{d|n}f(d)g(\frac nd)\)。
我们发现狄利克雷卷积满足交换律和结合律（显然）。
接下来，我们来认识几个简单的数论函数。
\(I(n)=1\)，就是无论 \(n\) 取啥值函数值都为 \(1\)，仅仅只是为了写成函数形式，叫做恒等函数。
\(\epsilon(n)\)，也写作 \(e(n)\)，\(\epsilon(n)=\begin{cases}1&n=1\\0&\text{else}\end{cases}\)，叫做元函数，因为它是卷积的单位元。
\(id(n)=n\)，叫做单位函数， \(id^x(n)=n^x\) 叫做幂函数。\(id(n)\) 是 \(id^x(n)\) 的特殊情况。
我们在前面讲过很多次积性函数，但是并没有给出严谨的定义和证明（笔者的谢罪），现在我们来给出。
定义 5.1.3：积性函数：对于一个数论函数 \(f(n)\)，如果对于 \(\gcd(a,b)=1\) 时，都有 \(f(ab)=f(a)f(b)\)，那么该函数是积性函数。
完全积性函数：对于一个数论函数 \(f(n)\)，如果对于任意整数 \(a,b\) 都有 \(f(ab)=f(a)f(b)\)，那么该函数是完全积性函数。
我们发现 \(I(n),\epsilon(n),id(n)\) 都是完全积性函数，并且同时也是积性函数。前面讲到的 \(\varphi(n)\) 是积性函数，但并不是完全积性函数函数。
显然，完全积性函数属于积性函数。
积性函数有很多优美的性质，我们马上就会来研究。
性质 5.1.4：对于一个积性函数 \(f(n)\)，有 \(f(1)=1\)。
根据积性函数的定义，有 \(f(1)=f(1)f(1)\)，所以有 \(f(1)=1\)。（暂不考虑 \(f(1)=0\) 的情况）
性质 5.1.5：对于一个积性函数 \(f(n)\)，如果对于 \(x=\prod_{i=1}^kp_i^{e_i}(p_i\in prime)\)，那么 \(f(x)=\prod_{i=1}^k f(p_i^{e_i})\)。
这个也很好证明，\(p_i^{e_i}\) 两两互质，根据积性函数的性质，显然满足这个条件。
性质 5.1.6：两个积性函数的卷积也是积性函数。
假设两个积性函数 \(f_1(n),f_2(n)\)，它们的狄利克雷卷积 \(g(n)=\sum_{d|n}f_1(d)f_2(\frac nd)\)。
假设两个数 \(a,b\) 满足 \(\gcd(a,b)=1\)，那么：
所以 \(g(n)\) 也是一个积性函数。
接下来的性质和函数的逆有关，我们先来定义函数的逆。
定义 5.1.7：定义狄利克雷卷积意义下的逆：若两个数论函数 \(f,g\) 满足 \(e=f\times g\)，则称 \(f,g\) 互逆。
这个定义可类比 0.2 中的乘法逆元来理解。
性质 5.1.7：一个积性函数的逆也是积性函数。
假设一个积性函数 \(f(n)\)，设数论函数 \(g\) 满足 \(e=f\times g\)，我们需要证明对于 \(\gcd(a,b)=1\)，\(g(ab)=g(a)g(b)\)。
- 若 \(a,b\) 中有一个等于 \(1\) 时：
如果 \(a=b=1\)，则 \(e(1)=1=f(1)g(1)\)，由于积性，根据性质 5.1.4 得 \(f(1)=1\)，所以 \(g(1)=1\)，所以当 \(a,b\) 中有一个等于 \(1\) 的题时，结论显然成立。
- 当 \(ab>1\) 时：
考虑数学归纳法，假设我们已经证明了所有 \(a'<a,b'<b\) 时的结论成立。
引理得证。
5.2 莫比乌斯函数
聊了这么久的狄利克雷卷积，那么这和我们要讲的莫比乌斯函数有设么关系呢？
我们先假设两个数论函数：\(f,g\) 满足一下关系：\(f(n)=\sum_{d|n}g(d)\)，
那么我们也可以把这个式子写成这样：\(f(n)=\sum_{d|n}I(d)g(\frac nd)\)。
然后，我们就发现：\(f=I\times g\)，
我们可以很轻松的通过 \(g\) 在 \(O(n\log n)\) 的时间内求出 \(f\)，但是如果知道 \(f\) 如何快速求 \(g\) 呢？
我们可以在等式两边同时乘以 \(I^{-1}\) 得：\(g=I^{-1}\times f\)。
我们把 \(I^{-1}\) 这个函数取一个好听一点的名字，叫做莫比乌斯函数 \(\mu\)。
我们发现 \(\mu\) 是一个积性函数，根据我们之前的引理 5.1.7 可以得知。
我们来研究一下这个函数的性质。
问题 5.2.1：探究 \(\mu\) 的函数值。
- 当 \(n=1\) 时，显然 \(\mu(n)=1\)，
- 当 \(n \in prime\) 时，根据 \(\sum_{d|n}I(d)\mu(\frac nd)=\mu(n)+1=e(n)=0\)，可以得知 \(\mu(n)=-1\)，
- 当 \(n=p^k(k>1)(p\in prime)\) 时，我们发现 \(\sum_{d|n}I(d)\mu(\frac nd)=\sum_{i=0}^k\mu(p^i)=e(n)=0\)，我们发现，要保证这个式子始终为 \(0\)，除了 \(\mu(1)=1,\mu(p)=-1\) 外，\(\mu(p^k)=0\)。
- 当 \(n=\prod_{i=1}^k p_i^{e_i}(p_i\in prime\and e_i\ge 1)\) 时，根据 \(\mu\) 是积性函数，我们能够得知，\(\mu(n)=\prod_{i=1}^k \mu(p_i^{e_i})\)，如果这之中有一个的函数值为 \(0\)，那么 \(\mu(n)\) 就为 \(0\)，否则的话就看 \(k\) 的奇偶性。因为此时 \(e_i\ge 1\and e_i\le 1\)，所以这里的 \(e_i=1\)，所以 \(\mu(n)=\prod_{i=1}^k \mu(p_i)=(-1)^k\)。
至此，我们得出了莫比乌斯函数的定义式：
引理 5.2.2：\(\sum_{d|n}\mu(d)=[n=1]\)。
这个式子根据最初 \(\mu\) 的定义就能够知道了。
\(\sum_{d|n}\mu(d)=\sum_{d|n}\mu(d)I(\frac nd)=e(n)=[n=1]\)。
根据这个式子，我们就可以进行最基础的莫比乌斯反演了。
问题 5.2.3：如何线性求莫比乌斯函数的值。
在问题 0.1.1 的时候（不会真的有人去翻把？）就有说过，欧拉筛可以线性筛出积性函数的值，莫比乌斯函数也不例外。
在欧拉筛的过程中，如果 \(i\% pri_j=0\)，我们发现 \(i\cdot pri_j\) 种肯定含一个平方因子：\(pri_j^2\)，所以直接不处理（赋值为 \(0\)），否则的或相当于在 \(i\) 的基础上增加一个质数。
我们可以轻松地打出这个筛。
bool flg[N];
int pri[N],tot,mu[N];
void get_mu(){
mu[1]=1;
for(int i=1;i<=n;i++){
if(!flg[i])pri[++tot]=i,mu[i]=-1;//p是质数的情况
for(int j=1;j<=tot&&pri[j]*i<=n;j++){
flh[i*pri[j]]=1;
if(i%pri[j]==0)break;//肯定含平方因子
mu[i*pri[j]]=-mu[i]; //如果i中含平方因子，那么两个的函数值都为0，否则1变-1，-1变1
}
}
}
但是当我们使用莫比乌斯反演时，我们大概率需要整除分块的帮忙，否则单次查询时间复杂度还是 \(O(n)\)，并不是很优。
5.3 整除分块
问题 5.3（整除分块）：快速求 \(\sum_{i=1}^n \lfloor\frac ni\rfloor\)。
直接求 \(O(n)\) 显然不够优，我们考虑优化。
我们发现，这些 \(\lfloor\frac ni\rfloor\) 中，有很多值是相同的，我们考虑能不能同时统计。
我们在求这个式子之前，我们也需要做一些准备活动。（老套路了
引理 5.3.1：\(\lfloor \frac ni\rfloor\) 有 \(O(\sqrt n)\) 种取值，具体地，最多只有 \(2\sqrt n\) 种取值。
对于 \(i\le \sqrt n\)，有 \(\le \sqrt n\) 个数，所以取值也只有 \(\le \sqrt n\) 种；对于 \(i>\sqrt n\)，有 \(\frac ni<\sqrt n\)，也最多只有 \(\sqrt n\) 种，所以总共不会超过 \(2\sqrt n\)。
引理 5.3.2：对于一个 \(i\in [1,n]\)，所有的 \(j\in[i,\lfloor\frac n{\lfloor\frac ni\rfloor}\rfloor]\) 都有 \(\lfloor \frac ni\rfloor=\lfloor \frac nj\rfloor\)，且对于所有 \(j>\lfloor\frac n{\lfloor\frac ni\rfloor}\rfloor\)，都不满足 \(\lfloor\frac ni\rfloor=\lfloor\frac nj\rfloor\)。
换句话说，也就是 \(\lfloor\frac n{\lfloor\frac ni\rfloor}\rfloor\) 是 \(\lfloor\frac ni\rfloor\) 相等的上界。
我们假设 \(k=\lfloor\frac ni\rfloor,ki+r=n\)，假设 \(\lfloor\frac n{i+d}\rfloor=k\)，那么存在一个 \(0\le r'<i+d\) 满足 \(k(i+d)+r'=n\)。
我们得到：\(r'=r-kd\)，我们发现 \(d_{\max}=\lfloor\frac rk\rfloor\)。
引理得证。
至此，我们就可以用 \(O(\sqrt n)\) 的时间内求出问题 5.3 的答案了。
具体地，我们可以开两个指针 \(l,r\)，从低位往高位扫，每次统一统计，更新两个指针即可。
int sum=0;
for(int l=1,r;l<=n;l=r+1){
int x=n/l;r=n/x;
sum+=(r-l+1)*x;
}
这个结论在向上取整的时候同样适用，具体式子可以用类似方法推一下，整除分块也可以解决许多其他计数问题，这里不细讲了。
5.5 莫比乌斯函数的运用
我们给几道例题：
问题 5.5.1（问题 4.3.2’）：给定 \(n,m\)，求 \(\sum_{i=1}^n\sum_{j=1}^m\gcd(i,j)\)。
回到那一道无数倍经验的题，我们把 \(i,j\) 的上限改成不一样的，那应该怎么做呢？
我们发现此时欧拉函数的做法已经不适用了。
我们考虑莫比乌斯函数。
在下面的过程中，我们不妨假设 \(n\le m\)。
枚举 \(\gcd\)：
还是同时除以 \(d\)，让式子中出现 \(=1\) 的形式：
接下来，我们进行最关键的一步：把 \(\sum_{d|n}\mu(d)=[n=1]\) 带进去：
别看现在这么多求和，其实一会就没了。
我们发现：\(t|\gcd(i,j)\Leftrightarrow t|i\and t|j\)，证明的话……显然。
接下来的一步，可能会让很多初学者头疼：调换求和顺序，当同时出现很多个求和的时候，我们就可以考虑这个方法。
为了简洁一点，我们假设 \(x=\lfloor\frac nd\rfloor,y=\lfloor\frac md\rfloor\)。
我们考虑枚举 \(t\) 代替枚举 \(i,j\)，我们发现要满足 \(i,j\) 都是 \(t\) 的倍数，就能产生 \(1\) 的贡献。
我们发现，\(i\) 是 \(t\) 的倍数和 \(j\) 是 \(t\) 的倍数这两个条件是独立的，所以我们可以把它拆开。
我们发现，\([1,x]\) 里的 \(t\) 的倍数总共有 \(\lfloor\frac xt\rfloor\) 个，所以我们可以把这个带进去。
这个时候，事情就差不多了。
我们可以发现，后面这个式子可以用到整除分块！
我们预处理除 \(\mu(t)\) 的前缀和，然后对于 \(\lfloor\frac xt\rfloor,\lfloor\frac yt\rfloor\) 都不变的区间 \([L,R]\)，我们可以直接求出右边式子的答案，这样可以做到时间复杂度 \(O(n\sqrt n)\)。
我们发现还不够优，把 \(x,y\) 带回来，\(\lfloor\frac nd\rfloor,\lfloor\frac md\rfloor\) 同样可以整除分块，这两个值一样的区间，对应右边的值一样。
对于一个区间 \(d\) 的求和，不要告诉我你不会！
这样就可以做到总时间复杂度 \(O(\sqrt n\cdot\sqrt n)=O(n)\) 了。
看到这里，你可能觉得这个好难，但其实做多了你会发现这个事情非常套路，最后变得人均都会了。
我们再来一题：
问题 5.5.2（P1829 [国家集训队]Crash的数字表格 / JZPTAB）：给定 \(n,m\)，求 \(\sum_{i=1}^n\sum_{j=1}^m \text{lcm}(i,j)\)。
熟悉的套路，熟悉的配方，让我们用熟悉的方式过掉这道题吧：
我们同样假设 \(n\le m\)。
接下来无脑上求和公式和整除分块就可以了。
#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=1e7+10,p=20101009;
int n,m,mu[N],pri[N],tot;
bool flg[N];
int sum(int l,int r){
return ((l+r)*(r-l+1ll)/2ll)%p;
}
int gts(int n,int m){
int s=0;
for(int t=1,q=0;t<=n;t=q+1){
int x=n/t,y=m/t;
q=min(n/x,m/y);
// printf("%lld %lld %lld %lld\n",n,m,t,q);
s=(s+(mu[q]-mu[t-1]+p)%p*sum(1,x)%p*sum(1,y)%p)%p;
}
return s;
}
signed main(){
scanf("%lld%lld",&n,&m);
if(n>m)swap(n,m);
mu[1]=1;
for(int i=2;i<=n;i++){
if(!flg[i])pri[++tot]=i,mu[i]=p-1;
for(int j=1;j<=tot&&pri[j]*i<=n;j++){
flg[i*pri[j]]=1;
if(i%pri[j]==0)break;
mu[i*pri[j]]=p-mu[i];
}
mu[i]=(mu[i]*i%p*i%p+mu[i-1])%p;
}
int ans=0;
for(int d=1,q=0;d<=n;d=q+1){
int x=n/d,y=m/d;
q=min(n/x,m/y);
// printf("sao:%lld %lld %lld %lld\n",n,m,d,q);
ans=(ans+sum(d,q)*gts(x,y))%p;
}
printf("%lld\n",ans%p);
}
几道练习：
P5518 [MtOI2019]幽灵乐团 / 莫比乌斯反演基础练习题
6 二项式反演
经过刚刚莫比乌斯反演的学习，我们应该对反演有了一个基础的认识。
比如我们知道一个数列 \(\{f_i\}\) 和 \(\{g_i\}\) 的关系满足：\(g_n=\sum_{i=0}^na_{n,i}f_i\)，我们希望通过 \(\{g_i\}\) 得到 \(\{f_i\}\)，也就是说，我们希望通过数组 \(\{b_{n,i}\}\) 满足 \(f_n\sum_{i=0}^nb_{n,i}g_i\)，这样我们就能够通过 \(g_i\) 求出 \(f_i\) 了。
6.1 二项式反演的认识及证明
二项式反演大概讲的是这样的一个东西：
我们发现这个式子两边长得非常像，而且它还有一个等价式比较常用：
证明的话我们只需要把左边的式子带到右边去：
我们把组合数暴力展开：
交换求和顺序：
我们发现后面的求和长得非常丑陋，所以考虑转化一下，把所有的 \(i\) 变成 \(n-i\)：
然后我们发现这个东西其实就是性质 1.3.6，所以转化成 \([n-j=0]\)，从而转化成 \([n=j]\)，在求和中只有一次满足条件，所以整个式子的值为 \(1\)。
让我们回到原式：
证毕。
6.2 二项式反演的运用
大部分二项式反演的题都可以用容斥去推，具体可以看个人熟练度、喜好之类的决定用什么方法，最后推出来的式子一般都一样。
问题 6.2.1（错排问题）：求有多少个长度为 \(n\) 的排列 \(\{p_i\}\)，满足 \(\forall i\in[1,n],p_i\not=i\)。
问题其实就是求所有数字的装错位置的排列数。
我们记错排数为 \(\{D_i\}\)，根据前面的知识，我们可以发现随便装的方案为 \(n!\)。
然后对于这 \(n!\) 种方案，可能会装错 \(0,1,2,\dots,n\) 个。
我们枚举装错了几个，每次钦定若干个位置装错，统计起来，答案应等于 \(n!\)。
发现它就是一个二项式反演，直接套公式：
拆开组合数，化简一下：
这就是我们常说的错排公式了。
当然，错排公式的推导过程还有很多，生成函数和容斥都可以推导，这里不在赘叙。
问题 6.2.2：有一排长度为 \(n\) 的格子，用 \(m(m\ge2)\) 种颜色涂色，要求相邻两个不相同，且每种颜色都被使用，求方案数。
这题直接容斥会很简单，直接枚举几种颜色没使用即可，由于我们现在在学二项式反演，所以还是用二项式反演来推。
我们记题目答案为 \(h(m,n)\)，如果我们不管所有颜色都被使用，我们会得到异色方案数涂色为 \(m(m-1)^{n-1}\)，含义就是第一个格子乱填，后面每一个都要排除前面的颜色，所以方案数为 \(m(m-1)^{n-1}\)。
我们同样枚举集中颜色没使用，每次钦定若干个颜色不使用，全部加起来也应该等于 \(m(m-1)^{n-1}\)。
我们发现这个和二项式反演略有区别，我们可以把边界保留然后按照刚才的证明思路再推一遍，会得到：
注意特判 \(n\) 比较小的情况。
问题 6.2.3（CF1228E Another Filling the Grid）：有一个 \(n\times n\) 的矩阵，每个格子填 \(1\) 到 \(k\) 的一个整数，要求每行每列的最小值都为 \(1\) ,求方案数对 \(1000000009\) 取模的结果。
如果 \(k=1\) 直接特判，我们考虑 \(k\ge2\) 的情况。
我们记 \(f(i,j)\) 表示恰有 \(i\) 行 \(j\) 列不合法的方案数，总方案数为 \(k^{n\times n}\) 中方案数，所以能得出式子：
这个是二维的二项式反演，反演方法也是一样的：
如果考虑容斥含义，那么后面那个 \(ni+nj-ij\) 其实就是 \(i\) 行 \(j\) 列覆盖的范围，这些地方强制不合法，要减去 \(1\) 的情况，每一个有 \((k-1)\) 种，其他乱选。
#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=1010,p=1e9+7;
int n,k,c[N][N];
int qpow(int x,int y=p-2){
int m=1;
for(;y;y>>=1,x=x*x%p)if(y&1)m=m*x%p;
return m;
}
signed main(){
scanf("%lld%lld",&n,&k);
c[0][0]=1;
for(int i=1;i<=n;i++){
c[i][0]=1;
for(int j=1;j<=i;j++)c[i][j]=(c[i-1][j-1]+c[i-1][j])%p;
}
int ans=0;
for(int i=0;i<=n;i++){
for(int j=0;j<=n;j++){
ans=(ans+c[n][i]*c[n][j]%p*((i+j&1)?p-1:1)%p*qpow(k-1,n*i+n*j-i*j)%p*qpow(k,n*n-n*i-n*j+i*j))%p;
}
}
printf("%lld\n",ans);
}
