这题做的历程堪称惊心动魄
刚刚学了莫比乌斯反演的我高高兴兴的和cbx一起反演式子
期间有突破,有停滞,有否定
然后苟蒻的我背着cbx偷偷打开了题解
看到了
我。。。。。。
去你的有个性质啊(当然还是自己知识储备不足)
具体证明
(其实当时主要是想的方向偏了,不然这个定理自己也能想出来)
然后就可以愉快的反演了
Σ(i∈[1,n])Σ(j∈[1.m])d(x,y)
=Σ(i=1)Σ(j=1)Σ(x|i)Σ(y|j)[gcd(x,y)==1]
=Σ(i=1)Σ(j=1)((n/i)*(m/j))Σ(d|i&&d|j)μ(d)
=Σ(d=1)μ(d)Σ(i=1) (n/(d*i)) Σ(j=1)(m/(d*j))
然后我们观察Σ(n/(d*i))
根据性质 (n/(d*i))==((n/d)/i)
我们发现这个东西可以用数论分块O(sqrt(n))预处理,设为f[i]
则原式= Σ(d=1)(μ(d)f[n/d]*f[m/d])
再用数论分块就好了
复杂度O(n*sqrt(n)+T*sqrt(n))
1 #include<iostream>
2 #include<cstdio>
3 #include<cmath>
4 #define ll long long
5 using namespace std;
6 int mu[50100],p[50010],top;ll tot[50100],f[50100];bool v[50010];
7 int main(){
8 f[1]=1;tot[1]=1;
9 for(int i=2;i<=50000;i++){
10 if(!v[i]){
11 p[++top]=i;
12 mu[i]=-1;
13 }
14 for(int j=1;j<=top&&i*p[j]<=50000;j++){
15 if(!(i%p[j])){
16 v[i*p[j]]=1;
17 break;
18 }
19 mu[i*p[j]]=-mu[i];
20 v[i*p[j]]=1;
21 }
22 tot[i]=tot[i-1]+mu[i];
23 int x;
24 for(int j=1;j<=i;j=x+1){
25 x=(i/(i/j));
26 f[i]+=(x-j+1)*(i/j);
27 }
28 }
29 int j,n,m,t;ll ans;
30 scanf("%d",&t);
31 while(t--){
32 scanf("%d%d",&n,&m);
33 if(n>m) swap(n,m);ans=0;
34 for(int i=1;i<=n;i=j+1){
35 j=min((n/(n/i)),(m/(m/i)));
36 ans+=(tot[j]-tot[i-1])*f[n/i]*f[m/i];
37 }
38 printf("%lld
",ans);
39 }
40 }