莫队略解

command_block · 2018-12-09 11:50:20 · 个人记录

1.经典莫队

抽象概论 : 经典莫队适用于这样的一类问题 :

给出一个长度为 n 的序列，对 q 个区间进行询问。

由区间 [l,r] 的答案可以快速得到区间 [l,r+1],[l,r-1],[l-1,r],[l+1,r] 的答案。即 : 端点可以快速移动。

允许离线。(下面的所有问题若未加说明均默认允许离线)
**题意** : 给出一个序列和参数 $k$ ，每次询问一个区间内有多少子区间的异或和为 $k$。 $n,m,k,a_i\leq 10^5$。

首先做异或前缀和，问题就变为了区间内选取一个对子，使得异或值为 k 的方案数。

对于计算“对子”贡献的问题，可以每次考虑新加入元素和旧元素之间的贡献。

记录 s[x] 为目前区间内 x 的出现次数，假设新加入了 y 那么能够异或出 k 的 x=k\ {\rm xor}\ y ，贡献就是 s[k\ {\rm xor}\ y]。然后将 s[y] 加一。

当去除某个元素时，也容易找到逆操作。

我们发现，通过加删操作，我们可以在 O(1) 内从区间 [l,r] 的答案得到 [l-1,r],[l,r-1],[l+1,r],[l,r+1] ，将区间端点移动。

若上一次回答 [L,R] ，这一次回答 [l,r] ，指针需要跳动 |L-l|+|R-r| 次。

我们发现，询问的回答顺序会影响复杂度。由于允许离线处理，我们可以任意安排回答顺序。

如何找到最优的回答顺序? 这可以转化成平面曼哈顿距离最短哈密顿路问题。

精确地解决这个问题较为困难，不过我们只需要找到一个近似的做法来保证复杂度：分块。

具体方法：对原序列分块，左端点在同一块里的，按照右端点排序，否则按照左端点排序。

例子:

排序前
              +-+
        +-----+
  +-+
    +---------+
+---------+
1 2 3|4 5 6|7 8 9|

排序后(三个数为一块)

  +-+
+---------+
    +---------+
        +-----+
              +-+
1 2 3|4 5 6|7 8 9|

肉眼可见，时间节省了很多。

复杂度证明 ：

设每块的大小为 B ，数列长为 n ，询问个数为 q 。

这些询问最多被分成 n/B 块，块内右端点有序，而且安排在一起回答。

r$ 指针不会往回跳，所以一个块复杂度 $O(n)$ ，共 $O(n/B)$ 块复杂度为 $O(n^2/B)

每个块内的左端点无序，但是它们相距不是很远，最多相距 O(B) 那么远。

所以每个询问左端点最多跳 O(B) 次，总的复杂度 O(qB)

综上所述,总共的复杂度是 O(qB+n^2/B)。

如果令 t=\sqrt{n} ,得复杂度是 O((n+q)\sqrt{n})
如果令 t=\dfrac{n}{\sqrt{q}} ,得复杂度是 O(n\sqrt{q}) (总是优于楼上)

What? 复杂度从 O(n^2) 一下变成了 O(n\sqrt{q}) ?

~~没错,离线真的可以为所欲为。~~

事实上，这个数量级也可以证明为平面最短哈密顿路的下界。

上代码，初学者要注意细节。

指针移动和贡献计算的先后顺序。
先执行将区间变大的操作，避免出现 l>r 的情况(在某些题中可能导致错误)。

#include<algorithm>
#include<cstdio>
#include<cmath>
#define ll long long
#define MaxN 100500
using namespace std;
struct Data{int l,r,p;}b[MaxN];
int BS;
bool cmp(const Data &A,const Data &B)
{return A.l/BS==B.l/BS ? A.r<B.r : A.l<B.l;}
int n,m,k,a[MaxN],s[MaxN<<1];
ll ret[MaxN];
int main()
{
  scanf("%d%d%d",&n,&m,&k);
  BS=n/sqrt(m)+1;
  for (int i=1;i<=n;i++){
    scanf("%d",&a[i]);
    a[i]^=a[i-1];
  }
  for (int i=1;i<=m;i++){
    scanf("%d%d",&b[i].l,&b[i].r);
    b[i].l--;b[i].p=i;
  }sort(b+1,b+m+1,cmp);
  int l=1,r=0;ll ans=0;
  for (int i=1;i<=m;i++){
    while(l>b[i].l){
      l--;
      ans+=s[k^a[l]];
      s[a[l]]++;
    }
    while(r<b[i].r){
      r++;
      ans+=s[k^a[r]];
      s[a[r]]++;
    }
    while(l<b[i].l){
      s[a[l]]--;
      ans-=s[k^a[l]];
      l++;
    }
    while(r>b[i].r){
      s[a[r]]--;
      ans-=s[k^a[r]];
      r--;
    }ret[b[i].p]=ans;
  }
  for (int i=1;i<=m;i++)
   printf("%lld\n",ret[i]);
  return 0;
}

