【bzoj2820】GCD

Time Limit: 3000 ms Memory Limit: 256 MB

description

​ 神犇GJS虐完数论后给zzHGR出了一个数论题。
​ 给定n,m，求1≤x≤n,1≤y≤m，且gcd(x,y)为质数的(x,y)有多少对。
​ zzHGR必然不会了，于是向你来请教……
​ 多组输入。

input

​ 第一行一个整数T，表述数据组数。
​ 接下来T行，每行两个正整数，表示n,m。

output

​ T行，每行一个整数表示第i组数据的结果

sample input

2
10 10
100 100

sample output

30
2791

HINT

​ (T=10000 ，n,m≤10^7)

---

反演的第一题！！啊哈哈哈哈哈哈哈超级高兴！

咳咳

首先写出表达式

相当于求（这里默认(n<m) ）

[egin{aligned} ans &= sumlimits_{i=1}^{n}sumlimits_{j=1}^{n} [ gcd(i , j) = prime ]\ &=sumlimits_{p=prime}^{}sumlimits_{i=1}^{n/p}sumlimits_{j=1}^{n/p} [ gcd(i , j) = 1 ]\ &=sumlimits_{p=prime}^{}sumlimits_{i=1}^{n/p}sumlimits_{j=1}^{m/p} sumlimits_{d|i, d|j} mu(d)&(sumlimits_{d|gcd(i,j)}mu(d) [gcd(i,j) = 1])\ &=sumlimits_{p=prime}^{}sumlimits_{p|k}^{k<=n} mu(frac{k}{p}) * lfloor frac{n}{k} floor * lfloorfrac{m}{k} floor&(k = pd)\ \ &然后。。观察一下前两个sigma我们先枚举了p然后再枚举了p在范围内的倍数。。\&那。。。不就是枚举了所有范围内的数吗哈！哈！哈！\ \ &=sumlimits_{k<=n} g(k) * lfloor frac{n}{k} floor * lfloorfrac{m}{k} floor& (g(k)=sumlimits_{p=prime}^{p|k}mu(frac{k}{p}))\ end{aligned} ]

费劲千辛万苦搞到这条看起来十分友善的式子，那么剩下的就是……

**考虑(g(k))怎么求 **

((g(k) = sumlimits_{p=prime}^{p|k}mu(frac{k}{p}))，我们设(k = p_0 * x)，其中(p_0)为质数)

那么就有两种情况：

1.如果(p_0|x)

​ 如果说(p_0=p)，那么(mu(frac{k}{p}) = mu(k))

​ 如果说(p_0 e p)，那么(p_0 * x)就有质数平方因子了，而除以(p)又不能将其消掉，所以(mu)值为0

​ 综上就是(g(k) = mu(k))

2.如果(p0)不是(x)的因子

​ 如果说(p_0 = p)，和上面一样

​ 如果说(p_0 e p)，

​ 那么由于(mu)是积性函数，且(p_0)与(frac{x}{p})互质，所以(mu(frac{k}{p}) = mu(frac{p_0 * x}{p}) = mu(p_0) * mu(frac{x}{p}))

​ 然后把(mu(p_0))提前我们就可以发现后面的式子其实就是(g(x))，而(mu(p_0) = -1)

​ 所以此时(mu(frac{k}{p}) = -g(x))

​ 综上就是(g(k) = mu(k) - g(x))，其中(x * p_0 = k)且(p_0)为质数

然后我们就可以十分愉快地将(g(k))筛出来啦

然而是多组数据

接下来再来枚举……稳稳的T啊……

咋办嘞。。显然要让求(f(x))变快啊，然后就发现因为后面两个东西是下取整

下取整。。那就说明有很多数其实弄到最后是一样的

和n/i下取整一样的数最大是n/(n/i)，下一个区间的开始刚好就是最大的那个数+1

然后就考虑到可以用前缀和来搞一下，瞬间舒服ovo

然后就十分愉快地做完啦ovo

---

#include<iostream>
#include<cstdio>
#include<cstring>
#define ll long long
using namespace std;
const int MAXN=1e7+10;
int miu[MAXN],p[MAXN],g[MAXN],sum[MAXN];
bool vis[MAXN];
ll ans;
int n,m,T,cnt,pos;
int prework(int n);

int main(){
#ifndef ONLINE_JUDGE
	freopen("a.in","r",stdin);
#endif
	scanf("%d",&T);
	prework(10000000);
	for (int o=1;o<=T;++o){
		scanf("%d%d",&n,&m);
		if (n>m) swap(n,m);
		ans=0;
		for (int i=1;i<=n;i=pos+1){
			pos=min(n/(n/i),m/(m/i));
			ans+=(ll)(sum[pos]-sum[i-1])*(n/i)*(m/i);
		}
		printf("%lld
",ans);
	}
}

int prework(int n){
	memset(vis,0,sizeof(vis));
	miu[1]=1;
	cnt=0;
	for (int i=2;i<=n;++i){
		if (!vis[i]){
			p[++cnt]=i;
			miu[i]=-1; g[i]=1;
		}
		for (int j=1;j<=cnt&&i*p[j]<=n;++j){
			vis[i*p[j]]=1;
			if (i%p[j]){
				miu[i*p[j]]=-miu[i];
				g[i*p[j]]=miu[i]-g[i];
			}
			else{
				miu[i*p[j]]=0,g[i*p[j]]=miu[i];
				break;
			}
		}
	}
	for (int i=1;i<=n;++i)
		sum[i]=sum[i-1]+g[i];
}