BZOJ 3529: [Sdoi2014]数表
3529: [Sdoi2014]数表
Time Limit: 10 Sec Memory Limit: 512 MBSubmit: 1407 Solved: 699
[Submit][Status][Discuss]
Description
有一张N×m的数表,其第i行第j列(1 < =i < =礼,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
4 4 3
10 10 5
Sample Output
20
148
148
HINT
1 < =N.m < =10^5 , 1 < =Q < =2×10^4
Source
分析:
论我的心情变化:劳资不陪你玩了,我还不信推不出来了--->TAT你赢了--->%PoPoQQQ...
大家都膜拜了Po姐,那我们也来膜拜一发...
首先公式恐惧症的我看到公式就GG,所以先把a的限制去掉...
定义F(i)代表i的约数之和...g(i)代表gcd(x,y)=i的数对个数...根据BZOJ1101我们可以得出g(i)的做法...
Σ(1<=i<=n) Σ(1<=j<=n) F(gcd(i,j)) (还是规定n<m)
=Σ(1<=i<=n) F(i)*g(i)
=Σ(1<=i<=n) F(i)*Σ(1<=x<=n/i) μ(x)*(n/i/x)*(m/i/x)
=Σ(1<=i<=n) F(i)*Σ(i|x&&x<=n) μ(x/i)*(n/x)*(m/x)
=Σ(1<=x<=n) (n/x)*(m/x)*Σ(i|x)μ(x/i)*F(i)
所以我们还是可以采用分段的思想...
先暴力的预处理出F函数和μ函数,然后枚举i,再枚举i的倍数x,用树状数组维护前缀和...
现在加上了a的限制...那就离线,按照a从小到大排序...
代码:
1 #include<algorithm> 2 #include<iostream> 3 #include<cstring> 4 #include<cstdio> 5 //by NeighThorn 6 using namespace std; 7 //大鹏一日同风起,扶摇直上九万里 8 9 const int maxn=100000+5,maxm=20000+5; 10 11 int cas,cnt,Max,Mod=2147483647,tr[maxn],miu[maxn],vis[maxn],prime[maxn]; 12 13 struct M{ 14 int n,m,a,id; 15 unsigned int ans; 16 }q[maxm]; 17 18 struct N{ 19 int pos,val; 20 }f[maxn]; 21 22 inline bool cmp1(M x,M y){ 23 return x.a<y.a; 24 } 25 26 inline bool cmp2(M x,M y){ 27 return x.id<y.id; 28 } 29 30 inline bool cmp3(N x,N y){ 31 return x.val<y.val; 32 } 33 34 inline void insert(int x,int y){ 35 for(;x<=Max;x+=x&(-x)) 36 tr[x]+=y; 37 } 38 39 inline int query(int x){ 40 int sum=0; 41 for(;x;x-=x&(-x)) 42 sum+=tr[x]; 43 return sum; 44 } 45 46 signed main(void){ 47 scanf("%d",&cas);Max=0; 48 for(int i=1;i<=cas;i++){ 49 scanf("%d%d%d",&q[i].n,&q[i].m,&q[i].a),q[i].id=i; 50 if(q[i].n>q[i].m) 51 swap(q[i].n,q[i].m); 52 Max=max(Max,q[i].m); 53 }miu[1]=1;cnt=0; 54 memset(vis,0,sizeof(vis)); 55 for(int i=2;i<=Max;i++){ 56 if(!vis[i]) 57 vis[i]=1,prime[++cnt]=i,miu[i]=-1; 58 for(int j=1;j<=cnt&&prime[j]*i<=Max;j++){ 59 vis[i*prime[j]]=1; 60 if(i%prime[j]==0){ 61 miu[i*prime[j]]=0;break; 62 } 63 miu[i*prime[j]]=-miu[i]; 64 } 65 } 66 for(int i=1;i<=Max;i++){ 67 f[i].pos=i; 68 for(int j=i;j<=Max;j+=i) 69 f[j].val+=i; 70 } 71 sort(f+1,f+Max+1,cmp3); 72 sort(q+1,q+cas+1,cmp1); 73 int lala=0; 74 for(int i=1;i<=cas;i++){ 75 while(lala+1<=Max&&f[lala+1].val<=q[i].a){ 76 lala++; 77 for(int j=f[lala].pos;j<=Max;j+=f[lala].pos) 78 insert(j,f[lala].val*miu[j/f[lala].pos]); 79 } 80 for(int j=1,r;j<=q[i].n;j=r+1){ 81 r=min(q[i].n/(q[i].n/j),q[i].m/(q[i].m/j)); 82 q[i].ans+=(query(r)-query(j-1))*(q[i].n/j)*(q[i].m/j); 83 } 84 q[i].ans&=Mod; 85 } 86 sort(q+1,q+cas+1,cmp2); 87 for(int i=1;i<=cas;i++) 88 printf("%u\n",q[i].ans); 89 return 0; 90 }
By NeighThorn