莫比乌斯反演学习笔记

引入

对于数论问题中的一些函数 \(f(n)\)，如果很难直接求出它的值，却容易求出其倍数和或约数和 \(g(n)\)，那么可以通过莫比乌斯反演化简运算，求得 \(f(n)\) 的值。

定义

\(\mu\) 为莫比乌斯函数，定义如下：

\[\mu(n)=\begin{cases} 1，n=1 \\ 0，n含有平方因子 \\ (-1)^k，k为 n的本质不同质因子个数 \end{cases}\]

令 \(n=\prod_{i=1}^{k}p_i^{c_i}\)，其中 \(p_i\) 为质因子，\(c_i \geq 1\)。上述定义表示：

1.当 \(n=1\) 时，\(\mu(n)=1\)；

2.当 \(n \ne 1\) 时：

• 当存在 \(i \in [1,k]\)，使得 \(c_i>1\) 时，\(\mu(n)=0\)，也就是只要有某个质因子出现次数超过一次，\(\mu(u)\) 就等于 \(0\)。

• 当任意 \(i \in [1,k]\)，都满足 \(c_i=1\) 时，\(\mu(n)=(-1)^k\)，\(k\) 即为出现的质因子的总个数。

代码实现

根据定义，\(\mu(1) \sim \mu(n)\) 可以在线性筛的过程中同时得到，时间复杂度为 \(O(n)\)：

mu[1]=1;
for(int i=2;i<=n;i++)
{
	if(!vis[i]) prime[++tot]=i,mu[i]=-1;
	for(int j=1;j<=tot&&i*prime[j]<=n;j++)
	{
		vis[i*prime[j]]=1;
		if(i%prime[j]==0)
		{
			mu[i*prime[j]]=0;
			break;
		} 
		mu[i*prime[j]]=-mu[i];
	}
}

性质

\[\sum_{d|n} \mu(d)=\left\{\begin{matrix}1，n=1 \\ 0，n \ne 1 \end{matrix}\right.\]

证明

令 \(n=\prod_{i=1}^{k}p_i^{c_i},n'=\prod_{i=1}^{k}p_i\)

那么 \(\sum_{d|n} \mu(d)=\sum_{d|n'} \mu(d)=\sum_{i=0}^{k} \binom{k}{i} (-1)^i=(1+(-1))^k\)

根据二项式定理，只有当 \(k=0\) 即 \(n=1\) 时原式值为 \(1\)，反之则为 \(0\)。

结论

\(\sum_{d|\gcd(i,j))}\mu(d)=[\gcd(i,j)=1]\)

这个等式的证明比较显然，根据前面提到的性质，只有当 \(\gcd(i,j)=1\) 时，才满足\(\mu(d)=1\)，反之 \(\mu(d)=0\)。

莫比乌斯变换

形式一

\[f(n)=\sum_{d|n} g(d)\Longleftrightarrow g(n)=\sum_{d|n}\mu(d)f(\frac{n}{d}) \]

形式二

\[f(n)=\sum_{n|d} g(d)\Longleftrightarrow g(n)=\sum_{n|d}\mu(\frac{d}{n})f(d) \]

「HAOI 2011」Problem b

题意

\(T\) 组询问，每次询问求值：

\[\sum_{i=x}^{n}\sum_{i=y}^{m}[gcd(i,j)=k] \]

\(1 \leq T,x,y,n,m,k \leq 5 \times 10^4\)。

思路

根据差分的思路，可以将每次询问拆成四个相同形式的式子：

\[\sum_{i=1}^{n}\sum_{i=1}^{m}[gcd(i,j)=k] \]

化简该式子：

\[\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor }\sum_{i=1}^{\left \lfloor \frac{m}{k} \right \rfloor }[gcd(i,j)=1] \]

利用莫比乌斯函数的性质，进一步化简，得到：

\[\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor }\sum_{i=1}^{\left \lfloor \frac{m}{k} \right \rfloor } \sum_{d|\gcd(i,j)} \mu(d) \]

变换求和顺序，先枚举 \(d | \gcd(i,j)\) 可得：

\[\sum_{d=1} \mu(d)\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor } [d|i]\sum_{i=1}^{\left \lfloor \frac{m}{k} \right \rfloor } [d|j] \]

而 \(1 \sim \left \lfloor \frac{n}{k} \right \rfloor\) 中，\(d\) 的倍数有 \(\left \lfloor \frac{n}{kd} \right \rfloor\) 个，故原式化为：

\[\sum_{d=1}^{\min(\left \lfloor \frac{n}{k} \right \rfloor,\left \lfloor \frac{m}{k} \right \rfloor)} \mu(d)\left \lfloor \frac{n}{kd} \right \rfloor \left \lfloor \frac{m}{kd} \right \rfloor \]

