莫比乌斯反演
重新学习了一遍莫比乌斯反演,整理一下。
莫比乌斯函数
莫比乌斯函数\(\mu\)是一个积性函数。
即对于一个数\(x\)的莫比乌斯函数分三种情况讨论
\(1.\)当\(x=1\)时,\(\mu(x)=1\)
\(2.\)当\(x\)为互异素数的乘积时\(\mu(x)=(-1)^k\),\(k\)是\(x\)分解后的素数个数。
\(3.\)其他情况\(\mu(x)=0\)
莫比乌斯函数是个积性函数,所以可以线性筛。
性质
\(1.\) $$\sum_{d|n}\mu(d)=[n=1]$$
\(2.\) $$\sum_{d|n}\frac{\mu(d)}{d}=\frac{\phi(n)}{n}$$
线性筛莫比乌斯函数代码
for(int i = 2;i <= END;++i) {
if(!vis[i]) {
mu[i] = -1;
dis[++js] = i;
}
for(int j = 1;j <= js && dis[j] * i <= END;++j) {
vis[i * dis[j]] = 1;
if(i % dis[j]) mu[i * dis[j]] = -mu[i];
else {
mu[i * dis[j]] = 0;
break;
}
}
}
莫比乌斯反演
其实莫比乌斯反演就是两个公式
设\(F(n)\)和\(f(n)\)是定义在非负整数上的两个函数
莫比乌斯反演仅仅靠这两个公式是很难做题的,做完后面两道例题应该就算入门了
公式一
若\(F(n)\)和\(f(n)\)满足
则有$$f(n)=\sum\limits_{d|n}\mu(d)F(\lfloor\frac{n}{d}\rfloor)$$
公式二
若\(F(n)\)和\(f(n)\)满足
则有
证明
留坑待填吧\(233\)
题目
例题1
求出有多少\(x,y\)满足\(1\le x\le B\quad 1\le y \le D\)且\(gcd(x,y) = k\)
评测戳这里
这是非常基础的一道题目。为了便于讨论,我们下面都假定\(B\le D\)
我们设\(f(n)\)表示\(gcd(x,y)=n\)的数量,\(F(n)\)表示\(gcd(x,y)\% n= 0\)的数量
显然有下面的式子
然后根据公式二,可以反演出下面的式子
现在\(f(k)\)就是答案,如果我们让\(B\)和\(D\)都先除以\(k\),那么\(f(1)\)就是答案。
所以下面我们假定\(B=\frac{B}{k},D=\frac{D}{k}\)
现在就是求$$f(1)=\sum\limits_{d=1}^B\mu(d)F(d)$$
然后现在只要能够快速的计算出\(F(d)\)就能够\(O(B)\)的计算了。
显然有$$F(i)=\lfloor\frac{B}{i}\rfloor\lfloor\frac{D}{i}\rfloor$$
所以$$f(1)=\sum\limits_{d=1}^B\mu(d)\lfloor\frac{B}{d}\rfloor\lfloor\frac{D}{d}\rfloor$$
然后计算就行了。因为题目中说要去重,发现被计算两遍的那些数字是\(x,y\) 都小于\(B\)的情况,所以最后在减去就行了。
代码如下
/*
* @Author: wxyww
* @Date: 2019-02-22 11:00:34
* @Last Modified time: 2019-02-22 11:24:13
*/
#include<cstdio>
#include<iostream>
#include<cstdlib>
#include<cstring>
#include<algorithm>
#include<queue>
#include<vector>
#include<ctime>
using namespace std;
typedef long long ll;
#define int ll
const int N = 100000 + 100;
ll read() {
ll x=0,f=1;char c=getchar();
while(c<'0'||c>'9') {
if(c=='-') f=-1;
c=getchar();
}
while(c>='0'&&c<='9') {
x=x*10+c-'0';
c=getchar();
}
return x*f;
}
int mu[N],dis[N],vis[N];
void pre() {
mu[1] = 1;
int js = 0;
int END = 100000;
for(int i = 2;i <= END;++i) {
if(!vis[i]) {
mu[i] = -1;
dis[++js] = i;
}
for(int j = 1;j <= js && dis[j] * i <= END;++j) {
vis[i * dis[j]] = 1;
if(i % dis[j]) mu[i * dis[j]] = -mu[i];
else {
mu[i * dis[j]] = 0;
break;
}
}
}
}
int solve(int x,int y) {
int ans = 0;
for(int i = 1;i <= min(x,y);++i)
ans += mu[i] * (x / i) * (y / i);
return ans;
}
signed main() {
int T = read();
pre();
for(int i = 1;i <= T;++i) {
read();int n = read();read();int m = read();int K = read();
int z = min(n,m);
if(K != 0)
printf("Case %lld: %lld\n",i,solve(n / K,m / K) - solve(z / K,z / K) / 2);
else printf("Case %lld: 0\n",i);
}
return 0;
}
例题2
直接看原题面吧
首先一个明显的转化就是,转化成对于一个数字\(k\),有多少\(1\le x \le n,1\le y\le m\)满足\(gcd(x,y)=k\)。然后对于每个\(k\),计算答案就行了。
发现上面这个式子不就是上面那道题么。同样设\(n \le m\)
我们设\(f(i)\)表示\(gcd(x,y)=i\)时的答案。\(fb_i\)斐波那契数列第\(i\)项那么这道题最终的答案就是。$$\prod\limits_{i=1}nfb_i$$
根据上一题推出的式子,我们可以推出来
原式=$$\prod_{i=1}nfb_i^{n/i}\mu(d)\lfloor\frac{n}{id}\rfloor\lfloor\frac{m}{id}\rfloor}$$
设\(T=id\)
原式=
到了这里,对于外面的\(\prod\limits_{T=1}^n...^{\lfloor\frac{n}{T}\rfloor\lfloor\frac{m}{T}\rfloor}\)可以数论分块做到\(\sqrt{m}\)
对于里面\(\prod\limits_{d|T}fb_d^{\mu(T/d)}\),暴力预处理是\(nlog\)的(调和级数)。
代码如下
/*
* @Author: wxyww
* @Date: 2019-02-22 15:09:25
* @Last Modified time: 2019-02-24 08:45:10
*/
#include<cstdio>
#include<iostream>
#include<cstdlib>
#include<cstring>
#include<algorithm>
#include<queue>
#include<vector>
#include<ctime>
using namespace std;
typedef long long ll;
#define int ll
const int N = 1000000 + 10,mod = 1e9 + 7;
ll read() {
ll x=0,f=1;char c=getchar();
while(c<'0'||c>'9') {
if(c=='-') f=-1;
c=getchar();
}
while(c>='0'&&c<='9') {
x=x*10+c-'0';
c=getchar();
}
return x*f;
}
int f[N],mu[N],dis[N],vis[N],g[N];
int qm(ll x,int y) {
int ret = 1;
if(y < 0) y += mod - 1;
for(;y;y >>= 1,x = 1ll * x * x % mod)
if(y & 1) ret = 1ll * ret * x % mod;
return ret;
}
int tmp[5];
void pre() {
int END = 1000000;
mu[1] = 1;
f[1] = 1,g[0] = 1;
g[1] = 1;
int js = 0;
for(int i = 2;i <= END;++i) {
f[i] = f[i - 1] + f[i - 2];
g[i] = 1;
f[i] >= mod ? f[i] -= mod : 0;
if(!vis[i]) {
dis[++js] = i;
mu[i] = -1;
}
for(int j = 1;j <= js && dis[j] * i <= END;++j) {
vis[i * dis[j]] = 1;
if(i % dis[j]) mu[i * dis[j]] = -mu[i];
else break;
}
}
for(int i = 1;i <= END;++i) {
tmp[0] = qm(f[i],-1);
tmp[1] = 1;
tmp[2] = f[i];
for(int j = i;j <= END;j += i)
g[j] = 1ll * g[j] * tmp[mu[j / i] + 1] % mod;
}
for(int i = 1;i <= END;++i) g[i] = 1ll * g[i - 1] * g[i] % mod;
}
signed main() {
pre();
int T = read();
while(T--) {
int n = read(),m = read();
if(n > m) swap(n,m);
int r;
ll ans = 1;
for(int l = 1;l <= n;l = r + 1) {
r = min(n / (n / l),m / (m / l));
ans = ans * qm(g[r] * qm(g[l - 1],mod - 2) % mod , (n / l) * (m / l) % (mod - 1)) % mod;
}
printf("%lld\n",ans);
}
return 0;
}