ZigZagK的博客
[莫比乌斯函数+线性筛+离线+除法分块+调和级数]BZOJ3529(Sdoi2014)【数表】题解

题目概述

有一张 $n\times m$ 的数表,其第 $i$ 行第 $j$ 列的数值为能同时整除 $i$ 和 $j$ 的所有自然数之和。给定 $a$ , 计算数表中不大于 $a$ 的数之和。

解题报告

我太菜了根本不会做。如果没有限制的话,就是很斯波的莫比乌斯函数题:

$$ \sum_{i=1}^{n}\sum_{j=1}^{m}\sigma[(i,j)]=\sum_{d=1}^{min\{n,m\}}\sum_{i=1}^{n}\sum_{j=1}^{m}\sigma(d)[(i,j)=d]\\ =\sum_{d=1}^{min\{n,m\}}\sigma(d)\sum_{i=1}^{\lfloor{n\over d}\rfloor}\sum_{j=1}^{\lfloor{m\over d}\rfloor}e[(i,j)]=\sum_{d=1}^{min\{n,m\}}\sigma(d)\sum_{k=1}^{min\{\lfloor{n\over d}\rfloor,\lfloor{m\over d}\rfloor\}}\mu(k)\lfloor{n\over dk}\rfloor\lfloor{m\over dk}\rfloor\\ =\sum_{T=1}^{min\{n,m\}}\lfloor{n\over T}\rfloor\lfloor{m\over T}\rfloor\sum_{d|T}\sigma(d)\mu({T\over d}) $$

令 $sum(n)=\sum_{T=1}^{n}\sum_{d|T}\sigma(d)\mu({T\over d})$ ,筛出 $sum(n)$ 之后就可以用除法分块做。但是有限制,我们要想办法去掉限制:离线。

询问按照 $a$ 从大到小排序,$d$ 根据 $\sigma(d)$ 从大到小排序,用树状数组维护 $sum$ ,这样就不用管限制了。

示例程序

#include<cstdio>
#include<cctype>
#include<algorithm>
using namespace std;
typedef long long LL;
const int maxn=100000,maxq=20000;const LL MOD=(LL)1<<31;

int Q,ID[maxn+5],S[maxn+5],u[maxn+5],p[maxn+5];bool pri[maxn+5];
struct data {int n,m,A,ID;} q[maxq+5];LL c[maxn+5],ans[maxq+5];

#define Eoln(x) ((x)==10||(x)==13||(x)==EOF)
inline char readc(){
    static char buf[100000],*l=buf,*r=buf;
    return (l==r&&(r=(l=buf)+fread(buf,1,100000,stdin),l==r))?EOF:*l++;
}
template<typename T> inline int readi(T &x){
    T tot=0;char ch=readc(),lst='+';
    while (!isdigit(ch)) {if (ch==EOF) return EOF;lst=ch;ch=readc();}
    while (isdigit(ch)) tot=(tot<<3)+(tot<<1)+(ch^48),ch=readc();
    return lst=='-'?x=-tot:x=tot,Eoln(ch);
}
inline bool Scmp(const int &a,const int &b) {return S[a]<S[b];}
inline void Make(int n){
    u[1]=1;S[1]=1;ID[1]=1;
    for (int i=2;i<=n;ID[i]=i,i++){
        if (!pri[i]) p[++p[0]]=i,u[i]=-1,S[i]=i+1;
        for (int j=1,t;j<=p[0]&&(t=i*p[j])<=maxn;j++)
            if (i%p[j]) u[t]=-u[i],S[t]=S[i]*S[p[j]],pri[t]=true;
            else {u[t]=0;S[t]=S[i]*S[p[j]]-S[i/p[j]]*p[j];pri[t]=true;break;}
    }sort(ID+1,ID+1+n,Scmp);
}
inline bool Acmp(const data &a,const data &b) {return a.A<b.A;}
inline void AMOD(LL &x,LL y) {if ((x+=y)>=MOD) x-=MOD;}
inline void Update(int x,LL y) {for (int i=x;i<=maxn;i+=i&-i) AMOD(c[i],y);}
inline LL Sum(int x) {LL sum=0;for (int i=x;i;i-=i&-i) AMOD(sum,c[i]);return sum;}
int main(){
    freopen("program.in","r",stdin);freopen("program.out","w",stdout);
    readi(Q);for (int i=1;i<=Q;i++) readi(q[i].n),readi(q[i].m),readi(q[i].A),q[i].ID=i;
    sort(q+1,q+1+Q,Acmp);Make(maxn);
    for (int t=1,j=1;t<=Q;t++){
        int n=q[t].n,m=q[t].m,A=q[t].A,i=q[t].ID;
        for (int d=ID[j];j<=maxn&&S[d]<=A;d=ID[++j])
            for (int T=d;T<=maxn;T+=d) Update(T,(MOD+S[d]*u[T/d])%MOD);
        for (int l=1,r;l<=n&&l<=m;l=r+1){
            r=min(n/(n/l),m/(m/l));LL sum=MOD+Sum(r)-Sum(l-1);
            AMOD(ans[i],sum*(n/l)%MOD*(m/l)%MOD);
        }
    }for (int i=1;i<=Q;i++) printf("%lld\n",ans[i]);return 0;
}
版权声明:本博客所有文章除特别声明外,均采用 CC BY 4.0 CN协议 许可协议。转载请注明出处!