于是这个式子就可以用数论分块求解。时间复杂度为 \(O(T \sqrt{n})\)

code：

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
const int N=5e4+10;
int prime[N],mu[N],tot;bool vis[N];
void init()
{
	mu[1]=1;
	for(int i=2;i<N;i++)
	{
		if(!vis[i]) prime[++tot]=i,mu[i]=-1;
		for(int j=1;j<=tot&&prime[j]*i<N;j++)
		{
			vis[prime[j]*i]=true;
			if(i%prime[j]==0) {mu[i*prime[j]]=0;break;}
			mu[i*prime[j]]=-mu[i];
		}
	}
	for(int i=2;i<N;i++) mu[i]+=mu[i-1];
}
int solve(int n,int m)
{
	int res=0;
	for(int l=1,r;l<=min(n,m);l=r+1)
	{
		r=min(n/(n/l),m/(m/l));
		res+=(mu[r]-mu[l-1])*(n/l)*(m/l);
	}
	return res;
}
int main()
{
	int T,n,m,x,y,k;scanf("%d",&T);init();
	while(T--)
	{
		scanf("%d%d%d%d%d",&x,&n,&y,&m,&k);
	    printf("%d\n",solve(n/k,m/k)-solve((x-1)/k,m/k)-solve(n/k,(y-1)/k)+solve((x-1)/k,(y-1)/k));
	}
	return 0;
}

Crash的数字表格

求值：

\[\sum_{i=1}^{n}\sum_{j=1}^{m}\mathrm{lcm} (i,j) \]

对 20101009$ 取模。

\(1 \leq n,m \leq 10^7\)。

思路

不难发现，原式等价于：

\[\sum_{i=1}^{n}\sum_{j=1}^{m}\frac{i \cdot j}{\mathrm{gcd} (i,j)} \]

枚举 \(d=\gcd(i,j)\)：

\[\sum_{d=1}^{n}\sum_{i=1}^{\left \lfloor \frac{n}{d} \right \rfloor }\sum_{j=1}^{\left \lfloor \frac{m}{d} \right \rfloor }[\gcd(i,j)=1]\ i \cdot j \cdot d \]

把 \(d\) 提到前面来，套用莫比乌斯函数：

\[\sum_{d=1}^{n}d\sum_{i=1}^{\left \lfloor \frac{n}{d} \right \rfloor }\sum_{j=1}^{\left \lfloor \frac{m}{d} \right \rfloor }\sum_{d|\gcd(i,j)} \mu(d) i \cdot j \]

单独计算后面两项的值，记

\[\mathrm{sum}(n,m)=\sum_{i=1}^{n}\sum_{j=1}^{m}\sum_{d|\gcd(i,j)} \mu(d) i \cdot j \]

变换求和顺序，先枚举 \(d | \gcd(i,j)\) 可得：

\[\sum_{d=1} \sum_{d|i}^{n}\sum_{d|j}^{m}\mu(d) i \cdot j \]

令 \(i=i' \cdot d,j=j' \cdot d\)，带入式子可得：

\[\sum_{d=1} \mu(d) d^2 \sum_{i=1}^{\left \lfloor \frac{n}{d} \right \rfloor }\sum_{j=1}^{\left \lfloor \frac{m}{d} \right \rfloor } i \cdot j \]

注意到这个式子前半部分可以预处理前缀和用整除分块做，后半部分则可以直接计算，那么原式就是：

\[\sum_{i=1}^{n} d \times \mathrm{sum}(\left \lfloor \frac{n}{d} \right \rfloor ,\left \lfloor \frac{m}{d} \right \rfloor) \]

于是总的时间复杂度就为 \(O(n+m)\)。

code：

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
const int N=1e7+10,mod=20101009;
int mu[N],sum[N],prime[N],tot;bool vis[N];
void init()
{
	mu[1]=1;
	for(int i=2;i<N;i++)
	{
		if(!vis[i]) prime[++tot]=i,mu[i]=-1;
		for(int j=1;j<=tot&&prime[j]*i<N;j++)
		{
			vis[prime[j]*i]=true;
			if(i%prime[j]==0) {mu[i*prime[j]]=0;break;}
			mu[i*prime[j]]=-mu[i];
		}
	}
	for(int i=1;i<N;++i) sum[i]=(sum[i-1]+1ll*i*i%mod*(mu[i]+mod))%mod;
}
int Sum(int n,int m){return (1ll*n*(n+1)/2%mod)*(1ll*m*(m+1)/2%mod)%mod;}
int func(int n,int m)
{
    int res=0;
    for(int l=1,r;l<=min(n,m);l=r+1)
	{
        r=min(n/(n/l),m/(m/l));
        res=(res+1ll*(sum[r]-sum[l-1]+mod)*Sum(n/l,m/l)%mod)%mod;
    }
    return res;
}
int solve(int n,int m)
{
	int res=0;
    for(int l=1,r;l<=min(n,m);l=r+1)
	{
        r=min(n/(n/l),m/(m/l));
        res=(res+1ll*(r-l+1)*(l+r)/2%mod*func(n/l,m/l)%mod)%mod;
    }
    return res;
}
int main()
{
	int n,m;scanf("%d%d",&n,&m);init();printf("%d\n",solve(n,m));
	return 0;
}

