莫比乌斯反演 习题

$Tips$：本文默认 $n\le m,a\le b$

[POI2007] ZAP-Queries

（题目传送门）

题意

给出 $a,b,d$，求满足 $1\le x\le a$，$1\le y\le b$，且 $gcd(x,y)=d$ 的二元组 $(x,y)$ 的数量。

$1\le n\le 5\times {10}^4$，$1\le d\le a,b\le 5\times {10}^4$。

解题思路

• 原题要求的是

  $$\sum _{i=1}^a\sum _{j=1}^b[gcd(i,j)=d]$$

• 将 $i,j$ 除以 $d$ 得

  $$=\sum _{i=1}^{\lfloor \frac{a}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{b}{d}\rfloor }[gcd(i,j)=1]$$

• 由莫比乌斯反演得

  $$=\sum _{i=1}^{\lfloor \frac{a}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{b}{d}\rfloor }\sum _{k|gcd(i,j)}\mu (k)$$

• 再将 $i,j$ 除以 $k$，这样就不关心 $k$ 是否为 $i,j$ 的约数，可以直接枚举 $k$ ：

  $$=\sum _{k=1}^{\lfloor \frac{a}{d}\rfloor }\mu (k)\lfloor \frac{a}{kd}\rfloor \lfloor \frac{b}{kd}\rfloor$$

• 这个式子的后半部分可以整除分块来优化，前半部分预处理前缀和即可。在代码实现时，可以预先让 $a,b$ 同除以 $d$，方便实现

代码

#include<bits/stdc++.h>
#define LL long long
using namespace std;

const int N=50010;

int T,a,b,d,n,m;
int prime[N],v[N],tot,mu[N],sum[N];
LL ans;

void primes()
{
	mu[1]=1;
	for(int i=2; i<=50000; i++)
	{
		if(!v[i])
		{
			v[i]=i;
			prime[++tot]=i;
			mu[i]=-1;
		}
			
		for(int j=1; j<=tot; j++)
		{
			if(prime[j]>v[i] || prime[j]>50000/i)
				break;
			v[i*prime[j]]=prime[j];
			if(i%prime[j]==0)
				break;
			mu[prime[j]*i]=-mu[i];
		}
	}
	
	for(int i=1; i<=50000; i++)
		sum[i]=sum[i-1]+mu[i];
}
 	
int main()
{
	primes();
	
    scanf("%d",&T);
    while(T--)
    {
    	scanf("%d%d%d",&a,&b,&d);
    	n=a/d;  m=b/d;
    	if(n>m)
    		swap(n,m);
    	
    	ans=0;
    	for(int l=1,r; l<=n; l=r+1)
    	{
    		r=min(n/(n/l),m/(m/l));
    		ans+=1LL*(sum[r]-sum[l-1])*1LL*(n/l)*1LL*(m/l);
		}
		
		printf("%lld\n",ans);
	}
    
    return 0;
}

[HAOI2011] Problem b

（题目传送门）

题意

对于给出的 $n$ 个询问，每次求有多少个数对 $(x,y)$，满足 $a\le x\le b$，$c\le y\le d$，且 $gcd(x,y)=k$，$gcd(x,y)$ 函数为 $x$ 和 $y$ 的最大公约数。

$1\le n,k\le 5\times {10}^4$，$1\le a\le b\le 5\times {10}^4$，$1\le c\le d\le 5\times {10}^4$。

解题思路

• 类比上面一题 ZAP，设 $f(n,m)$ 表示满足 $1\le x\le n$，$1\le y\le m$，且 $gcd(x,y)=k$ 的二元组 $(x,y)$ 的数量

• 类似二维前缀和，容易得出

  $$ans=f(b,d)-f(a-1,d)-f(b,c-1)+f(a-1,c-1)$$

• 同样，可以一开始先令 $a,b,c,d$ 除以 $k$，方便代码实现

代码

#include<bits/stdc++.h>
#define LL long long
using namespace std;

const int N=50010;

int T,a,b,c,d,k;
int prime[N],v[N],tot,mu[N],sum[N];
LL ans;

void primes()
{
	mu[1]=1;
	for(int i=2; i<=50000; i++)
	{
		if(!v[i])
		{
			v[i]=i;
			prime[++tot]=i;
			mu[i]=-1;
		}
			
		for(int j=1; j<=tot; j++)
		{
			if(prime[j]>v[i] || prime[j]>50000/i)
				break;
			v[i*prime[j]]=prime[j];
			if(i%prime[j]==0)
				break;
			mu[prime[j]*i]=-mu[i];
		}
	}
	
	for(int i=1; i<=50000; i++)
		sum[i]=sum[i-1]+mu[i];
}

