题目描述
有一张N*m的数表,其第i行第j列(1 < =i < =n,1 < =j < =m)的数值为能同时整除i和j的所有自然数之和。给定a,计算数表中不大于a的数之和。
输入输出格式
输入格式:
输入包含多组数据。
输入的第一行一个整数Q表示测试点内的数据组数
接下来Q行,每行三个整数n,m,a(|a| < =10^9)描述一组数据。
输出格式:
对每组数据,输出一行一个整数,表示答案模2^31的值。
输入输出样例
输入样例#1:
2
4 4 3
10 10 5
输出样例#1:
20
148
说明
1 < =N.m < =10^5 , 1 < =Q < =2*10^4
分析:i,j的公约数可以看作是gcd(i,j)的约数,设为F(gcd(i,j)),其中F(gcd(i,j))<=a。
设
,
则
前面显然可以整除分块,关键是后面了。
设
显然 ,就是一个积性函数,线筛一波,顺便把 也给筛出来。
我们发现 有贡献时, 。
于是我们可以把询问离线,把各个询问的a从小到大排序,也就是说,每次把ai-1到ai的F对g的贡献给算出来加入到g,也就是要个F排序了,那么相当于动态维护小于等于a时的g前缀和,显然这个时候就是树状数组一波了。
代码:
// luogu-judger-enable-o2
#include <iostream>
#include <cstdio>
#include <cmath>
#include <algorithm>
const int maxc=2e4+7;
const int maxn=1e5+7;
using namespace std;
int test,cnt;
int not_prime[maxn],prime[maxn],mul[maxn],low[maxn];
int t[maxn];
struct rec{
int x,y;
}f[maxn];
struct node{
int n,m,k,ans,num;
}q[maxc];
bool cmp1(node x,node y)
{
return x.k<y.k;
}
bool cmp2(rec x,rec y)
{
return x.y<y.y;
}
bool cmp3(node x,node y)
{
return x.num<y.num;
}
void getmul(int n)
{
f[1]=(rec){1,1};
mul[1]=1;
for (int i=2;i<=n;i++)
{
if (!not_prime[i])
{
prime[++cnt]=i;
f[i]=(rec){i,i+1};
low[i]=i;
mul[i]=-1;
}
for (int j=1;j<=cnt;j++)
{
if (i*prime[j]>n) break;
not_prime[i*prime[j]]=1;
if (low[i]%prime[j]==0) low[i*prime[j]]=low[i]*prime[j];
else low[i*prime[j]]=prime[j];
if (i*prime[j]==low[i*prime[j]])
{
f[i*prime[j]]=(rec){i*prime[j],f[i].y+i*prime[j]};
mul[i*prime[j]]=0;
}
else
{
int x=low[i*prime[j]],y=i*prime[j]/x;
f[i*prime[j]]=(rec){i*prime[j],f[x].y*f[y].y};
mul[i*prime[j]]=mul[x]*mul[y];
}
}
}
sort(f+1,f+n+1,cmp2);
}
int lowbit(int x)
{
return x&(-x);
}
void ins(int x,int k)
{
while (x<=1e5)
{
t[x]+=k;
x+=lowbit(x);
}
}
int getsum(int x)
{
int ans=0;
while (x>0)
{
ans+=t[x];
x-=lowbit(x);
}
return ans;
}
int main()
{
scanf("%d",&test);
getmul(1e5);
for (int i=1;i<=test;i++)
{
scanf("%d%d%d",&q[i].n,&q[i].m,&q[i].k);
q[i].num=i;
if (q[i].n>q[i].m) swap(q[i].n,q[i].m);
}
sort(q+1,q+test+1,cmp1);
int d=1;
for (int i=1;i<=test;i++)
{
while ((f[d].y<=q[i].k) && (d<=1e5))
{
for (int j=f[d].x;j<=1e5;j+=f[d].x)
{
ins(j,f[d].y*mul[j/f[d].x]);
}
d++;
}
int x=0,y=0;
int n=q[i].n,m=q[i].m;
for (int j=1,last;j<=n;j=last+1)
{
last=min(n/(n/j),m/(m/j));
y=getsum(last);
q[i].ans+=(n/j)*(m/j)*(y-x);
x=y;
}
}
sort(q+1,q+test+1,cmp3);
for (int i=1;i<=test;i++) printf("%d\n",q[i].ans&0x7fffffff);
}