【Ynoi2018】GOSICK

Author Avatar
mrsrz 2月 04, 2020

二次离线莫队。

题目大意：

给定一个序列 $a_1,a_2,a_3,\ldots,a_n$，每次询问给出 $l,r$，查询有多少对二元组 $(i,j)$ 满足 $l\leq i,j\leq r$ 且 $a_i$ 是 $a_j$ 的倍数。

题解：

第十四分块。

这里只讲如何从 $[1,i]$ 递推到 $[1,i+1]$，关于二次离线莫队的具体写法请参见第十四分块(前体)的题解。

考虑用 $b_i$ 表示 $i$ 的因数的出现次数，$c_i$ 表示 $i$ 的倍数的出现次数。那么一个数 $x$ 当前的贡献就为 $b_x+c_x$。

考虑新加入 $x$ 这个数，所有 $x$ 的因数 $y$ 对应的 $c_y$ 就要加 $1$，所有 $x$ 的倍数 $z$ 对应的 $b_z$ 也要加 $1$。

$5\times 10^5$ 内的数的因数个数不多，单次枚举因数部分的复杂度可以看做 $O(\sqrt n)$。

但是当 $x$ 比较小的时候，枚举 $x$ 的倍数就是 $O(n)$ 的。

考虑根号分治，对 $x>\sqrt n$ 的部分直接暴力枚举，单次复杂度 $O(\sqrt n)$，对于 $x\leq \sqrt n$ 的部分，使用另外的方法。

对于 $x\leq\sqrt n$ 的部分，我们考虑解决二次离线莫队中两个需要维护的东西的计算。

  • 计算 $[1,i-1]$ 对 $i$ 的贡献，即计算 $a_ i$ 的倍数在 $[1,i-1]$ 中的出现次数。

    $a_j$ 是 $a_i$（$j\lt i$）的倍数，相当于 $a_i$ 是 $a_j$ 的因数，所以我们对每个数统计它作为因数的出现次数即可。

  • 计算 $[1,i]$ 对 $[l,r]$ 的贡献。我们还是要解决询问单点的倍数/因数个数的问题。但这里总询问次数为 $O (n\sqrt n)$，所以我们必须 $O(1)$ 查询。

    现在的问题在于，我们无法统计 $x$ 的因数中，小于等于 $\sqrt n$ 的因数的出现次数。对于大于等于 $\sqrt n$ 的因数，直接枚举其倍数并加上贡献的复杂度是正确的，但小于等于时是不正确的。

    我们反过来考虑。我们对每个因数 $x$，如果知道了 $l,r$ 内 $x$ 的倍数的出现次数，那么这个 $x$ 作为因数的贡献就得到了。我们对每个 $\leq \sqrt n$ 的因数 $x$ 都计算贡献即可。

    考虑如何快速计算区间 $[l,r]$ 内 $x$ 的倍数的出现次数。我们对每个 $x$，记录 $pre_{x,i}$ 表示前 $i$ 个数中，$x$ 的倍数的出现次数。然后在计算 $[1,i]$ 对 $[l,r]$ 的贡献的时候，我们记录 $[1,i]$ 中每个数的出现次数，然后 $pre_{x,r}-pre_{x,l-1}$ 乘上 $x$ 在 $[1,i]$ 的出现次数，就是 $[l,r]$ 中的数作为倍数，$[1,i]$ 中的数作为因数（不超过 $\sqrt n$） 时的贡献。

这样各部分维护的复杂度都为 $O(n\sqrt n)$ 或 $O(m\sqrt n)$。

时间复杂度 $O((n+m)\sqrt n)$。

关于空间复杂度，直接对 $O(\sqrt n)$ 个因数开前缀和数组存是 $O(n\sqrt n)$ 的。可以只用一个前缀和数组，对每个因数分别重新计算，可以做到 $O(n)$。

卡常数部分：

  1. 线性空间的做法会被卡 Cache（Update：目前不会被卡了，但是根号空间会被卡，下面的代码仅供参考），所以要写空间 $O(n\sqrt n)$ 的做法。由于空间限制，根号分治的阈值只能开到 $100\sim 200$。
  2. 每次都 $O(\sqrt n)$ 计算一个数的因数效率非常低，可以开 vector 存每个数的因数，一开始预处理每个数的因数，然后就可以直接遍历 vector。时空复杂度 $O(n\log n)$。
  3. 存储二次离线的询问时，STL 的 vector 较慢，需要手动分配连续内存。
  4. IO 优化。

Code：