LL f(int n,int m)
{
	if(n>m)
		swap(n,m);
	LL res=0;
	for(int l=1,r; l<=n; l=r+1)
   	{
    	r=min(n/(n/l),m/(m/l));
    	res+=1LL*(sum[r]-sum[l-1])*1LL*(n/l)*1LL*(m/l);
	}
	
	return res;
}
 	
int main()
{
	primes();
	
    scanf("%d",&T);
    while(T--)
    {
    	scanf("%d%d%d%d%d",&a,&b,&c,&d,&k);
    	
    	ans=0;
    	ans=f(b/k,d/k)-f((a-1)/k,d/k)-f(b/k,(c-1)/k)+f((a-1)/k,(c-1)/k);
		
		printf("%lld\n",ans);
	}
    
    return 0;
}

YY的GCD

（题目传送门）

题意

$T$ 组数据，每组数据给定 $N,M$，求 $1\le x\le N$，$1\le y\le M$ 且 $gcd(x,y)$ 为质数的 $(x,y)$ 有多少对。

$T={10}^4$，$N,M\le {10}^7$。

解题思路

• 设 $primes$ 表示质数集合

• 原题要求的是

  $$ans=\sum _{k\in primes}\sum _{i=1}^n\sum _{j=1}^m[gcd(i,j)=k]$$

• 同除以 $k$ 得：

  $$=\sum _{k\in primes}\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{k}\rfloor }[gcd(i,j)=k]$$

• 由莫比乌斯反演得

  $$=\sum _{k\in primes}\sum _{i=1}^{\lfloor \frac{n}{k}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{k}\rfloor }\sum _{d|gcd(i,j)}\mu (d)$$

• 我们枚举 $d$ 得：

  $$=\sum _{k\in primes}\sum _{d=1}^{\lfloor \frac{n}{k}\rfloor }\mu (d)\lfloor \frac{n}{kd}\rfloor \lfloor \frac{m}{kd}\rfloor$$

• 设 $T=kd$，把 $T$ 提到前面去

  $$=\sum _{T=1}^n\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor \sum _{k|T,k\in primes}\mu (\frac{T}{k})$$

• 对于后半部分式子，我们可以在线性筛后预处理出来，再做一次整除分块即可求出答案

代码

#include<bits/stdc++.h>
#define LL long long
using namespace std;

const int N=10000010;

int T,n,m;
int prime[N],v[N],tot,mu[N];
LL sum[N],g[N];

void primes()
{
	mu[1]=1;
	for(int i=2; i<=1e7; i++)
	{
		if(!v[i])
		{
			v[i]=i;
			prime[++tot]=i;
			mu[i]=-1;
		}
			
		for(int j=1; j<=tot; j++)
		{
			if(prime[j]>v[i] || prime[j]>1e7/i)
				break;
			v[i*prime[j]]=prime[j];
			if(i%prime[j]==0)
				break;
			mu[prime[j]*i]=-mu[i];
		}
	}
	
	for(int i=1; i<=tot; i++)
	{
		for(int j=1; j*prime[i]<=1e7; j++)
			g[prime[i]*j]+=mu[j]; 
	}
	
	for(int i=1; i<=1e7; i++)
		sum[i]=sum[i-1]+g[i];
}

LL cacl(int n,int m)
{
	if(n>m)
		swap(n,m);
	
	LL res=0; 
	for(int l=1,r; l<=n; l=r+1)
	{
		r=min(n/(n/l),m/(m/l));
		res+=1LL*(sum[r]-sum[l-1])*(n/l)*(m/l);
	}
	
	return res;
}
 	
int main()
{
	primes();
	
    scanf("%d",&T);
    while(T--)
    {
    	scanf("%d%d",&n,&m);
		printf("%lld\n",cacl(n,m));
	}
    
    return 0;
}

于神之怒加强版

（题目传送门）

题意

给定 $n,m,k$，计算

$$\sum _{i=1}^n\sum _{j=1}^{m}gcd(i,j{)}^k$$

对 ${10}^9+7$ 取模的结果。

对于全部的测试点，保证 $1\le T\le 2\times {10}^3$，$1\le n,m,k\le 5\times {10}^6$。

解题思路

• 直接开始推式子

  $$ans=\sum _{d=1}^n\sum _{i=1}^n\sum _{j=1}^m{d}^k[gcd(i,j)=d]$$

• 除以 $d$ 得：

  $$=\sum _{d=1}^n\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }{d}^k[gcd(i,j)=1]$$
  $$=\sum _{d=1}^n{d}^k\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }\sum _{x|gcd(i,j)}\mu (x)$$

