Description
有一张N×m的数表,其第i行第j列(1 < =i < =n,1 < =j < =m)的数值为
能同时整除i和j的所有自然数之和。给定a,计算数表中不大于a的数之和。
Input
输入包含多组数据。
输入的第一行一个整数Q表示测试点内的数据组数,接下来Q行,每行三个整数n,m,a(|a| < =10^9)描述一组数据。
Output
对每组数据,输出一行一个整数,表示答案模2^31的值。
Sample Input
2
4 4 3
10 10 5
Sample Output
20
148
HINT
1 < =N.m < =10^5 , 1 < =Q < =2×10^4
这个太卡时间了,搞了我好久,不过最后跑了8s,利用了他取模的数,C++正好可以直接爆int,然后小于0就加2的31次方就行了,我就用longint强行截取了后半部分相当于爆int
首先a[i,j]上的数我们可以看成是F[gcd(i,j)],F[i]我们都预处理出来
然后我们要求的就是 ΣF[i]*g[i] (F[i]<=a,g[i]是gcd=i的个数)
所以我们要求的就是 ΣF[i]*Σtrunc(n/d)*trunc(m/d)*μ(d/i) (F[i]<=a,i|d)
所以我们要求的就是 Σtrunc(n/d)*trunc(m/d)*ΣF[i]*μ(d/i) (F[i]<=a,i|d)
因为trunc(n/d)*trunc(m/d)只有根号n级别的个数,所以我们要处理ΣF[i]*μ(d/i)的前缀和,我用的是树状数组
然后每次询问都可以根号n*logn回答了
但是还有一个条件,就是F[i]<=a
前缀和还是动态的,怎么办
首先把F[i]排序是肯定的,对于每一个i,影响到的是他的倍数,这个很麻烦啊
于是就离线做了,把a排序,然后暴力修改前缀和,即枚举倍数
因为排序了,所以我们最多每个i都枚举倍数一次,应该是nlognlogn的
然后就写完了
const
maxn=;
h=<<;
type
node=record
n,m,a,id:longint;
end;
var
q:array[..maxn]of node;
flag:array[..maxn]of boolean;
p,u,k,f,ans,c:array[..maxn]of longint;
t,tot,max:longint; function min(x,y:longint):longint;
begin
if x<y then exit(x);
exit(y);
end; procedure swap(var x,y:longint);
var
t:longint;
begin
t:=x;x:=y;y:=t;
end; procedure swap(var x,y:node);
var
t:node;
begin
t:=x;x:=y;y:=t;
end; procedure sort(l,r:longint);
var
i,j:longint;
y:int64;
begin
i:=l;
j:=r;
y:=f[k[(l+r)>>]];
repeat
while f[k[i]]<y do
inc(i);
while f[k[j]]>y do
dec(j);
if i<=j then
begin
swap(k[i],k[j]);
inc(i);
dec(j);
end;
until i>j;
if i<r then sort(i,r);
if j>l then sort(l,j);
end; procedure pre;
var
i,j,s:longint;
begin
f[]:=;
u[]:=;
for i:= to max do
begin
if flag[i]=false then
begin
inc(tot);
p[tot]:=i;
f[i]:=i+;
u[i]:=-;
end;
for j:= to tot do
begin
if int64(i)*p[j]>max then break;
flag[i*p[j]]:=true;
if i mod p[j]= then
begin
s:=p[j];
while i mod (int64(s)*p[j])= do
s:=s*p[j];
if s=i then f[i*p[j]]:=(s*p[j]*p[j]-)div(p[j]-)
else f[i*p[j]]:=f[i div s]*f[s*p[j]];
break;
end
else
begin
f[i*p[j]]:=f[i]*(p[j]+);
u[i*p[j]]:=-u[i];
end;
end;
end;
for i:= to max do
k[i]:=i;
sort(,max);
end; procedure sort2(l,r:longint);
var
i,j,y:longint;
begin
i:=l;
j:=r;
y:=q[(l+r)>>].a;
repeat
while q[i].a<y do
inc(i);
while q[j].a>y do
dec(j);
if i<=j then
begin
swap(q[i],q[j]);
inc(i);
dec(j);
end;
until i>j;
if i<r then sort2(i,r);
if j>l then sort2(l,j);
end; procedure init;
var
i:longint;
begin
read(t);
for i:= to t do
begin
read(q[i].n,q[i].m,q[i].a);
q[i].id:=i;
if q[i].n>q[i].m then swap(q[i].n,q[i].m);
if max<q[i].n then max:=q[i].n;
end;
sort2(,t);
end; procedure add(x,y:longint);
begin
while x<= do
begin
c[x]:=longint(int64(c[x])+y);
x:=x+(x and (-x));
end;
end; function sum(x:longint):longint;
begin
sum:=;
while x> do
begin
sum:=longint(int64(sum)+c[x]);
x:=x-(x and (-x));
end;
end; procedure main;
var
i,j,kk,s,lasta,s1,s2:longint;
begin
lasta:=;
for i:= to t do
begin
while (lasta<=max) and (f[k[lasta]]<=q[i].a) do
begin
s:=k[lasta];
j:=;
while s<=max do
begin
if u[j]<> then add(s,f[k[lasta]]*u[j]);
inc(s,k[lasta]);
inc(j);
end;
inc(lasta);
end;
kk:=;
while kk<=q[i].n do
begin
s1:=q[i].n div kk;
s2:=q[i].m div kk;
s:=min(trunc(q[i].n/s1),trunc(q[i].m/s2));
ans[q[i].id]:=longint(int64(ans[q[i].id])+int64(longint(int64(s1)*s2))*longint(int64(sum(s))-sum(kk-)));
kk:=s+;
end;
end;
for i:= to t do
if ans[i]< then writeln(ans[i]+h)
else writeln(ans[i]);
end; begin
init;
pre;
main;
end.