可以看到，莫队算法的实现是较为简洁的。

[P1494 [国家集训队]小Z的袜子](https://www.luogu.org/problemnew/show/P1494) [P2709 小B的询问](https://www.luogu.org/problemnew/show/P2709)
卡常技巧 : 奇偶块排序

仔细分析两个指针的行为。

在某一块内，r 指针从左到右跳一整遍，l 指针每次在块内乱动，而在两个块间切换时，r 指针又要回到序列开头。

我们可以在相邻的两块中，前一块按照 r 升序排序，后一块按照 r 降序排序，这样 r 指针就不需要回头了。

定义一个可重集合的权值为 : 其中的同色无序对数。给出三个数列 $A,B,C$ ，每次给出三个数 $l_A,l_B,l_C$ ，询问 $A[1,l_a],B[1,l_b],C[1,l_C]$ 三个前缀组成的可重集的权值。 $n\leq 5\times 10^4,q\leq 10^5$。

子问题 : 维护一个可重集，支持插入删除，回答权值。

类似 P1494 ，使用桶记下每个元素的出现次数，不难实现 O(1) 插入删除。

如果把询问 (l_A,l_B,l_C) 看做三维空间内的点，我们就可以利用上述算法每次移动一步。

若要寻找耗时最少的回答顺序，则是三维曼哈顿距离最短哈密顿路问题。仍然可以使用分块的方法。

将整个 n\times n\times n 的空间均匀地分成若干立方体区域，假设边长为 B。

那么，在任意两个在同一块内的询问之间跳转的代价，和在相邻两块的询问之间跳转的代价都是 O(B) 的。

如果按照某种顺序依次访问各个块(这显然能够做到)，复杂度则是 O(\text{块数}*B) ，块数是 (n/B)^3。

总复杂度则是 O(n^3/B^2+qB) ，令 n^3/B^2=qB 解得 B=n/\sqrt[3]{q} ，最优复杂度为 O(nq^{2/3}) 。

上面只是复杂度的理论分析，具体实现可见下文“带修莫队”。

扩展 : k 维莫队的复杂度可以做到 O(nq^{1-1/k})。

**题意** : 给出一个序列和若干形为 $(l_1,r_2,l_2,r_2)$ 的询问，求下列式子的值 $$\sum\limits_{x}\text{get}(l_1,r_1,x)\times \text{get}(l_2,r_2,x)$$ $\text{get}(l,r,x)$ 表示计算区间 $[l,r]$ 中，数字 $x$ 出现了多少次。 $n,q\leq 5\times 10^4$。

不难发现，当 l_1,r_1,l_2,r_2 之一移动一步时，只有一个 x 的贡献会改变，且利用桶容易计算。

那么，我们就得到了一个四维莫队的 O(nq^{3/4}) 搞笑做法。但这也告诉我们四个指针同时对付是没有出路的。

考虑差分，设 c(l,x)=\text{get}(1,l,x) ，则有 get(l,r,x)=c(1,r,x)-c(1,l-1,x)。

将欲求式子改写成 :

\sum\limits_{x}\big(c(r_1,x)-c(l_1-1,x)\big)\big(c(r_2,x)-c(l_2-1,x)\big)

拆开括号后可以得到 4 个形如下式的子问题。

\sum\limits_{x}c(r_1,x)c(r_2,x)

这就是二维莫队可以解决的了，复杂度 O(n\sqrt{m})。

#include<algorithm>
#include<cstdio>
#include<cmath>
#define ll long long 
#define MaxN 50500
using namespace std;
int read(){}
struct Data{int p1,p2,p,op;}b[MaxN*4];
int BS;
bool cmp(const Data &A,const Data &B)
{return A.p1/BS==B.p1/BS ? A.p2<B.p2 : A.p1<B.p1;}
int c1[MaxN],c2[MaxN],a[MaxN];
ll ret[MaxN];
int n,m,q;
int main()
{
  n=read();
  for (int i=1;i<=n;i++)a[i]=read();
  q=read();
  for (int i=1,l1,l2,r1,r2;i<=q;i++){
    l1=read();r1=read();
    l2=read();r2=read();
    b[++m]=(Data){l1-1,l2-1,i,1};
    b[++m]=(Data){r1,l2-1,i,-1};
    b[++m]=(Data){l1-1,r2,i,-1};
    b[++m]=(Data){r1,r2,i,1};
  }BS=n/sqrt(m)+1;
  sort(b+1,b+m+1,cmp);
  int p1=0,p2=0;ll ans=0;
  for (int i=1;i<=m;i++){
    while(p1<b[i].p1){
      ans+=c2[a[++p1]];
      c1[a[p1]]++;
    }while(p1>b[i].p1){
      c1[a[p1]]--;
      ans-=c2[a[p1--]];
    }
    while(p2<b[i].p2){
      ans+=c1[a[++p2]];
      c2[a[p2]]++;
    }while(p2>b[i].p2){
        c2[a[p2]]--;
      ans-=c1[a[p2--]];
    }ret[b[i].p]+=b[i].op*ans;
  }
  for (int i=1;i<=q;i++)
    printf("%lld\n",ret[i]);
  return 0;
}

类似的题 : P4689 [Ynoi2016]这是我自己的发明

-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-

2.带修莫队

前面我们用莫队解决的是一类无修改的“静态”问题，现在我们考虑修改。

**题意** : 带修区间数颜色问题。

我们把一个询问描述为 (l,r,t) : 前 t 个修改生效后，区间 [l,r] 的颜色数。

考虑如何在这个三维空间内移动， l,r 的移动方法我们是熟悉的，下面来研究 t。

考虑一个修改对当前状态 [l,r] 的贡献，如果不是改了 [l,r] 内的位置，显然无影响。否则相当于一次加元素和一次删元素。

接下来就是三维莫队的活了。

设块大小为 B。

当左端点不在一个块内，则比较左端点。

反之则考虑右端点，如果不在一个块内，则比较左端点。反之比较时间。

(意思是，若不在一个块内，优先向左走，次之向右走，最末沿时间轴走)

根据前文的分析，令 q=n/q^{1/3} 可得复杂度为 O(nq^{2/3})。

#include<algorithm>
#include<cstdio>
#define MaxN 140000
#define MaxS 1000500
using namespace std;
int read(){
  int X=0;char ch=0;
  while(ch<48||ch>57)ch=getchar();
  while(ch>=48&&ch<=57)X=X*10+(ch^48),ch=getchar();
  return X;
}
bool readc(){
  char ch=getchar();
  while('Q'!=ch&&'R'!=ch)ch=getchar();
  return 'Q'==ch;
}
struct Data
{int l,r,pos,t;}b[MaxN];
struct Ord
{int pos,x,las;}c[MaxN];
int a[MaxN],bl[MaxN],o[MaxS],ret[MaxN];
bool cmp(const Data &A,const Data &B){
  return bl[A.l]==bl[B.l] ? 
        (bl[A.r]==bl[B.r] ? A.t<B.t : A.r<B.r)
         : A.l<B.l;
}
#define add(p) if (++o[a[p]]==1)ans++;
#define del(p) if (--o[a[p]]==0)ans--;
char ord;int n,m,BS=10,tb,tc;
int main()
{
  n=read();m=read();
  while(1ll*BS*BS*BS<=1ll*n*n*5)BS++;
  for (int i=1;i<=n;i++)
    {a[i]=read();bl[i]=i/BS;}
  for (int i=1;i<=m;i++){
    if (readc()){
      b[++tb].l=read();b[tb].r=read();
      b[tb].pos=tb;b[tb].t=tc;
    }else{
      c[++tc].pos=read();
      c[tc].las=a[c[tc].pos];
      a[c[tc].pos]=c[tc].x=read();
    }
  }sort(b+1,b+tb+1,cmp);
  int l=1,r=0,tim=tc,ans=0;
  for (int i=1,p;i<=tb;i++){
    while(tim<b[i].t){
      p=c[++tim].pos;
      if (l<=p&&p<=r)del(p);
      a[p]=c[tim].x;
      if (l<=p&&p<=r)add(p);
    }while(tim>b[i].t){
      p=c[tim].pos;
        if (l<=p&&p<=r)del(p);
      a[p]=c[tim--].las;
      if (l<=p&&p<=r)add(p);
    }while(l<b[i].l)del(l++);
    while(l>b[i].l)add(--l);
    while(r>b[i].r)del(r--);
    while(r<b[i].r)add(++r);
    ret[b[i].pos]=ans;
  }for (int i=1;i<=tb;i++)
     printf("%d\n",ret[i]);
  return 0;
}

由于修改的常数较大，建议把块大小略微调大。加强了数据导致卡常，不开 O_2 过不了。

-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-

3.回滚莫队

回滚莫队用于解决这样一类问题 : 我们向区间中加入元素相对容易，而删除元素较为困难。

也就是说，区间 [l,r] 容易快速转移到 [l-1,r] 或者 [l,r+1]。

比如此题中，我们只需对每个颜色记录最靠左和最靠右的出现，就能资瓷加入元素。

考虑巧妙地加入元素来规避删除操作。

先暴力处理两个端点在同一块内的操作。

观察到同一块中， r 指针只会向右跳(加入)，这部分不会产生删除操作。

问题就在于 l 指针乱跳可能出现删除。

我们可以一开始将 l 置为本块的右边界。

回答每个询问时，先扩展 r 指针，然后再移动 l ，最后将 l 回滚（撤销操作），回到右边界。

(本质上是需要资瓷快速计算临时的 O(\sqrt{n}) 个加入的贡献)

这样子 l 指针常数大了 2~3 倍，而且跨越两块需要清空，写得不好复杂度可能变为O((n+m)\sqrt{n})。

#include<algorithm>
#include<cstdio>
#define MaxN 200500
using namespace std;
inline int read(){}
int BS,bl[MaxN];
struct Data
{int l,r,tl,tr,p;}b[MaxN];
bool cmp(const Data &A,const Data &B)
{return bl[A.l]==bl[B.l] ? A.r<B.r : A.l<B.l;}
int x[MaxN],a[MaxN],l[MaxN],r[MaxN],sl[MaxN],sr[MaxN]
   ,n,q,m,ret[MaxN],ans;
void clr()
{for (int k=1;k<=m;k++){l[k]=n;r[k]=0;}}
int main()
{
  n=read();for (BS=5;BS*BS*2<=n;BS++);
  for (int i=1;i<=n;i++){
    x[i]=a[i]=read();
    bl[i]=i/BS;
  }sort(x+1,x+n+1);
  m=unique(x+1,x+n+1)-x-1;
  for (int i=1;i<=n;i++)
    a[i]=lower_bound(x+1,x+m+1,a[i])-x;
  q=read();
  for (int i=1;i<=q;i++){
    b[i].l=read();b[i].r=read();
    b[i].p=i;
  }sort(b+1,b+q+1,cmp);clr();
  int tl=(1/BS+1)*BS-1,tr=tl-1;
  for (int i=1;i<=q;i++){
    if (bl[b[i].l]!=bl[b[i-1].l]){
      tl=(b[i].l/BS+1)*BS-1;
      tr=tl-1;ans=0;clr();
    }
    if (bl[b[i].l]==bl[b[i].r]){
      int sav=0;
      for (int k=b[i].l;k<=b[i].r;k++)sr[a[k]]=k;
      for (int k=b[i].l;k<=b[i].r;k++)
        sav=max(sav,sr[a[k]]-k);
      ret[b[i].p]=sav;
      continue;
    }
    while(tr<b[i].r){
      ++tr;
      if (l[a[tr]]==n)l[a[tr]]=tr;
      r[a[tr]]=tr;
      ans=max(ans,tr-l[a[tr]]);
    }
    int sav=ans,sp=tl;
    for (int k=b[i].l;k<=tl;k++)
      {sl[a[k]]=l[a[k]];sr[a[k]]=r[a[k]];}
    while(b[i].l<tl){
      --tl;
      l[a[tl]]=tl;
      if (r[a[tl]]==0)r[a[tl]]=tl;
      ans=max(ans,r[a[tl]]-tl);
    }ret[b[i].p]=ans;
    tl=sp;ans=sav;
    for (int k=b[i].l;k<=tl;k++)
      {l[a[k]]=sl[a[k]];r[a[k]]=sr[a[k]];}
  }
  for (int i=1;i<=q;i++)
    printf("%d\n",ret[i]);
  return 0;
}

4.树上莫队

特指解决树链问题的莫队。显然，这个问题强于区间莫队。

**题意** : 树链数颜色。

全dfs序 : 在 dfs 每次到达某个点时(入栈)，在序列中加入该点。离开某个点时，再加入一次，得到长为 2n 的序列。

称 dfn[u] 为 u 第一次被写入欧拉序的位置，out[u] 为第二次被写入欧拉序的位置。
结论 : 路径 u,v （不妨设 dfn[u]\leq dfn[v] ）中的点集为 :

全 dfs 序列 \big[out[u],dfn[v]\big] 之间恰好出现一次的点并上 {\rm lca}(u,v) 。
证明 : S=\big[out[u],dfn[v]\big] 中的全 dfs 序其实就是记录了在 u,v 之间 dfs 的过程。

令 t={\rm lca}(u,v) ，那么 u,t 之间的点的入栈操作均不在 S 中，而出栈操作都在 S 中。

对于 v,t 之间的点，出栈操作均不在 S 中，而入栈操作都在 S 中。

而之间经过的分支不属于 u,v 之间的路径，相应的，其出入栈都在 S 内。

对于 t ，则需要特殊考虑。

由于这个性质，我们只需将第二次加入看做删除即可。这样，我们就将树链问题转化为了区间问题。

例题② : 在路上了

我们可以用树上回滚莫队将其解决。但回滚莫队不支持删除，无法使用前文中全 dfs 序的做法。

于是，我们考虑另一种转化。

综合练习 ： P4074 [WC2013] 糖果公园

5.二次离线莫队

例题： P5047 [Ynoi2019 模拟赛] Yuno loves sqrt technology II、

题意：查询区间逆序对数。允许离线。

首先有一个显然的 O(n\sqrt{n}\log n) 暴力。

用莫队维护，并用树状数组维护当前集合。每次加入时计算当前数与其他数之间的贡献，只需查询当前集合中 \leq x 的数的个数。

该算法需要 O(n\sqrt{n}) 次单点修改与前缀查询，故无法使用简单的根号平衡来解决。

注意到，我们在指针移动时所需的询问形如： “在区间 [l,r] 内有多少个数 \leq k”，而这可以差分为“在前缀 [1,r] 内有多少个数 \leq k”。

我们可以用可持久化分块（O(\sqrt{n}) 修改 O(1) 查询）来支持 O(1) 求解这个问题。这样可以做到时空复杂度 O(n\sqrt{n}) ，但常数较大。

由于莫队处理的是离线问题，我们可以将需要的 O(n\sqrt{n}) 个询问再次离线下来，然后按照 r 排序，用普通值域分块一次处理（扫描线），这样就不需要可持久化了。很遗憾，记下询问所需的空间复杂度仍然是 O(n\sqrt{n}) 的。

为了优化空间复杂度（以及减小时间常数），我们进一步观察莫队产生的询问的特殊性质：

记 (r,x)= “在前缀 [1,r] 内 \leq x 的数的数目”

l$ 指针向左移动到 $l'$ ： $[l,r],[l-1,r],[l-2,r]\dots[l',r]
其中的贡献为 \sum\limits_{i\in[l',l)}(r,A_i-1)-(i,A_i-1) 。

向右移动只需将贡献取反。
r$ 指针向右移动到 $r'$ ： $[l,r],[l,r+1],[l,r+2]\dots[l,r']
其中的贡献为 \sum\limits_{i\in(r,r']}(i-l+1)-\big((i-1,A_i)-(l-1,A_i)\big) 。

对于 \sum\limits_{i\in(r,r']}(i-l+1) 部分容易计算。

向左移动只需将贡献取反。

观察其中的贡献形式，可以分类为以下四种：

A:(l,r,x)$ 表示 $\sum\limits_{i\in[l,r]}(x,A_i-1)
D:(l,r,x)$ 表示 $\sum\limits_{i\in[l,r]}(x,A_i)

只有 A,D 两种需要扫描线。

我们只需要将 (l,r,x) 的三元组挂在 l 处，动态维护目前生效的三元组，若期限 r 到了，就将其删除。

于是空间复杂度优化到了 O(n) ，时间常数也大为减小。

#include<algorithm>
#include<cstdio>
#include<vector>
#include<cmath>
#define ll long long
#define pb push_back
#define MaxN 100500
using namespace std;
int n;
struct SumDS
{
  #define lbit(x) (x&-x)
  int t[MaxN];
  void add(int p)
  {while(p<=n){t[p]++;p+=lbit(p);}}
  int qry(int p){
    int ret=0;
    while(p){ret+=t[p];p^=lbit(p);}
    return ret;
  }
}T;
struct SumDS2
{
  int ob[505],o[MaxN],BS;
  void add(int p){
    int bp=p/BS+1,br=bp*BS;
    for (int i=bp;i<=BS;i++)ob[i]++;
    for (int i=p;i<br;i++)o[i]++;
  }
  inline int qry(int p)
  {return ob[p/BS]+o[p];}
}T2;
ll ans[MaxN];
struct Data2{int l,r,c,d,p;};
vector<Data2> s[MaxN];
struct Data{int l,r,p;}b[MaxN];
ll s1[MaxN],s2[MaxN];
int q,BS;
bool cmp(const Data &A,const Data &B)
{return A.l/BS==B.l/BS ? ((A.l/BS)&1)^(A.r<B.r) : A.l<B.l;}
void MoQ()
{
  sort(b+1,b+q+1,cmp);
  int l=1,r=0;
  for (int i=1;i<=q;i++){
    if (b[i].l<l){
      ans[b[i].p]-=s1[l-1]-s1[b[i].l-1];
      s[r].pb((Data2){b[i].l,l-1,1,-1,b[i].p});
      l=b[i].l;
    }
    if(r<b[i].r){
      ans[b[i].p]+=1ll*(r+b[i].r+1-2*l)*(b[i].r-r)/2-(s2[b[i].r]-s2[r]);
      s[l-1].pb((Data2){r+1,b[i].r,1,0,b[i].p});
      r=b[i].r;
    }
    if(l<b[i].l){
      ans[b[i].p]+=s1[b[i].l-1]-s1[l-1];
      s[r].pb((Data2){l,b[i].l-1,-1,-1,b[i].p});
      l=b[i].l;
    }
    if(b[i].r<r){
      ans[b[i].p]-=1ll*(r+b[i].r+1-2*l)*(r-b[i].r)/2-(s2[r]-s2[b[i].r]);
      s[l-1].pb((Data2){b[i].r+1,r,-1,0,b[i].p});
      r=b[i].r;
    }
  }
}
int a[MaxN],x[MaxN];
int main()
{
  scanf("%d%d",&n,&q);
  for (int i=1;i<=n;i++)
    {scanf("%d",&a[i]);x[i]=a[i];}
  sort(x+1,x+n+1);
  for (int i=1;i<=n;i++)
    a[i]=lower_bound(x+1,x+n+1,a[i])-x;
  for (int i=1;i<=n;i++){
    s2[i]=s2[i-1]+T.qry(a[i]);
    T.add(a[i]);
    s1[i]=s1[i-1]+T.qry(a[i]-1);
  }
  for (int i=1;i<=q;i++){
    scanf("%d%d",&b[i].l,&b[i].r);
    b[i].p=i;
  }
  BS=1.0*n/sqrt(q)+5;
  MoQ();
  T2.BS=sqrt(n)+1;
  for (int i=1;i<=n;i++){
    T2.add(a[i]);
    for (int k=0;k<s[i].size();k++){
      Data2 u=s[i][k];
      for (int j=u.l;j<=u.r;j++)
        ans[u.p]+=u.c*T2.qry(a[j]+u.d);
    }
  }
  for (int i=1;i<=q;i++)ans[b[i].p]+=ans[b[i-1].p];
  for (int i=1;i<=q;i++)printf("%lld\n",ans[i]);
  return 0;
}