• 枚举 $x$ 得：

  $$\sum _{d=1}^n{d}^k\sum _{x=1}^{\lfloor \frac{n}{d}\rfloor }\mu (x)\lfloor \frac{n}{dx}\rfloor \lfloor \frac{m}{dx}\rfloor$$

• 类似上一题的优化，设 $T=dx$ 得：

  $$\sum _{T=1}^n\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor \sum _{d|T}{d}^k\mu (\frac{T}{d})$$

• 观察到，$\sum _{d|T}{d}^k\mu (\frac{T}{d})$ 是狄利克雷卷积的形式，而 $f(d)=d^k$ 和 $\mu (\frac{T}{d})$ 又是 积性函数，所以 $g(T)=\sum _{d|T}{d}^k\mu (\frac{T}{d})$ 也是积性函数。因此，我们可以考虑线性筛预处理出 $g$ 的值

• 设 $T$ 质因数分解为 $T=\prod _{i=1}^q{p}_i^{c_i}$，则

  $$g(T)=\prod _{i=1}^{q}g(p_i^{c_i})$$
  $$=\prod _{i=1}^q\sum _{d|p_i^{c_i}}{d}^k\mu (\frac{p_i^{c_i}}{d})$$

• 观察到，当 $d=1,p_i,p{i}^2\dots p_i^{c_i-2}$ 时，$\mu (\frac{p_i^{c_i}}{d})=0$，所以

  $$g(T)=\prod _{i=1}^q{p}_i^{(c_i-1)k}\times \mu (p_i)+p_i^{c_{i}k}\times \mu (1)$$
  $$=\prod _{i=1}^q{p}_i^{(c_i-1)k}\times (p_i^k-1)$$

• 这个式子是可以预处理出来的

• 所以

  $$ans=\sum _{T=1}^n\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor g(T)$$

• 利用整除分块即可求解

代码

#include<bits/stdc++.h>
#define LL long long
using namespace std;

const int N=5000010;
const int MOD=1e9+7;

int T,k,n,m;
int prime[N],v[N],tot;
LL g[N],mi[N],sum[N];

LL ksm(int a,int b)
{
	if(b==0)
		return 1;
	if(b==1)
		return a;
		
	LL tmp=ksm(a,b/2);
	if(b%2==0)
		return (tmp*tmp)%MOD;
	else
		return ((tmp*tmp)%MOD*(LL)a)%MOD;
}

void primes()
{
	g[1]=1;
	for(int i=2; i<=5e6; i++)
	{
		if(!v[i])
		{
			v[i]=i;
			prime[++tot]=i;
			mi[i]=ksm(i,k);
			g[i]=(mi[i]-1+MOD)%MOD; //i为质数时g(i)=i^k-1
		}
			
		for(int j=1; j<=tot; j++)
		{
			if(prime[j]>v[i] || prime[j]>5e6/i)
				break;
			v[i*prime[j]]=prime[j];
			
			if(i%prime[j]==0)
				g[i*prime[j]]=(g[i]*mi[prime[j]])%MOD; //当i和prime[j]不互质时，读者可自行证明该结论正确性
			else
				g[i*prime[j]]=(g[i]*g[prime[j]])%MOD;
		}
	}
	
	for(int i=1; i<=5e6; i++)
		sum[i]=sum[i-1]+g[i];
}

LL cacl(int n,int m)
{
	if(n>m)
		swap(n,m);
	
	LL res=0; 
	for(int l=1,r; l<=n; l=r+1)
	{
		r=min(n/(n/l),m/(m/l));
		res=(res+((sum[r]-sum[l-1]+MOD)%MOD)*(n/l)%MOD*(m/l)%MOD)%MOD;
	}
	
	return res;
}
 	
int main()
{
	scanf("%d%d",&T,&k);
	
	primes();
	
    while(T--)
    {
    	scanf("%d%d",&n,&m);
		printf("%lld\n",cacl(n,m));
	}
    
    return 0;
}

[国家集训队] Crash的数字表格 / JZPTAB

（题目传送门）

题意

给定 $n,m$，计算

$$\sum _{i=1}^n\sum _{j=1}^m\mathrm{lcm}(i,j)$$

$1\le n,m\le {10}^7$

解题思路

• 对 $\mathrm{lcm}(i,j)$ 变换一下得：

  $$ans=\sum _{i=1}^n\sum _{j=1}^m\frac{ij}{gcd(i,j)}$$

