Description
Input
Output
对每组数据,输出一行一个整数,表示答案模2^31的值。
Sample Input
4 4 3
10 10 5
Sample Output
148
题解Here!
莫比乌斯反演的标准模板题。
首先,把那个$a$的约束与取模运算先丢一边去。
设$f(x)$为题目要求的约束和,即$f(x)=sum_{d|x}d$。
先保证$n<=m$。
那么题目要求的就是:$Ans=sum_{i=1}^nsum_{j=1}^mf(gcd(i,j))$。
套路,更换枚举项:
$$Ans=sum_{d=1}^nsum_{i=1}^nsum_{j=1}^mf(d)[gcd(i,j)==d]$$
把那个$f(d)$提到前面来:
$$Ans=sum_{d=1}^nf(d)sum_{i=1}^nsum_{j=1}^m[gcd(i,j)==d]$$
后面这个式子好熟悉啊:$sum_{i=1}^nsum_{j=1}^m[gcd(i,j)==d]$
板子了,不会请看这里:洛谷P3455 [POI2007]ZAP-Queries
于是:
$$Ans=sum_{d=1}^nf(d)sum_{i=1}^{lfloorfrac{n}{d}
floor}mu(i)lfloorfrac{n}{id}
floorlfloorfrac{m}{id}
floor$$
那个$id$很烦人,设$D=id$,则:
$$Ans=sum_{D=1}^nsum_{d|D}f(d)mu(frac{D}{d})lfloorfrac{n}{D}
floorlfloorfrac{m}{D}
floor$$
后面那一堆与$d$无关,提到前面来:
$$Ans=sum_{D=1}^nlfloorfrac{n}{D}
floorlfloorfrac{m}{D}
floorsum_{d|D}f(d)mu(frac{D}{d})$$
很显然嘛,前面的数论分块。
那后面的$f(d)$呢?又不能线性筛。。。
等等,不能线性筛就暴力算啊!
每个数暴力算到它的倍数里去。
复杂度?
首先有个式子:(别问我怎么证明的,我也不知道。。。)
$$frac{n}{1}+frac{n}{2}+frac{n}{3}+...+frac{n}{n}=nlog_2n$$
于是复杂度就是$O(nlog_2n)$的。
再看$a$的限制。
我们可以离线处理,将询问按$a$从小到大排序,$i$按$f(i)$从小到大排序。
每次询问将$f(i)<=a$的$i$插入树状数组,维护前缀和。
那,取模怎么办?
因为模数是$2^{31}$,正好是$int$型的二进制上限。
我们可以令其自然溢出,最后输出的时候按位与$2^{31}-1$即可。
复杂度是$O(nlog_2n+qsqrt nlog_2n)$。
附代码:
#include<iostream> #include<algorithm> #include<cstdio> #define MAXN 100010 using namespace std; int n=0,q; int bit[MAXN],ans[MAXN]; int k=0,prime[MAXN],mu[MAXN],sum[MAXN],pos[MAXN]; bool np[MAXN]; struct Question{ int n,m,a,id; }que[MAXN]; inline int read(){ int date=0,w=1;char c=0; while(c<'0'||c>'9'){if(c=='-')w=-1;c=getchar();} while(c>='0'&&c<='9'){date=date*10+c-'0';c=getchar();} return date*w; } inline bool cmp1(const int &p,const int &q){ return sum[p]<sum[q]; } inline bool cmp2(const Question &p,const Question &q){ return p.a<q.a; } inline int lowbit(int x){return x&(-x);} inline void update(int x,int v){for(;x<=n;x+=lowbit(x))bit[x]+=v;} inline int query(int x){int s=0;for(;x;x-=lowbit(x))s+=bit[x];return s;} void make(){ int m=n; mu[1]=1; for(int i=2;i<=m;i++){ if(!np[i]){ prime[++k]=i; mu[i]=-1; } for(int j=1;j<=k&&prime[j]*i<=m;j++){ np[prime[j]*i]=true; if(i%prime[j]==0)break; mu[prime[j]*i]=-mu[i]; } } for(int i=1;i<=m;i++){ pos[i]=i; for(int j=i;j<=m;j+=i)sum[j]+=i; } sort(pos+1,pos+m+1,cmp1); } int solve(int n,int m){ int ans=0; for(int i=1,last=1;i<=n;i=last+1){ last=min(n/(n/i),m/(m/i)); ans+=(n/i)*(m/i)*(query(last)-query(i-1)); } return ans; } void work(){ for(int i=1,now=1;i<=q;i++){ for(;sum[pos[now]]<=que[i].a&&now<=n;now++) for(int k=pos[now];k<=n;k+=pos[now]) update(k,sum[pos[now]]*mu[k/pos[now]]); ans[que[i].id]=solve(que[i].n,que[i].m); } for(int i=1;i<=q;i++)printf("%d ",ans[i]&2147483647); } void init(){ q=read(); for(int i=1;i<=q;i++){ que[i].n=read();que[i].m=read();que[i].a=read(); if(que[i].n>que[i].m)swap(que[i].n,que[i].m); que[i].a=max(0,que[i].a); que[i].id=i; n=max(n,que[i].m); } sort(que+1,que+q+1,cmp2); } int main(){ init(); make(); work(); return 0; }