题目大意:n次询问,每一次询问有多少数对(x,y)满足a<=x<=b,c<=y<=d,gcd(x,y)=k。a,b,c,d,k已经给出。 a,b,c,d,k,n<=50000。
首先容易想到,实际上问题可以转化成求 [1,b]和[1,d],[1,a]和[1,c],[1,a]和[1,d],[1,b]和[1,c]四个区间gcd(x,y)=k的,然后容斥原理直接做即可。0.0
然后gcd(x,y)=k可以转化成gcd(x/k,y/k)=1,在n/k和m/k的范围内。
那么咱们的问题实际上就是转化成了求 1<=x<=m,1<=y<=n 的gcd(x,y)=k的个数即可。做四次然后容斥搞。
这个貌似咱们会做了额...莫比乌斯反演搞一搞就可以了...但是!,那个的复杂度是O(max(a,b,c,d)/k)的!,如果都取极限数据的话50000个询问显然会O(n^2)从而TLE....所以咱们这里需要考虑一下如何优化额 0.0
咱反正是想了半天不知道怎么优化额...【太弱了咩= =】
神犇A君:我有特别的优化技巧!
神犇B君:这道题是原题啊2333!
2333333....
优化的关键在于这里,最后的答案显然都是这样的形式:ans=∑miu[i]*(m/i)*(n/i),实际上这里存在一个冗余...咱们首先证明 m/i(i取遍1~m) 最多只有 O(sqrt(m)) 种取法,其中/为除,取整除法的那种。
一开始没想出来...跑去问叉姐一行神犇= =然后被虐傻了...
其实非常单纯的证明额...当d<=sqrt(m)的时候,最多有sqrt(m)种组合【显然的吧...】,然后当d>sqrt(m)的时候,m/d就在sqrt(m)内部了,所以还是sqrt(m)的...
所以总的就是O(sqrt(m))的了= =。
明明是非常单纯的证法...听完就觉得自己太弱了...
然后由于是连续的区间,所以显然 m/i,n/i所有的相同的连续区间最多有O(sqrt(m)+sqrt(n))种。
然后咱们处理下莫比乌斯函数的前缀和,就可以将a,b按连续区间分块,然后计算的时候按块来,复杂度就是O(sqrt(m)+sqrt(n))
但是---------咱们还是得处理出那个a,b的所有的分块的区间来额...不然复杂度还是不会变。咱们只需要保证得到了某一个区间能够在O(1)求出下一个区间即可,这样子就可以保证复杂度是O(sqrt(m)+sqrt(n))不变。这样子最后的复杂度就是O(n(sqrt(a)+sqrt(b))),对于5*10^4还是可以过的。但是...咱们还是不知道怎么O(1)的球区间...实际上这也好考虑...咱们先来考虑对于单个m,怎么求出连续的区间。假设某个值j的连续的区间的开头是i,即m/i=j,i是所有满足条件的i的第一个,那么咱们就考虑最后一个是多少,它是可以求出来的,就是 m/(m/i),这里的/为整除。至于为什么是这个,咱们可以证明出来,咱们设m/i=j,则 j*i+k=m,当咱们增加i的时候,j*(i+1)+k'=m,则k变成k-j,然后一直到k变成负数之前都是可行的,而达到负数的时候实际上是k'<j的时候,这个时候的i实际上就是 m/j,即 m/(m/i)
所以咱们从1,1开始,每一次挑出min(n/(n/i),m/(m/i)),结合当前的区间就是答案。复杂度为O(N(sqrt(a)+sqrt(b))),这里预处理出miu函数的前缀和,只需要O(n)即可。
另外...这题学习到了又一个东西= =那就是能用longint的地方不要用int64,尤其是在卡时间的时候最重要的= =这道题调了半天常数终于过去了。还有一个优化是BLADEVIL神犇提出来的,咱们可以在中途使用longint计算,然后赋值给int64的时候可能会越界,这个时候咱们使用一种强制转换的函数,即 int64(),括弧里面是一个变量,能将某个变量在运算中强制变成某个类型,这样子的话速度貌似可以更快!
=================================
program bzoj2301;
const maxn=50010;
var p,miu:array[0..maxn]of longint;
v:array[0..maxn]of boolean;
sum:array[0..maxn]of int64;
a,b,c,d,k,n:longint;
ans:int64;
function min(a,b:longint):longint;
begin
if a<b then exit(a) else exit(b);
end;
procedure prepare;
var i,j:longint;
begin
fillchar(v,sizeof(v),true);
miu[1]:=1;p[0]:=0;
for i:=2 to 50000 do
begin
if v[i] then
begin
inc(p[0]);p[p[0]]:=i;
miu[i]:=-1;
end;
for j:=1 to p[0] do
begin
if i*p[j]<=50000 then v[i*p[j]]:=false;
if i*p[j]>50000 then break;
if i mod p[j]=0 then
begin
miu[i*p[j]]:=0;
break;
end
else miu[i*p[j]]:=-miu[i];
end;
end;
fillchar(sum,sizeof(sum),0);
for i:=1 to 50000 do
sum[i]:=miu[i]+sum[i-1];
end;
function calc(n,m:longint):longint;
var tt,i,ans,tmp:longint;
begin
if n>m then
begin
tt:=n;n:=m;m:=tt;
end;
i:=1;ans:=0;
while i<=n do
begin
tmp:=min(n div (n div i),m div (m div i));
ans:=ans+(sum[tmp]-sum[i-1])*(n div i)*(m div i);
i:=tmp+1;
end;
exit(ans);
end;
procedure init;
var i:longint;
begin
prepare;
read(n);
for i:=1 to n do
begin
read(a,b,c,d,k);
a:=(a-1)div k;
b:=b div k;
c:=(c-1)div k;
d:=d div k;
ans:=int64(calc(b,d))+int64(calc(a,c))-int64(calc(b,c))-int64(calc(a,d));
writeln(ans);
end;
end;
begin
init;
end.
评论