[bzoj3529]数表

Posted 2021-03-25 pywbktda

记f(n)表示n的约数和，先不考虑a的限制，那么即求$sum_{i=1}^{n}sum_{j=1}^{m}f(gcd(i,j))$
枚举$d=gcd(i,j)$，即$sum_{d=1}^{n}f(d)sum_{g|d}(n/gd)(m/gd)mu(g)$（后面就是指公约数为d的数对个数）
令$t=gd$，即$sum_{t=1}^{n}(n/t)(m/t)sum_{d|t}f(d)cdot mu(t/d)$，对后半部分预处理（先线性筛出f和$mu$，然后$o(nln(n))$求）+前半部分数论分块即可快速求出
离线，将询问按照a排序，然后不断插入f(i)（即修改所有是i约数的t位置，暴力），用树状数组来支持数论分块即可

 1 #include<bits/stdc++.h>
 2 using namespace std;
 3 #define N 100005
 4 #define L (k<<1)
 5 #define R (L+1)
 6 #define mid (l+r>>1)
 7 #define ui unsigned int
 8 struct ji{
 9     int n,m,a,id;
10     bool operator < (const ji &k)const{
11         return a<k.a;
12     }
13 }q[N];
14 int t,n,m,a,vis[N],p[N],mu[N],id[N];
15 ui f[N],ans[N],tr[N<<2];
16 bool cmp(int x,int y){
17     return f[x]<f[y];
18 }
19 void update(int k,int l,int r,int x,int y){
20     if (l==r){
21         tr[k]+=y;
22         return;
23     }
24     if (x<=mid)update(L,l,mid,x,y);
25     else update(R,mid+1,r,x,y);
26     tr[k]=tr[L]+tr[R];
27 }
28 ui query(int k,int l,int r,int x,int y){
29     if ((l>y)||(x>r))return 0;
30     if ((x<=l)&&(r<=y))return tr[k];
31     return query(L,l,mid,x,y)+query(R,mid+1,r,x,y);
32 }
33 int main(){
34     mu[1]=f[1]=1;
35     for(int i=2;i<N-4;i++){
36         if (!vis[i]){
37             p[++p[0]]=i;
38             mu[i]=-1;
39             f[i]=i+1;
40         }
41         for(int j=1;(j<=p[0])&&(i*p[j]<N-4);j++){
42             vis[i*p[j]]=1;
43             for(int k=i;;k/=p[j])
44                 if (k%p[j]){
45                     if (k>1)f[i*p[j]]=f[k]*f[i*p[j]/k];
46                     else f[i*p[j]]=f[i]+i*p[j];
47                     break;
48                 }
49             if (i%p[j]==0)break;
50             mu[i*p[j]]=-mu[i];
51         }
52     }
53     for(int i=1;i<N-4;i++)id[i]=i;
54     sort(id+1,id+N-4,cmp);
55     scanf("%d",&t);
56     for(int i=1;i<=t;i++){
57         scanf("%d%d%d",&n,&m,&a);
58         if (n>m)swap(n,m); 
59         q[i]=ji{n,m,a,i};
60     }
61     sort(q+1,q+t+1);
62     for(int i=1,j=1;i<=t;i++){
63         for(;(j<N-4)&&(f[id[j]]<=q[i].a);j++)
64             for(int k=1;k<=(N-5)/id[j];k++)update(1,1,N-5,k*id[j],f[id[j]]*mu[k]);
65         for(int l=1,r;l<=q[i].n;l=r+1){
66             r=min(q[i].n/(q[i].n/l),q[i].m/(q[i].m/l));
67             ans[q[i].id]+=(q[i].n/l)*(q[i].m/l)*query(1,1,N-5,l,r);
68         }
69     }
70     for(int i=1;i<=t;i++)printf("%u
",ans[i]&2147483647);
71 }

View Code