#include<algorithm>
#include<cstring>
#include<cstdio>
#include<cmath>
#include<cctype>
#include<vector>
#define memset0(x) memset(x,0,sizeof(x))
using namespace std;
typedef long long LL;
const int siz=144,N=5e5+3;
int B;
#define bel(x)((x-1)/B)
struct node{
    int l,r,id;
    inline int operator<(const node&rhs)const{return(bel(l)!=bel(rhs.l)?l<rhs.l:((bel(l)&1)?r>rhs.r:r<rhs.r));}
}q[N],Memory_of_winter[N*4],*it=Memory_of_winter;
#define gc pa==pb&&(pb=(pa=buf)+fread(buf,1,500000,stdin),pa==pb)?EOF:*pa++
static char buf[500000],*pa(buf),*pb(buf);
char bwf[9000005],*ww=bwf;
int xx[24];
inline void print(LL x){
    static int*q=xx;
    if(!x)*ww++='0';else{
        for(;x;x/=10)*q++=(x%10)^'0';
        while(q!=xx)*ww++=*--q;
    }
    *ww++='\n';
}
inline int readint() {
    int x=0;char c=gc;
    while(c<'0'||c>'9')c=gc;
    for(;c>='0'&&c<='9';c=gc)x=x*10+(c&15);
    return x;
}
struct Vector{
    node*begin,*end;
    inline void malloc(int len){
        begin=end=it;
        it+=len;
    }
    inline void push_back(const node&x){*end++=x;}
}L[N],R[N];
int totL[N],totR[N];
int n,m,a[N];
int ys[N],bs[N],cnt[N];
LL ans[N],L_R[N],R_L[N];
int pre[N][siz+2];
vector<int>fc[N];
int main(){
    n=readint(),m=readint();
    B=n/(sqrt(m*2/3)+1)+1;
    for(int i=1;i<=n;++i)
    for(int j=i;j<=n;j+=i)
    fc[j].push_back(i);
    for(int i=1;i<=n;++i){
        const int x=a[i]=readint();
        for(int j:fc[x])
        if(j<=siz)++pre[i][j];else break;
        for(int j=1;j<=siz;++j)pre[i][j]+=pre[i-1][j];
    }
    for(int i=1;i<=m;++i)
    q[i].l=readint(),q[q[i].id=i].r=readint();
    for(int i=1;i<=n;++i){
        const int x=a[i];
        L_R[i]=L_R[i-1]+bs[x];
        for(int j:fc[x])++bs[j],L_R[i]+=cnt[j];
        ++cnt[x];
    }
    memset0(cnt),memset0(bs);
    for(int i=n;i;--i){
        const int x=a[i];
        R_L[i]=R_L[i+1]+bs[x];
        for(int j:fc[x])++bs[j],R_L[i]+=cnt[j];
        ++cnt[x];
    }
    memset0(cnt),memset0(bs);
    sort(q+1,q+m+1);
    q[0].l=1,q[0].r=0;
    for(int i=1;i<=m;++i){
        const node&now=q[i],&pre=q[i-1];
        ans[i]=L_R[now.r]-L_R[pre.r]+R_L[now.l]-R_L[pre.l];
        if(now.r>pre.r)++totR[pre.l-1];
        if(now.r<pre.r)++totR[pre.l-1];
        if(now.l<pre.l)++totL[now.r+1];
        if(now.l>pre.l)++totL[now.r+1];
    }
    for(int i=0;i<=n+1;++i)
    L[i].malloc(totL[i]),R[i].malloc(totR[i]);
    for(int i=1;i<=m;++i){
        const node&now=q[i],&pre=q[i-1];
        if(now.r>pre.r)
        R[pre.l-1].push_back((node){pre.r+1,now.r,-i});
        if(now.r<pre.r)
        R[pre.l-1].push_back((node){now.r+1,pre.r,i});
        if(now.l<pre.l)
        L[now.r+1].push_back((node){now.l,pre.l-1,-i});
        if(now.l>pre.l)
        L[now.r+1].push_back((node){pre.l,now.l-1,i});
    }
    for(int i=1;i<=n;++i){
        const int x=a[i];
        for(int j:fc[x])++bs[j];
        if(x>siz)
        for(int j=x;j<=n;j+=x)++bs[j];else ++cnt[x];
        for(node*vc=R[i].begin;vc!=R[i].end;++vc)
        if(vc->id>0){
            LL&t=ans[vc->id];
            for(int k=vc->l;k<=vc->r;++k)t+=bs[a[k]];
            const int*R=pre[vc->r],*L=pre[vc->l-1];
            for(int k=1;k<=siz;++k)t+=(LL)cnt[k]*(R[k]-L[k]);
        }else{
            LL&t=ans[-vc->id];
            for(int k=vc->l;k<=vc->r;++k)t-=bs[a[k]];
            const int*R=pre[vc->r],*L=pre[vc->l-1];
            for(int k=1;k<=siz;++k)t-=(LL)cnt[k]*(R[k]-L[k]);
        }
    }
    memset0(cnt),memset0(bs);
    for(int i=n;i;--i){
        const int x=a[i];
        for(int j:fc[x])++bs[j];
        if(x>siz)
        for(int j=x;j<=n;j+=x)++bs[j];else ++cnt[x];
        for(node*vc=L[i].begin;vc!=L[i].end;++vc)
        if(vc->id>0){
            LL&t=ans[vc->id];
            for(int k=vc->l;k<=vc->r;++k)t+=bs[a[k]];
            const int*R=pre[vc->r],*L=pre[vc->l-1];
            for(int k=1;k<=siz;++k)t+=(LL)cnt[k]*(R[k]-L[k]);
        }else{
            LL&t=ans[-vc->id];
            for(int k=vc->l;k<=vc->r;++k)t-=bs[a[k]];
            const int*R=pre[vc->r],*L=pre[vc->l-1];
            for(int k=1;k<=siz;++k)t-=(LL)cnt[k]*(R[k]-L[k]);
        }
    }
    static LL out[N];
    for(int i=1;i<=m;++i)ans[i]+=ans[i-1],out[q[i].id]=ans[i]+q[i].r-q[i].l+1;
    for(int i=1;i<=m;++i)print(out[i]);
    fwrite(bwf,ww-bwf,1,stdout);
    return 0;
}