约数个数和

\(T\) 组询问。每次询问给定 \(n,m\)，求：

\[\sum_{i=1}^n\sum_{j=1}^md(i \cdot j) \]

其中 \(d(x)\) 为 \(x\) 的约数个数。

\(1\le T,n,m \le 50000\)

思路：

首先需要用到 \(d\) 函数的一个特殊性质，即：

\[d(i \cdot j)=\sum_{x|i}\sum_{y|i} [\gcd(x,y)=1] \]

考虑将这个式子化简：

\[d(i \cdot j)=\sum_{x|i}\sum_{y|i} [\gcd(x,y)=1] \]

\[=\sum_{x|i}\sum_{y|i} \sum_{p|\gcd(x,y)}\mu(p) \]

\[=\sum_{p=1}\mu(p) \sum_{x|i}\sum_{y|i} [p|\gcd(x,y)] \]

\[=\sum_{p|i , p|j}\mu(p) \sum_{x|i}\sum_{y|i} [p|\gcd(x,y)] \]

\[=\sum_{p|i , p|j}\mu(p) \sum_{x|{\frac{i}{p}}}\sum_{y|{\frac{j}{p}}} 1 \]

\[=\sum_{p|i , p|j}\mu(p) d(\frac{i}{p})d(\frac{j}{p}) \]

带回原式，得：

\[\sum_{i=1}^n\sum_{j=1}^md(i \cdot j) \]

\[=\sum_{i=1}^n\sum_{j=1}^m\sum_{p|i , p|j}\mu(p) d(\frac{i}{p})d(\frac{j}{p}) \]

\[=\sum_{p=1}^{\min(n,m)}\mu(p)\sum_{i=1}^n\sum_{j=1}^m [p|i][p|j] d(\frac{i}{p})d(\frac{j}{p}) \]

\[=\sum_{p=1}^{\min(n,m)}\mu(p)\sum_{i=1}^{\left \lfloor \frac{n}{p} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{m}{p} \right \rfloor} d(i)d(j) \]

\[=\sum_{p=1}^{\min(n,m)}\mu(p)\sum_{i=1}^{\left \lfloor \frac{n}{p} \right \rfloor} d(i)\sum_{j=1}^{\left \lfloor \frac{m}{p} \right \rfloor} d(j) \]

对于 \(d\) 函数，也可以在线性筛的同时求出。那么就可以做到 \(O(n)\) 预处理，\(O(T\sqrt{n})\) 询问。

code：

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
const int N=5e4+10;
#define LL long long
LL d[N];int t[N],mu[N],prime[N],tot;bool vis[N];//t表示最小质因子的出现次数 
void init()
{
	mu[1]=d[1]=1;
	for(int i=2;i<N;i++)
	{
		if(!vis[i]) prime[++tot]=i,d[i]=2,mu[i]=-1,t[i]=1;
		for(int j=1;j<=tot&&prime[j]*i<N;j++)
		{
			vis[prime[j]*i]=true;
			if(i%prime[j]==0)
			{
				mu[i*prime[j]]=0;
				d[i*prime[j]]=d[i]/(t[i]+1)*(t[i]+2);
 			    t[i*prime[j]]=t[i]+1;
 			    break;
			}
			mu[i*prime[j]]=-mu[i];
			d[i*prime[j]]=d[i]*2;
		    t[i*prime[j]]=1;
		}
	}
	for(int i=1;i<N;i++) mu[i]+=mu[i-1],d[i]+=d[i-1];
}
LL solve(int n,int m)
{
	LL res=0;
	for(int l=1,r;l<=min(n,m);l=r+1)
	{
		r=min(n/(n/l),m/(m/l));
		res+=(mu[r]-mu[l-1])*d[n/l]*d[m/l];
	}
	return res;
}
int main()
{
	int T,n,m;scanf("%d",&T);init();
	while(T--) scanf("%d%d",&n,&m),printf("%lld\n",solve(n,m));
	return 0;
}