• 枚举 $gcd(i,j)$ 得：

  $$=\sum _{i=1}^n\sum _{j=1}^m\sum _{d=1}^n\frac{ij}{d}[gcd(i,j)=d]$$

• 除以 $d$ 得：

  $$=\sum _{d=1}^n\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }ij\times d[gcd(i,j)=1]$$
  $$=\sum _{d=1}^{n}d\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }ij[gcd(i,j)=1]$$

• 把后半部分式子单独提出来，设：

  $$f(a,b)=\sum _{i=1}^a\sum _{j=1}^{b}ij[gcd(i,j)=1]$$

• 化简得：

  $$f(a,b)=\sum _{i=1}^a\sum _{j=1}^{b}ij\sum _{k|gcd(i,j)}\mu (k)$$
  $$=\sum _{k=1}^a\mu (k)\sum _{i=1}^{\lfloor \frac{a}{k}\rfloor }\sum _{j=1}^{\lfloor \frac{b}{k}\rfloor }ij\times k^2$$
  $$=\sum _{k=1}^a\mu (k)\times k^2\sum _{i=1}^{\lfloor \frac{a}{k}\rfloor }\sum _{j=1}^{\lfloor \frac{b}{k}\rfloor }ij$$

• 设

  $$g(a,b)=\sum _{i=1}^a\sum _{j=1}^{b}ij$$

• 化简得：

  $$=\frac{(1+a)a}{2}\times \frac{(b+1)b}{2}$$

• 所以，求出 $g(a,b)$ 的时间复杂度为 $O(1)$

• 将 $g(a,b)$ 代回 $f(a,b)$ 中得：

  $$f(a,b)=\sum _{k=1}^a\mu (k)\times k^2\times g(\lfloor \frac{a}{k}\rfloor ,\lfloor \frac{b}{k}\rfloor )$$

• 我们利用整除分块即可求出 $f(a,b)$

• 将 $f(a,b)$ 代回 $ans$ 得：

  $$ans=\sum _{d=1}^{n}d\times f(\lfloor \frac{n}{d}\rfloor ,\lfloor \frac{m}{d}\rfloor )$$

• 再做一次整除分块即可求解

代码

#include<bits/stdc++.h>
#define LL long long
using namespace std;

const int N=10000010;
const int MOD=20101009;

int n,m;
int prime[N],v[N],tot,mu[N];
LL sum[N];

void primes()
{
	mu[1]=1;
	for(int i=2; i<=1e7; i++)
	{
		if(!v[i])
		{
			v[i]=i;
			prime[++tot]=i;
			mu[i]=-1;
		}
			
		for(int j=1; j<=tot; j++)
		{
			if(prime[j]>v[i] || prime[j]>1e7/i)
				break;
			v[i*prime[j]]=prime[j];
			if(i%prime[j]==0)
				break;
			mu[i*prime[j]]=-mu[i];
		}
	}
	
	for(int i=1; i<=1e7; i++)
		sum[i]=(sum[i-1]+1LL*(mu[i]+MOD)%MOD*i%MOD*i%MOD)%MOD;
}

LL g(int a,int b)
{
	return (1LL*a*(a+1)/2%MOD)*(1LL*b*(b+1)/2%MOD)%MOD;
}

LL f(int a,int b)
{
	if(a>b)
		swap(a,b);
		
	LL res=0;
	for(int l=1,r; l<=a; l=r+1)
	{
		r=min(a/(a/l),b/(b/l));
		res=(res+(sum[r]-sum[l-1]+MOD)%MOD*g(a/l,b/l)%MOD)%MOD;
	}
	
	return res;
}

LL calc(int a,int b)
{
	if(a>b)
		swap(a,b);
		
	LL res=0;
	for(int l=1,r; l<=a; l=r+1)
	{
		r=min(a/(a/l),b/(b/l));
		res=(res+1LL*(r-l+1)*(l+r)/2%MOD*f(a/l,b/l)%MOD)%MOD;
	}
	
	return res;
}
 		
int main()
{
	primes();
	
	scanf("%d%d",&n,&m);
	printf("%lld",calc(n,m));
    
    return 0;
}

[SDOI2014] 数表

（题目传送门）

题意

有一张 $n\times m$ 的数表，其第 $i$ 行第 $j$ 列（$1\le i\le n$，$1\le j\le m$）的数值为能同时整除 $i$ 和 $j$ 的所有自然数之和。给定 $a$，计算数表中不大于 $a$ 的数之和。

$1\le n,m\le {10}^5$，$1\le Q\le 2\times {10}^4$。

解题思路

• 设 $f(x)$ 表示 $x$ 的因数之和

