「Luogu2522」[HAOI2011]Problem b

problem

Solution

题目要求

x=aby=cd[gcd(x,y)=k]\sum_{x=a}^b\sum_{y=c}^d[gcd(x,y)=k]

设上式为Ans(a,b,c,d,k)Ans(a,b,c,d,k)

不妨来想一想a=1,c=1a=1,c=1的时候怎么做

其实就跟「Luogu3455」[POI2007]ZAP-Queries一模一样了

略过一大堆式子,当a=c=1a=c=1

Ans(1,b,1,d,k)=t=1min(b,d)kμ(t)bktdktAns(1,b,1,d,k)=\sum_{t=1}^{\frac{min(b,d)}{k}}\mu(t)\lfloor\frac{b}{kt}\rfloor\lfloor\frac{d}{kt}\rfloor

容斥一下可得

Ans(a,b,c,d,k)=Ans(1,b,1,d,k)Ans(1,a1,1,d,k)Ans(1,b,1,c1,k)+Ans(1,a1,1,c1,k)Ans(a,b,c,d,k)=Ans(1,b,1,d,k)-Ans(1,a-1,1,d,k)\\-Ans(1,b,1,c-1,k)+Ans(1,a-1,1,c-1,k)

愉快地调用就可以了

Code

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
#include <cstdio>
#include <iostream>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <algorithm>
#define maxn 50005
#define N 50000
using namespace std;
typedef long long ll;

template <typename T> void read(T &t)
{
t=0;int f=0;char c=getchar();
while(!isdigit(c)){f|=c=='-';c=getchar();}
while(isdigit(c)){t=t*10+c-'0';c=getchar();}
if(f)t=-t;
}

int n;
int pri[maxn],pcnt,nop[maxn],mu[maxn];
int k;

void GetPrime()
{
mu[1]=1,nop[1]=1;
for(register int i=2;i<=N;++i)
{
if(!nop[i])pri[++pcnt]=i,mu[i]=-1;
for(register int j=1;j<=pcnt && i*pri[j]<=N;++j)
{
nop[i*pri[j]]=1;
if(i%pri[j]==0)break;
else mu[i*pri[j]]=-mu[i];
}
}
for(register int i=1;i<=N;++i)mu[i]+=mu[i-1];
}

int Calc(int a,int b)
{
int re=0,up=min(a,b)/k;
for(register int l=1,r;l<=up;l=r+1)
{
r=min(a/(a/l),b/(b/l));
re+=(mu[r]-mu[l-1])*(a/(l*k))*(b/(l*k));
}
return re;
}

int main()
{
GetPrime();
read(n);
while(n--)
{
int a,b,c,d;
read(a),read(b),read(c),read(d),read(k);
printf("%d\n",Calc(b,d)-Calc(b,c-1)-Calc(a-1,d)+Calc(a-1,c-1));
}
return 0;
}