如果你想卡空间，那么需要用到一些技巧（Update：这些东西现在不加上也能通过，仅供参考）。

线性空间（除去 $O(n\log n)$ 存储因数的部分）的做法，我们需要对每个 $x$ 都重新求前缀和数组，并对每个询问计算贡献。由于内存访问非常随机，卡不进 Cache，导致常数无法接受。

所以我们的目标是卡进 Cache。

一个非常 naive 的思路就是，对 $\sqrt n$ 这个阈值再分块，即我们开一个 $n\sqrt[4]n$ 的二维数组，存储连续的 $\sqrt[4]n$ 个数作为因数的信息。这样我们在减小空间的基础上，又尽可能地卡进 Cache。

理论空间复杂度是 $O(n\sqrt[4]n)$。实际上这个阈值没法取太大（因为大于阈值的部分只是一个纯粹的循环，常数非常小，而小于等于阈值的部分则非常吃常数），所以实际使用的空间也比较可观。

Update：简单地卡了一下空间，目前能进 $128\texttt{MB}$。

Code：

#include<algorithm>
#include<cstring>
#include<cstdio>
#include<cmath>
#include<cctype>
#include<vector>
#define memset0(x) memset(x,0,sizeof(x))
using namespace std;
typedef long long LL;
const int siz=99,N=5e5+3;
int B;
#define bel(x)((x-1)/B)
struct node{
    int l,r,id;
    inline int operator<(const node&rhs)const{return(bel(l)!=bel(rhs.l)?l<rhs.l:((bel(l)&1)?r>rhs.r:r<rhs.r));}
}q[N],Memory_of_winter[N*4],*it=Memory_of_winter;
#define gc pa==pb&&(pb=(pa=buf)+fread(buf,1,500000,stdin),pa==pb)?EOF:*pa++
static char buf[500000],*pa(buf),*pb(buf);
char bwf[9000005],*ww=bwf;
int xx[24];
inline void print(LL x){
    static int*q=xx;
    if(!x)*ww++='0';else{
        for(;x;x/=10)*q++=(x%10)^'0';
        while(q!=xx)*ww++=*--q;
    }
    *ww++='\n';
}
inline int readint() {
    int x=0;char c=gc;
    while(c<'0'||c>'9')c=gc;
    for(;c>='0'&&c<='9';c=gc)x=x*10+(c&15);
    return x;
}
struct Vector{
    node*begin,*end;
    inline void malloc(int len){
        begin=end=it;
        it+=len;
    }
    inline void push_back(const node&x){*end++=x;}
}L[N],R[N];
int n,m,a[N];
int bs[N],cnt[N];
LL ans[N],L_R[N],R_L[N];
int pre[N][9],ttt[N];
vector<int>fc[N];
int main(){
    n=readint(),m=readint();
    B=n/(sqrt(m*2/3)+1)+1;
    for(int i=1;i<=n;++i)
    for(int j=i;j<=n;j+=i)++ttt[j];
    for(int i=1;i<=n;++i)fc[i].reserve(ttt[i]);
    for(int i=1;i<=n;++i)
    for(int j=i;j<=n;j+=i)
    fc[j].push_back(i);
    for(int i=1;i<=n;++i)a[i]=readint();
    for(int i=1;i<=m;++i)
    q[i].l=readint(),q[q[i].id=i].r=readint();
    for(int i=1;i<=n;++i){
        const int x=a[i];
        L_R[i]=L_R[i-1]+bs[x];
        for(int j:fc[x])++bs[j],L_R[i]+=cnt[j];
        ++cnt[x];
    }
    memset0(cnt),memset0(bs);
    for(int i=n;i;--i){
        const int x=a[i];
        R_L[i]=R_L[i+1]+bs[x];
        for(int j:fc[x])++bs[j],R_L[i]+=cnt[j];
        ++cnt[x];
    }
    memset0(cnt),memset0(bs);
    sort(q+1,q+m+1);
    q[0].l=1,q[0].r=0;
    for(int i=1;i<=m;++i){
        const node&now=q[i],&pre=q[i-1];
        ans[i]=L_R[now.r]-L_R[pre.r]+R_L[now.l]-R_L[pre.l];
        if(now.r>pre.r)++cnt[pre.l-1];
        if(now.r<pre.r)++cnt[pre.l-1];
        if(now.l<pre.l)++bs[now.r+1];
        if(now.l>pre.l)++bs[now.r+1];
    }
    for(int i=0;i<=n+1;++i)
    L[i].malloc(bs[i]),R[i].malloc(cnt[i]);
    memset0(cnt),memset0(bs);
    for(int i=1;i<=m;++i){
        const node&now=q[i],&pre=q[i-1];
        if(now.r>pre.r)
        R[pre.l-1].push_back((node){pre.r+1,now.r,-i});
        if(now.r<pre.r)
        R[pre.l-1].push_back((node){now.r+1,pre.r,i});
        if(now.l<pre.l)
        L[now.r+1].push_back((node){now.l,pre.l-1,-i});
        if(now.l>pre.l)
        L[now.r+1].push_back((node){pre.l,now.l-1,i});
    }
    for(int i=1;i<=n;++i){
        const int x=a[i];
        for(int j:fc[x])++bs[j];
        if(x>siz)
        for(int j=x;j<=n;j+=x)++bs[j];
        for(node*vc=R[i].begin;vc!=R[i].end;++vc)
        if(vc->id>0){
            LL&t=ans[vc->id];
            for(int k=vc->l;k<=vc->r;++k)t+=bs[a[k]];
        }else{
            LL&t=ans[-vc->id];
            for(int k=vc->l;k<=vc->r;++k)t-=bs[a[k]];
        }
    }
    memset0(cnt),memset0(bs);
    for(int i=n;i;--i){
        const int x=a[i];
        for(int j:fc[x])++bs[j];
        if(x>siz)
        for(int j=x;j<=n;j+=x)++bs[j];
        for(node*vc=L[i].begin;vc!=L[i].end;++vc)
        if(vc->id>0){
            LL&t=ans[vc->id];
            for(int k=vc->l;k<=vc->r;++k)t+=bs[a[k]];
        }else{
            LL&t=ans[-vc->id];
            for(int k=vc->l;k<=vc->r;++k)t-=bs[a[k]];
        }
    }
    for(int w=1;w<=11;++w){
        const int Lt=(w-1)*9+1,Rt=w*9;
        for(int i=1;i<=n;++i)
        for(int j=Lt;j<=Rt;++j)
        pre[i][j-Lt]=pre[i-1][j-Lt]+!(a[i]%j);
        int tot[9];
        memset0(tot);
        for(int i=1;i<=n;++i){
            if(Lt<=a[i]&&a[i]<=Rt)++tot[a[i]-Lt];
            for(node*vc=R[i].begin;vc!=R[i].end;++vc)
            if(vc->id>0)
            for(int j=0;j<9;++j)
            ans[vc->id]+=(LL)tot[j]*(pre[vc->r][j]-pre[vc->l-1][j]);else
            for(int j=0;j<9;++j)
            ans[-vc->id]-=(LL)tot[j]*(pre[vc->r][j]-pre[vc->l-1][j]);
        }
        memset0(tot);
        for(int i=n;i;--i){
            if(Lt<=a[i]&&a[i]<=Rt)++tot[a[i]-Lt];
            for(node*vc=L[i].begin;vc!=L[i].end;++vc)
            if(vc->id>0)
            for(int j=0;j<9;++j)
            ans[vc->id]+=(LL)tot[j]*(pre[vc->r][j]-pre[vc->l-1][j]);else
            for(int j=0;j<9;++j)
            ans[-vc->id]-=(LL)tot[j]*(pre[vc->r][j]-pre[vc->l-1][j]);
        }
    }
    static LL out[N];
    for(int i=1;i<=m;++i)ans[i]+=ans[i-1],out[q[i].id]=ans[i]+q[i].r-q[i].l+1;
    for(int i=1;i<=m;++i)print(out[i]);
    fwrite(bwf,ww-bwf,1,stdout);
    return 0;
}
Powered by Hexo
Theme - Material