• 则

  $$ans=\sum _{i=1}^n\sum _{j=1}^{m}f(gcd(i,j))\times [f(gcd(i,j))\le a]$$

• 枚举 $gcd(i,j)$ 得：

  $$=\sum _{d=1}^n\sum _{i=1}^n\sum _{j=1}^{m}f(d)\times [gcd(i,j)=d]\times [f(d)\le a]$$
  $$=\sum _{d=1}^{n}f(d)\times [f(d)\le a]\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }[gcd(i,j)=1]$$

• 由莫比乌斯反演得：

  $$=\sum _{d=1}^{n}f(d)\times [f(d)\le a]\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }\sum _{k|gcd(i,j)}\mu (k)$$
  $$=\sum _{d=1}^{n}f(d)\times [f(d)\le a]\sum _{k=1}^{\lfloor \frac{n}{d}\rfloor }\mu (k)\lfloor \frac{n}{kd}\rfloor \lfloor \frac{m}{kd}\rfloor$$

• 套路设 $T=kd$ 得：

  $$=\sum _{T=1}^n\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor \sum _{x|T}\mu (\frac{T}{x})\times f(x)\times [f(x)\le a]$$

• $f(x)$ 可以在 $O(n\ln n)$ 的时间预处理出来

• 因为有 $a$ 的限制，当 $a$ 变大时，后半部分式子也在改变。所以考虑将询问离线，按 $a$ 的值从小到大排序，每次 $a$ 变大时都类似插入操作，考虑用树状数组维护

• 设 $g(T)=\sum _{x|T}\mu (\frac{T}{x})\times f(x)$，对于所有 $f(x)\le a$ 来说，它会对 $g(x),g(2x)\dots$ 产生贡献，所以给每个 $g$ 都加上 $\mu (\frac{T}{x})\times f(x)$ 即可

代码

#include<bits/stdc++.h>
#define LL long long
using namespace std;

const int N=100010,Q=20010;
const LL MOD=2147483648;

struct node
{
	int dat;
	int id;
}f[N];
struct que
{
	int n,m,s,id;
}a[Q];

int T;
int prime[N],v[N],tot,mu[N];
LL c[N],ans[N];

bool cmp(que xx,que yy)
{
	return xx.s<yy.s;
}

bool cmp2(node xx,node yy)
{
	return xx.dat<yy.dat;
}

void add(int x,LL y)
{
	for(; x<=1e5; x+=(x&-x))
		c[x]=(c[x]+y)%MOD;
}

LL ask(int x)
{
	LL s=0;
	for(; x; x-=(x&-x))
		s=(c[x]+s)%MOD;
	return s;
}

void update(int x,LL y)
{
	for(int i=1; i*x<=1e5; i++)
		add(i*x,1LL*mu[i]*y);
}

void primes()
{
	mu[1]=1;
	for(int i=2; i<=1e5; i++)
	{
		if(!v[i])
		{
			v[i]=i;
			prime[++tot]=i;
			mu[i]=-1;
		}
			
		for(int j=1; j<=tot; j++)
		{
			if(prime[j]>v[i] || prime[j]>1e5/i)
				break;
			v[i*prime[j]]=prime[j];
			
			if(i%prime[j]==0)
				break;
			mu[i*prime[j]]=-mu[i];
		}
	}
	
	for(int i=1; i<=1e5; i++)
		for(int j=1; i*j<=1e5; j++)
			f[i*j].dat=(f[i*j].dat+i)%MOD;
		
	for(int i=1; i<=1e5; i++)
		f[i].id=i;
		
	sort(f+1,f+1+100000,cmp2);
}

LL calc(int a,int b)
{
	if(a>b)
		swap(a,b);
		
	LL res=0;
	for(int l=1,r; l<=a; l=r+1)
	{
		r=min(a/(a/l),b/(b/l));
		res=(res+(ask(r)-ask(l-1)+MOD)%MOD*1LL*(a/l)%MOD*(b/l)%MOD)%MOD;
	}
	
	return res;
}
 		
int main()
{
	primes();

	scanf("%d",&T);
	for(int i=1; i<=T; i++)
		scanf("%d%d%d",&a[i].n,&a[i].m,&a[i].s),a[i].id=i;
	
	sort(a+1,a+1+T,cmp);
	
	int p=1;
	for(int i=1; i<=T; i++)
	{
		while(f[p].dat<=a[i].s && p<=1e5)
			update(f[p].id,f[p].dat),p++;
		ans[a[i].id]=calc(a[i].n,a[i].m);
	}
	
	for(int i=1; i<=T; i++)
		printf("%lld\n",ans[i]);
    
    return 0;
}