后缀自动机(SAM)+广义后缀自动机(GSA)

 

经过一顿操作之后竟然疑似没退役0 0

你是XCPC选手吗?我觉得我是!

稍微补一点之前丢给队友的知识吧,除了数论以外都可以看看,为Dhaka和新队伍做点准备...

不错的零基础教程见 IO WIKI - 后缀自动机,这篇就从自己的角度总结一下吧,感觉思路总是和别人不一样...尽量用我比较能接受的语言和逻辑重新组织一遍。

注意:在本文中,字符串的下标从1开始。

目前的 SAM 模板:

const int N=100005;
const int C=26; //注意检查字符集大小!

//在结构题外开任何与SAM状态相关的数组,都需要N<<1
struct SuffixAutomaton
{
    int sz,lst; //状态数上限=2|S|-1 
    
    int len[N<<1],link[N<<1];
    int nxt[N<<1][C];
    //map<char,int> nxt[N<<1];
    //extend(char),并使用nxt[clone]=nxt[q]替换memcpy
    
    SuffixAutomaton()
    {
        len[0]=0,link[0]=-1;
        lst=sz=0;
    }
    
    void extend(int c)
    {
        int cur=++sz;
        len[cur]=len[lst]+1;
        
        int p=lst;
        while(p!=-1 && !nxt[p][c])
            nxt[p][c]=cur,p=link[p];
        
        if(p==-1)
            link[cur]=0;
        else
        {
            int q=nxt[p][c];
            if(len[p]+1==len[q])
                link[cur]=q;
            else
            {
                int clone=++sz;
                len[clone]=len[p]+1;
                link[clone]=link[q];
                memcpy(nxt[clone],nxt[q],C*4);
                
                while(p!=-1 && nxt[p][c]==q)
                    nxt[p][c]=clone,p=link[p];
                link[q]=link[cur]=clone;
            }
        }
        lst=cur;
    }
    
    int id[N<<1];
    
    //构建完成后,id顺序为len递增(逆拓扑序)【仅可排一次】 
    void sort()
    {
        static int bucket[N<<1];
        
        memset(bucket,0,sizeof(bucket));
        for(int i=1;i<=sz;i++)
            bucket[len[i]]++;
        for(int i=1;i<=sz;i++)
            bucket[i]+=bucket[i-1];
        for(int i=1;i<=sz;i++)
            id[bucket[len[i]]--]=i;
    }
}sam;
View Code

目前的 GSA 模板:

typedef pair<int,int> pii;

const int N=1000005;
const int C=26; //注意检查字符集大小!

//在结构题外开任何与SAM状态相关的数组,都需要N<<1
struct GeneralSuffixAutomaton
{
    int sz; //状态数上限=2|S|-1 
    
    int len[N<<1],link[N<<1];
    int nxt[N<<1][C];
    //map<char,int> nxt[N<<1];
    //extend(int,char)
    
    GeneralSuffixAutomaton()
    {
        len[0]=0,link[0]=-1;
        sz=0;
    }
    
    //向trie树上插入一个字符串 
    void insert(char *s,int n)
    {
        for(int i=1,cur=0;i<=n;i++)
        {
            int c=s[i]-'a';
            if(!nxt[cur][c])
                nxt[cur][c]=++sz;
            cur=nxt[cur][c];
        }
    }
    
    //扩展GSA,对于lst的c转移进行扩展 
    int extend(int lst,int c)
    {
        int cur=nxt[lst][c];
        len[cur]=len[lst]+1;
        
        int p=link[lst];
        while(p!=-1 && !nxt[p][c])
            nxt[p][c]=cur,p=link[p];
        
        if(p==-1)
            link[cur]=0;
        else
        {
            int q=nxt[p][c];
            if(len[p]+1==len[q])
                link[cur]=q;
            else
            {
                int clone=++sz;
                len[clone]=len[p]+1;
                link[clone]=link[q];
                for(int i=0;i<C;i++)
                    nxt[clone][i]=len[nxt[q][i]]?nxt[q][i]:0;
                
                while(p!=-1 && nxt[p][c]==q)
                    nxt[p][c]=clone,p=link[p];
                link[q]=link[cur]=clone;
            }
        }
        return cur;
    }
    
    //基于trie树建立GSA
    void build()
    {
        queue<pii> Q;
        for(int i=0;i<C;i++)
            if(nxt[0][i])
                Q.push(pii(0,i));
        
        while(!Q.empty())
        {
            pii tmp=Q.front();
            Q.pop();
            
            int cur=extend(tmp.first,tmp.second);
            for(int i=0;i<C;i++)
                if(nxt[cur][i])
                    Q.push(pii(cur,i));
        }
    } 
    
    int id[N<<1];
    
    //构建完成后,id顺序为len递增(逆拓扑序)【仅可排一次】 
    void sort()
    {
        static int bucket[N<<1];
        
        memset(bucket,0,sizeof(bucket));
        for(int i=1;i<=sz;i++)
            bucket[len[i]]++;
        for(int i=1;i<=sz;i++)
            bucket[i]+=bucket[i-1];
        for(int i=1;i<=sz;i++)
            id[bucket[len[i]]--]=i;
    }
}gsa;
View Code

 

简单的目录:

背景

概念:endpos 等价类

概念:SAM 的状态

概念:后缀链接 link

SAM 的建立

正确性证明(其实没写啥...)

SAM 的性质

经典应用(SAM)

从 SAM 到广义 SAM

广义 SAM(GSA)的建立

经典应用(GSA) 

一些练习题

 

 

~ 背景 ~ (SAM 要解决什么问题,为什么“暴力”方法不行)


很久以前就听说过 SAM,一直有一个不明觉厉的印象,导致一直没敢看。但其实 SAM 本质上来说是一个偏板子的知识点(甚至板子也不算长),其中比较困难的地方在于对于 SAM 结构的理解,以及再次基础上对于其性质的利用。网络上的各种材料,感觉大多从一个已经有一定理解的角度来进行讲述,一开始看下来感觉有点概念不清晰。

简单来说,SAM 想要解决的是一个这样的问题:首先对于字符串 $s$ 建立一个自动机,再将另外的字符串$t$扔到自动机里面逐字符进行匹配。假设目前已经匹配到了 $t[i]$,那么这个匹配到的状态表示的意思是“$t[1:i]$ 的 $i$ 个后缀中,最长的出现在 $s$ 中的连续子字符串”。

举个例子,$s=abb$、$t=abbcab$,那么匹配到 $t[6]=b$ 时,匹配到的状态应该是 $ab$。虽然 $s$ 与 $t$ 的最长公共连续子字符串为 $abb$,但是因为 $abb$ 不是 $t[1:6]$ 的后缀,所以不是需要求的结果。

对于这个问题有一种非常暴力的做法,就是用 $s$ 的 $\dfrac{|s|\cdot (|s|-1)}{2}$ 个子串来建立 AC 自动机。在这个自动机上我们能够匹配到 $s$ 的任意连续子字符串,于是能够满足我们的需求,但是状态数直接爆炸了。我们需要一个类似这个 AC 自动机的结构(包含状态、转移边和 fail 边的自动机,通过转移边进行一次匹配,失配的时候通过 fail 边往回跳),但是状态数是 $O(|s|)$ 的。

 

 

~ 概念:endpos 等价类 ~ (规定一种将多个子串归入到同一个状态的规则)

 

为了将状态由 $O(|s|^2)$ 压缩到 $O(|s|)$,我们需要一个状态能够表示多个 $s$ 的子串。我们规定,若两个 $s$ 的子串 endpos 等价,则其属于同一个 SAM 状态

那么我们首先需要了解 endpos 的概念。对于子串 $s'\in s$,其 endpos 等于其出现在 $s$ 中的右端位置的集合。举个例子来说,$s=abcabbacab$、$s'=ab$,我们可以通过肉眼发现 $s'$ 在 $s$ 中出现的位置有$[1:2]$、$[4:5]$、$[9:10]$,那么有 $\mathrm{endpos}(s')=\{2,5,10\}$。需要注意的是,$s'$ 在 $s$ 中出现的位置可能重叠,比如$s=aaaa$、$s'=aa$,那么 $s'$ 的出现位置为$[1:2]$、$[2:3]$、$[3:4]$,此时 $\mathrm{endpos}(s')=\{2,3,4\}$。

对于两个 $s$ 的子串$s_1'$、$s_2'$,如果有 $\mathrm{endpos}(s_1')=\mathrm{endpos}(s_2')$,那么我们称 $s_1'$ 与 $s_2'$ endpos 等价,其属于同一个 endpos 等价类

下面以字符串 $s=abab$ 为例,对于其所有子串进行等价类划分:

子串 endpos 子串 endpos
abab $\{4\}$ ba $\{3\}$
aba $\{3\}$ a $\{1,3\}$
bab $\{4\}$ b $\{2,4\}$
ab $\{2,4\}$ $\empty$ /

根据上表,一共可以划分出 $5$ 个等价类:$\mathrm{endpos}(ba,aba)=\{3\}$,$\mathrm{endpos}(bab,abab)=\{4\}$,$\mathrm{endpos}(a)=\{1,3\}$,$\mathrm{endpos}(b,ab)=\{2,4\}$,以及需要单列出来的空串 $\empty$(我们可以认为 $\mathrm{endpos}(\empty)=\{0,1,2,\cdots,|s|\}$,包含 $0$ 是为了区分形如 $s=aa\cdots a$ 时的 $\mathrm{endpos}(a)$)。

 

 

~ 概念:SAM 的状态 ~ (同一 / 不同状态中的子串具有什么性质与联系,这对于理解建立 SAM 的过程十分关键)

 

有了 endpos 等价类的定义,我们考虑属于同一等价类的子串具有什么样的性质。

首先,假如两个子串 $s_1'$ 与 $s_2'$ 具有相同的 $\mathrm{endpos}$,考虑一个具体的出现位置 $r\in \mathrm{endpos}$ 则有 $s_1'=s[l_1:r]$、$s_2'=s[l_2:r]$,那么显然两者必须存在后缀关系

接着,我们考虑字符串 $s$ 的两个子串 $s_{long}'$ 与 $s_{short}'$,其中 $s_{short}'$ 是 $s_{long}'$ 的一个后缀。对于这两个子串的 endpos,显然每个 $s_{long}'$ 出现的位置上均出现了$s_{short}'$,这也就是说有:

\[\mathrm{endpos}(s_{long}')\subseteq \mathrm{endpos}(s_{short}')\]

这个式子说明,对于某个子串 $s'$ 的所有后缀,其 $\mathrm{endpos}$ 集合的大小随着后缀长度的增加而单调不增,那么属于同一状态(即 endpos 相同)的所有子串应该是一些长度连续的、互为后缀的字符串。【这里说实话有点绕,需要努力理解一下】

引用 OI WIKI 的表述,也许能更好的帮助理解 “长度连续、互为后缀”:

对于每一个状态 $v$,一个或多个子串与之匹配。我们记 $longest(v)$ 为其中最长的一个字符串,记 $len(v)$ 为它的长度。类似地,记 $shortest(v)$ 为最短的子串,它的长度为 $minlen(v)$。那么对应这个状态的所有字符串都是字符串 $longest(v)$ 的不同的后缀,且所有字符串的长度恰好覆盖区间 $[minlen(v), len(v)]$ 中的每一个整数。

举一个例子,假如字符串为 ababcabcd,那么集合划分如下,可以验证一下上述性质:

字符串 endpos 字符串 endpos 字符串 endpos
a $\{1,4\}$ c $\{3,6\}$ ca $\{4\}$
b $\{2,5\}$ bc bca
ab abc abca
cab $\{5\}$ cabc $\{6\}$  
bcab bcabc
abcab abcabc

虽然在实际中我们并不会存下一个 SAM 状态中的全部字符串,但暂时可以将每个状态形象地理解成 endpos 等价的的字符串集合,从而帮助理解建立 SAM 的过程。

从自动机匹配的角度来看的话,如果我们在 SAM 上匹配到了某一状态,这就表示自动机接受了这个字符串集合中的其中一个字符串(其实我们有办法确定具体是哪一个,但要等到学完怎么建 SAM 后才能说清楚)。

 

 

~ 概念:后缀链接 link ~ (SAM 中的 “fail” 转移边)

 

除了状态以外,自动机的另一个重要的组成部分就是转移。与 AC 自动机十分相似的是,后缀自动机也有 $\mathrm{nxt}$(匹配成功时)与 $\mathrm{link}$(匹配失败时,对应着 AC 自动机中的 $\mathrm{fail}$)两类转移。其中 $\mathrm{nxt}$ 的意义与 AC 自动机中的 $\mathrm{nxt}$ 完全一致,就不做展开了;需要稍微解释一下的是后缀链接 $\mathrm{link}$。

后缀链接的定义为:对于 SAM 中的一个状态 $state$,其对应 endpos 等价类中的字符串均为 $\mathrm{longest}(state)$ 的后缀(这是之前的结论);我们取一个最长的 $\mathrm{longest}(state)$ 的后缀,使得其与 $\mathrm{longest}(state)$ 不 endpos 等价(其实这个后缀就是 $\mathrm{shortest}(state)$ 再去掉最左边的一个字符,长度为 $\mathrm{minlen}(state)-1$),我们就将这个最长后缀所在的状态赋给 $\mathrm{link}(state)$。于是我们得到一个挺常用的性质:

\[\mathrm{minlen}(state)=\mathrm{len}\left(\mathrm{link}(state)\right)+1\]

根据上述定义,我们可以对上面的例子得到后缀链接(红色箭头):

我们尝试理解一下为什么失配时要跳后缀链接。对于当前已匹配到的状态 $state$,有$\mathrm{endpos}(state)=\{p_1, p_2,\cdots, p_n\}$,那么在状态 $state$ 尝试匹配字符 $c$ 失配,就表示 $s[p_1+1], s[p_2+1], \cdots, s[p_n+1]$ 均不等于 $c$(否则就存在 $\mathrm{nxt}(c)$ 的转移边),那么我们就需要舍弃一部分已匹配串的前缀(即取已匹配串的一个后缀)后再次尝试匹配。当我们通过后缀链接跳到 $\mathrm{link}(state)$ 时,这个状态的 endpos 集合更大了一些,我们说不定能找到某个新增的 $p_i$ 使得 $s[p_i+1]=c$,从而完成匹配。从这个角度来看,转移到 $\mathrm{link}(state)$ 时我们只需要舍弃最小的已匹配长度,就会带来 endpos 集合的变大,而 endpos 集合的变大则有可能新增一些 $\mathrm{nxt}$ 转移使得待匹配的字符 $c$ 被接受,这个逻辑和 AC 自动机的 $\mathrm{fail}$ 边是完全一致的。

 

 

~ SAM 的建立 ~

 

SAM 的建立是一个在线的过程。假如我们需要对于字符串 $s$ 构建 SAM,我们依次往 SAM 中添加每一个字符 $s[i]$ 来将 SAM 进行一次扩展,接下来细说一次扩展中需要做的事情。

首先,在未添加任何字符时,有一个初始状态对应着空集 $\empty$ 的等价类,我们将其记作 $state_0$。

在进行扩展的过程中,我们维护了一个变量:$last$。其表示在扩展字符 $s[i]$ 之前,字符串 $s[1:i-1]$ 在部分建成的 SAM 中所属的状态。初始时 $last=state_0$。

接下来,我们尝试在当前串 $s[1:i-1]$ 的后面再插入一个字符 $s[i]=c$。我们首先创建一个新的节点 $cur$,并赋值 $\mathrm{len}(cur)=\mathrm{len}(last)+1$,此时 $last$ 的转移边及 $\mathrm{link}(cur)$ 仍未确定,需要进行如下的讨论:

1. 我们从 $last$ 开始,沿着后缀链接 $\mathrm{link}$ 一直走,边走边将途径节点(包含 $last$)的 $\mathrm{nxt}(c)$ 均赋值为 $cur$,直到某个节点 $p$ 存在 $\mathrm{nxt}(c)$ 的转移边才停下。假如这样的节点 $p$ 不存在,那么我们会一直走到 $state_0$,此时我们直接令 $\mathrm{link}(cur)=state_0$ 即可完成操作。

2. 假如我们在一个节点 $p$ 停下,那么我们记 $p$ 通过字符 $c$ 转移到的节点为 $q$,即有 $\mathrm{nxt}_p(c)=q$。现在考虑 $\mathrm{len}(p)$ 与 $\mathrm{len}(q)$ 之间的关系。假如 $\mathrm{len}(p)+1=\mathrm{len}(q)$ 成立,那么我们仅需令 $\mathrm{link}(cur)=q$ 即可完成操作。

3. 假如有 $\mathrm{len}(p)+1\neq \mathrm{len}(q)$,那么需要通过依次通过四步来扩展 SAM:

  (1) 新建一个节点 $clone$,设定 $\mathrm{len}(clone)=\mathrm{len}(p)+1$,并复制 $q$ 的所有 $\mathrm{nxt}$ 转移。例如,假如有 $\mathrm{nxt}_q(c)=tmp$,就同样令 $\mathrm{nxt}_{clone}(c)=tmp$;

  (2) 同样是对于 $clone$,复制 $q$ 的后缀链接(即令 $\mathrm{link}(clone)=\mathrm{link}(q)$),同时抹去 $\mathrm{link}(q)$

  (3) 从 $p$ 开始(包含 $p$)沿着后缀链接 $\mathrm{link}$ 一直走到 $state_0$。假设在这个过程中经过某个节点 $x$,且存在 $\mathrm{nxt}_x(c)=q$ 的转移,就将该转移重定向至 $clone$(即令 $\mathrm{nxt}_x(c)=clone$);

  (4) 将 $q$ 和 $cur$ 的后缀链接指向 $clone$,即令 $\mathrm{link}(q)=\mathrm{link}(cur)=clone$。
做完以上四步以后,最麻烦的 $\mathrm{len}(p)+1\neq \mathrm{len}(q)$ 的 case 也已经解决了。

(上图中实际上节点 $p$ 也是需要重定向的,但为了表示的方便,没有画出重定向 $\mathrm{nxt}_p(c)$ 的箭头)

最后的最后,我们令 $last=cur$,为扩展下一个字符做准备。

虽然说了那么多,其实代码实现相当简单,与上面的讨论一一对应,甚至没有什么注释的必要:

int sz,lst; //状态数上限=2|S|-1 
int len[N<<1],link[N<<1];
int nxt[N<<1][C];
//初始情况下,sz=lst=0,len[0]=0,link[0]=-1
//令link[0]=-1是为了标记后缀链接路径的终点 

void extend(int c)
{
    int cur=++sz;
    len[cur]=len[lst]+1;
    
    int p=lst;
    while(p!=-1 && !nxt[p][c])
        nxt[p][c]=cur,p=link[p];
    
    if(p==-1)
        //case #1
        link[cur]=0;
    else
    {
        int q=nxt[p][c];
        if(len[p]+1==len[q])
            //case #2
            link[cur]=q;
        else
        {
            //case #3
            int clone=++sz;
            len[clone]=len[p]+1;
            link[clone]=link[q];
            memcpy(nxt[clone],nxt[q],C*4);
            
            while(p!=-1 && nxt[p][c]==q)
                nxt[p][c]=clone,p=link[p];
            link[q]=link[cur]=clone;
        }
    }
    lst=cur;
}

 

 

~ 正确性证明 ~

 

细节有点多(懒),具体可以参考 OI WIKI - SAM - 正确性证明 以及 OI WIKI - SAM - 对操作数为线性的证明不过对于抄板子不是那么关键。

为了便于感性的理解,我们不妨举一个扩展 SAM 的过程中最复杂的 Case #3 的例子。

字符串 $s=abb$,当仅扩展了 $ab$ 与扩展完了 $abb$ 时,SAM 的结构分别如下:

可以看出,做一次 Case #3 的扩展相当于对于一个 endpos 等价类(如图中的 $\{b,ab\}$)做拆分。更深入的分析就不展开了。

 

 

~ SAM 的性质 ~

 

对于一个长度为 $n$ 的字符串,SAM 的状态数是 $2n$ 级别的,转移数是 $3n$ 级别的。这个性质没什么特殊的应用,记得数组别开小就行了。

比较关键的性质在于 SAM 的后缀链接 $\mathrm{link}$ 上。

 

1. SAM 的后缀链接构成一个树形结构

我们不妨回顾一下后缀链接的定义:对于 SAM 中的一个状态 $state$,我们取一个最长的 $\mathrm{longest}(state)$ 的后缀,使得其与 $\mathrm{longest}(state)$ 不 endpos 等价,$\mathrm{link}(state)$ 就指向该后缀所属于的状态。

这也就是说,对于任意状态 $state$ 的后缀链接所指向的状态 $\mathrm{link}(state)$,$state$ 所包含的子串长度均长于 $\mathrm{link}(state)$ 中的子串长度(即 $\mathrm{minlen}(state)>\mathrm{len}(\mathrm{link}(state))$)。利用这个性质,我们发现后缀链接不会构成环,因为沿着后缀链接走时,所经过状态中的子串长度严格单调减

同时,由于所有后缀链接的路径一定止于 $state_0$,所以所有后缀链接一定构成树形结构。有说法称其为 parent 树。

 

2. 在 parent 树上,祖先对应的字符串是子孙的后缀

上文已经回答过了,“沿着后缀链接走时,所经过状态中的子串长度严格单调减”。

稍微推广一下,我们可以求出多个字符串的最长公共后缀:对于字符串 $s_1,\cdots,s_n$,我们找到其在 SAM 上的位置 $x_1,\cdots,x_n$,那么最长公共后缀就是 $\mathrm{LCA}(x_1,\cdots,x_n)$。

 

3. 若某个状态存在 $\mathrm{nxt}(c)$ 转移,那么其在 parent 树上的祖先均存在 $\mathrm{nxt}(c)$ 转移

子孙存在 $\mathrm{nxt}(c)$ 转移,说明子孙出现的某个位置 $r_i\in \mathrm{endpos}(\text{子孙})$ 后面有一个字符 $c$。因为祖先是子孙的后缀,所以子孙出现的所有地方祖先均出现,那么显然有 $r_i\in \mathrm{endpos}(\text{祖先})$,则祖先也存在 $\mathrm{next}(c)$ 转移。

 

4. 每个状态的 $\mathrm{endpos}$ 集合为其在 parent 树中子树 $\mathrm{endpos}$ 集合的并集

我们在构建 SAM 的过程中,将 $last$ 所经过的状态提取出来,则这些状态分别包含子串 $s[1:i]$,$1\leq i\leq n$。那么我们为 $s[1:i]$ 所在的状态打上 $i$ 的标签,表示 $i$ 为该状态的一个出现位置。然后我们考虑该状态在 parent 树上的祖先们,由于祖先是子孙的后缀,那么 $i$ 也是祖先们的一个出现位置,那么我们可以将 $i$ 的标签打在 $s[1:i]$ 所在状态一直到 $state_0$ 的路径上。对于所有 $s[1:i]$ 均这样打标签,最终每个状态的 $\mathrm{endpos}$ 就是被打上的标签集合。

不过这个结论暂时还不够严谨,我们需要证明 $s[1:i]$ 所在状态的子孙均不出现在位置 $i$。这也并不困难,因为一个子串出现在位置 $i$ 当且仅当其为 $s[1:i]$ 的一个后缀,而 $s[1:i]$ 的后缀仅会存在于 $s[1:i]$ 所在状态到 $state_0$ 的路径上,一定不在 $s[1:i]$ 所在状态的子树中。

 

淦,虽然感觉这里有很多性质可以挖,但写的时候只能列出比较显然的几条...另外的一些 DP 相关的性质需要在实际应用中体现吧。

 

补充5. 在 SAM 上,对于每个节点暴力跳后缀链接到 $state_0$,总的时间复杂度是 $\mathrm{O}(n^{1.5})$

虽然在少量博客上看到这个性质,但是并没有找到证明,不是十分确定。不过有题目用到了这个性质(如 2019 ICPC 徐州的 L 题),所以大概是对的。从直观上确实比较难卡出很坏的情况。

 

 

~ 经典应用(SAM) ~

 

1. 所有本质不同的子串及其出现次数(SAM 性质 / parent 树上 DP)

由于 SAM 中所有状态所表示的字符串都是互不相同的,所以本质不同的子串数量就是$\sum_{state}\mathrm{len}(state)-\mathrm{minlen}(state)+1$,即把每个状态中包含的字符串数量全部加起来。

例题:LOJ2033 「SDOI2016」生成魔咒

#include <map>
#include <cstdio>
#include <vector>
#include <cstring>
#include <algorithm>
using namespace std;

typedef long long ll;

const int N=100005;
const int C=26; //注意检查字符集大小!

ll ans=0;

struct SuffixAutomaton
{
    int sz,lst; //状态数上限=2|S|-1 
    
    int len[N<<1],link[N<<1];
    //int nxt[N<<1][C];
    map<int,int> nxt[N<<1];
    //extend(char),并使用nxt[clone]=nxt[q]替换memcpy
    
    SuffixAutomaton()
    {
        len[0]=0,link[0]=-1;
        lst=sz=0;
    }
    
    void extend(int c)
    {
        int cur=++sz;
        len[cur]=len[lst]+1;
        
        int p=lst;
        while(p!=-1 && !nxt[p][c])
            nxt[p][c]=cur,p=link[p];
        
        if(p==-1)
            link[cur]=0;
        else
        {
            int q=nxt[p][c];
            if(len[p]+1==len[q])
                link[cur]=q;
            else
            {
                int clone=++sz;
                len[clone]=len[p]+1;
                link[clone]=link[q];
                //memcpy(nxt[clone],nxt[q],C*4);
                nxt[clone]=nxt[q];
                
                while(p!=-1 && nxt[p][c]==q)
                    nxt[p][c]=clone,p=link[p];
                link[q]=link[cur]=clone;
            }
        }
        lst=cur;
        ans+=len[cur]-len[link[cur]];
    }
}sam;

int n;
int a[N];

int main()
{
    scanf("%d",&n);
    for(int i=1;i<=n;i++)
        scanf("%d",&a[i]);
    
    for(int i=1;i<=n;i++)
        sam.extend(a[i]),printf("%lld\n",ans);
    return 0;
}
View Code

另一方面,某一子串的出现次数就等于其 $\mathrm{endpos}$ 集合的大小。根据性质“每个状态的 $\mathrm{endpos}$ 集合为其在 parent 树中子树 $\mathrm{endpos}$ 集合的并集”,我们对于 $s[1:i]$ 所在状态赋为 $1$,其余状态为 $0$,然后通过在 parent 树上 dfs 求子树中 $1$ 的个数即可。

例题:洛谷P3804 【模板】后缀自动机 (SAM)

#include <cstdio>
#include <vector>
#include <cstring>
#include <algorithm>
using namespace std;

typedef long long ll;

const int N=1000005;
const int C=26; //注意检查字符集大小!

struct SuffixAutomaton
{
    int sz,lst; //状态数上限=2|S|-1 
    
    int len[N<<1],link[N<<1];
    int nxt[N<<1][C];
    //map<char,int> nxt[N<<1];
    //extend(char),并使用nxt[clone]=nxt[q]替换memcpy
    
    SuffixAutomaton()
    {
        len[0]=0,link[0]=-1;
        lst=sz=0;
    }
    
    void extend(int c)
    {
        int cur=++sz;
        len[cur]=len[lst]+1;
        
        int p=lst;
        while(p!=-1 && !nxt[p][c])
            nxt[p][c]=cur,p=link[p];
        
        if(p==-1)
            link[cur]=0;
        else
        {
            int q=nxt[p][c];
            if(len[p]+1==len[q])
                link[cur]=q;
            else
            {
                int clone=++sz;
                len[clone]=len[p]+1;
                link[clone]=link[q];
                memcpy(nxt[clone],nxt[q],C*4);
                
                while(p!=-1 && nxt[p][c]==q)
                    nxt[p][c]=clone,p=link[p];
                link[q]=link[cur]=clone;
            }
        }
        lst=cur;
    }
    
    int id[N<<1];
    
    //构建完成后,id顺序为len递增(逆拓扑序)【仅可排一次】 
    void sort()
    {
        static int bucket[N<<1];
        
        memset(bucket,0,sizeof(bucket));
        for(int i=1;i<=sz;i++)
            bucket[len[i]]++;
        for(int i=1;i<=sz;i++)
            bucket[i]+=bucket[i-1];
        for(int i=1;i<=sz;i++)
            id[bucket[len[i]]--]=i;
    }
}sam;

int n;
char s[N];

int cnt[N<<1];

int main()
{
    scanf("%s",s+1);
    n=strlen(s+1);
    
    for(int i=1;i<=n;i++)
        sam.extend(s[i]-'a'),cnt[sam.lst]=1;
    
    sam.sort();
    
    ll ans=0;
    for(int i=sam.sz;i>=1;i--)
    {
        int id=sam.id[i];
        if(cnt[id]>1)
            ans=max(ans,1LL*cnt[id]*sam.len[id]);
        
        cnt[sam.link[id]]+=cnt[id];
    }
    printf("%lld\n",ans);
    return 0;
}
View Code

在以上的实现中,并没有显式地对于 parent 树进行 dfs,而是将所有状态拓扑排序后按序处理。进行拓扑排序的策略并不复杂,我们只需要按照 $\mathrm{len}(state)$ 进行一趟桶排,这样即可保证祖先($\mathrm{len}$ 较小)一定排在子孙($\mathrm{len}$ 较大)的前面,此时倒序即为一个合法的拓扑排序。

int id[N<<1];

//构建完成后,id顺序为len递增(逆拓扑序)【仅可排一次】 
void sort()
{
    static int bucket[N<<1];
    
    memset(bucket,0,sizeof(bucket));
    for(int i=1;i<=sz;i++)
        bucket[len[i]]++;
    for(int i=1;i<=sz;i++)
        bucket[i]+=bucket[i-1];
    for(int i=1;i<=sz;i++)
        id[bucket[len[i]]--]=i;
}

 

2. 字典序第 K 大子串(SAM 上 DP)

我们先考虑本质相同的子串仅计算一次的情况。刚刚在计算本质不同子串数量的时候,我们利用状态的性质解决了,但实际上还有一个在 SAM 上 DP 的做法。

有这样的一个结论:一个字符串本质不同的子串数量,相当于在 SAM 上通过 $\mathrm{nxt}$ 转移产生的路径数量。我们尝试理解这个结论。首先,这样的计数一定不会少算(有子串没被算到),因为每个本质不同的子串都会在 SAM 上存在唯一对应的路径;同时,这个计数也不会多算(算到多余的子串),因为 SAM 仅接受模式串的子串。

利用这个结论,我们直接利用 DP 计算从每个状态出发的路径数量即可(由于 $\mathrm{nxt}$ 转移构成的图是一个 DAG,我们仍然可以通过令 $\mathrm{len}$ 从小到大来获得一个正拓扑序)。DP 的转移方程为 $path[x]=1+\sum_{y=\mathrm{nxt}_x(c_i)}path[y]$,需要按照逆拓扑序依次计算。

对于本质相同的子串计算出现次数的情况,将上面转移方程中的 $1$ 替换成 $cnt[x]$ 即可。

那么求字典序第 $K$ 大子串,只需要从 $state_0$ 开始贪心,每次选取最小的前缀路径数量大于等于 $K$ 的转移即可(原理同线段树上二分,只不过 SAM 上每个状态有字符集大小的分叉)。

例题:LOJ2102 「TJOI2015」弦论

#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;

typedef long long ll;

const int N=500005;
const int C=26; //注意检查字符集大小!

struct SuffixAutomaton
{
    int sz,lst; //状态数上限=2|S|-1 
    
    int len[N<<1],link[N<<1];
    int nxt[N<<1][C];
    //map<char,int> nxt[N<<1];
    //extend(char),并使用nxt[clone]=nxt[q]替换memcpy
    
    SuffixAutomaton()
    {
        len[0]=0,link[0]=-1;
        lst=sz=0;
    }
    
    void extend(int c)
    {
        int cur=++sz;
        len[cur]=len[lst]+1;
        
        int p=lst;
        while(p!=-1 && !nxt[p][c])
            nxt[p][c]=cur,p=link[p];
        
        if(p==-1)
            link[cur]=0;
        else
        {
            int q=nxt[p][c];
            if(len[p]+1==len[q])
                link[cur]=q;
            else
            {
                int clone=++sz;
                len[clone]=len[p]+1;
                link[clone]=link[q];
                memcpy(nxt[clone],nxt[q],C*4);
                
                while(p!=-1 && nxt[p][c]==q)
                    nxt[p][c]=clone,p=link[p];
                link[q]=link[cur]=clone;
            }
        }
        lst=cur;
    }
    
    int id[N<<1];
    
    //构建完成后,id顺序为len递增(逆拓扑序)【仅可排一次】 
    void sort()
    {
        static int bucket[N<<1];
        
        memset(bucket,0,sizeof(bucket));
        for(int i=1;i<=sz;i++)
            bucket[len[i]]++;
        for(int i=1;i<=sz;i++)
            bucket[i]+=bucket[i-1];
        for(int i=1;i<=sz;i++)
            id[bucket[len[i]]--]=i;
    }
}sam;

int n;
char s[N];
int type,K;

int cnt[N<<1];
ll path[N<<1][2];

int main()
{
    scanf("%s",s+1);
    scanf("%d%d",&type,&K);
    n=strlen(s+1);
    
    for(int i=1;i<=n;i++)
        sam.extend(s[i]-'a'),cnt[sam.lst]=1;
    
    sam.sort();
    
    for(int i=sam.sz;i>=0;i--)
    {
        int id=sam.id[i];
        
        path[id][0]=1,path[id][1]=cnt[id];
        for(int j=0;j<C;j++)
        {
            int nxt=sam.nxt[id][j];
            if(nxt)
                for(int k=0;k<2;k++)
                    path[id][k]+=path[nxt][k];
        }
        if(id!=-1)
            cnt[sam.link[id]]+=cnt[id];
    }
    
    if(path[0][type]-cnt[0]<K)
        printf("-1\n");
    else
    {
        int cur=0;
        while(K)
        {
            for(int i=0;i<C;i++)
            {
                int nxt=sam.nxt[cur][i];
                if(!nxt)
                    continue;
                
                if(path[nxt][type]>=K)
                {
                    cur=nxt,putchar('a'+i);
                    K-=(type?cnt[cur]:1);
                    break;
                }
                K-=path[nxt][type];
            }
        }
    }
    return 0;
}
View Code

 

3. 两个字符串的最长公共子串(SAM 上匹配)

对于 $s,t$ 两个字符串,对于 $s$ 串建立 SAM,$t$ 在这个 SAM 上进行匹配,失配时跳 $\mathrm{link}$。实现上和 AC 自动机差不多。

例题:SPOJ-LCS Longest Common Substring

#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;

typedef long long ll;

const int N=250005;
const int C=26; //注意检查字符集大小!

struct SuffixAutomaton
{
    int sz,lst; //状态数上限=2|S|-1 
    
    int len[N<<1],link[N<<1];
    int nxt[N<<1][C];
    //map<char,int> nxt[N<<1];
    //extend(char),并使用nxt[clone]=nxt[q]替换memcpy
    
    SuffixAutomaton()
    {
        len[0]=0,link[0]=-1;
        lst=sz=0;
    }
    
    void extend(int c)
    {
        int cur=++sz;
        len[cur]=len[lst]+1;
        
        int p=lst;
        while(p!=-1 && !nxt[p][c])
            nxt[p][c]=cur,p=link[p];
        
        if(p==-1)
            link[cur]=0;
        else
        {
            int q=nxt[p][c];
            if(len[p]+1==len[q])
                link[cur]=q;
            else
            {
                int clone=++sz;
                len[clone]=len[p]+1;
                link[clone]=link[q];
                memcpy(nxt[clone],nxt[q],C*4);
                
                while(p!=-1 && nxt[p][c]==q)
                    nxt[p][c]=clone,p=link[p];
                link[q]=link[cur]=clone;
            }
        }
        lst=cur;
    }
    
    int id[N<<1];
    
    //构建完成后,id顺序为len递增(逆拓扑序)【仅可排一次】 
    void sort()
    {
        static int bucket[N<<1];
        
        memset(bucket,0,sizeof(bucket));
        for(int i=1;i<=sz;i++)
            bucket[len[i]]++;
        for(int i=1;i<=sz;i++)
            bucket[i]+=bucket[i-1];
        for(int i=1;i<=sz;i++)
            id[bucket[len[i]]--]=i;
    }
}sam;

int n,m;
char s[N],t[N];

int main()
{
    scanf("%s%s",s+1,t+1);
    n=strlen(s+1),m=strlen(t+1);
    
    for(int i=1;i<=n;i++)
        sam.extend(s[i]-'a');
    
    int ans=0,cur=0,len=0;
    for(int i=1;i<=m;i++)
    {
        int p=cur;
        while(p!=-1 && !sam.nxt[p][t[i]-'a'])
            p=sam.link[p],len=sam.len[p];
        
        if(p!=-1)
            cur=sam.nxt[p][t[i]-'a'],len++;
        else
            cur=len=0;
        
        ans=max(ans,len);
    }
    printf("%d\n",ans);
    return 0;
}
View Code

 

4. 多个字符串的最长公共子串(SAM 上匹配 + parent 树上 DP)

想同时处理所有字符串是比较困难的,我们考虑先对于一个串建立 SAM,然后对于其余的串依次求出“当匹配到状态 $state$ 时,最长的匹配长度”。

考虑对于具体某一个待匹配串 $t$,我们还是像两个字符串一样求 $t[1:i]$ 在 SAM 上匹配得到的状态以及匹配上的长度 $len$。但不同的是,我们还需要额外维护一个数组 $curlen$,表示字符串 $t$ 对于每个状态的 $\mathrm{longest}(state)$ 能匹配的最长后缀长度。更新的规则为 $curlen[x]=\mathrm{max}(curlen[x],len)$,因为 $t[1:i]$ 所对应的状态可能由多个状态转移而来,我们需要取其中最大的匹配长度。

需要注意到,我们在维护 $curlen$ 的过程中,只考虑到了 $t[1:i]$ 这不超过 $|t|$ 个状态,但实际上对于 SAM 上的任意状态均应计算 $curlen$。所以我们按照 $\mathrm{len}$ 从大到小的顺序计算没经过的状态的 $curlen$,并使用 $curlen[x]$ 来更新 $curlen[\mathrm{link}(x)]$,规则为 $curlen[\mathrm{link}(x)]=\mathrm{max}(curlen[\mathrm{link}(x)],curlen[x])$。

对于一个待匹配串 $t$ 我们能按照如上的方法求得每个状态的 $curlen$,那么所有字符串的在每个状态的最长匹配长度,就是 $ans[state]=\mathrm{min}_t(curlen_t[state])$。而最终答案为 $\mathrm{min}_{state}(ans[state])$。

例题:SPOJ-LCS2 Longest Common Substring II

#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;

typedef long long ll;

const int N=100005;
const int C=26; //注意检查字符集大小!

struct SuffixAutomaton
{
    int sz,lst; //状态数上限=2|S|-1 
    
    int len[N<<1],link[N<<1];
    int nxt[N<<1][C];
    //map<char,int> nxt[N<<1];
    //使用nxt[clone]=nxt[q]替换memcpy
    
    SuffixAutomaton()
    {
        len[0]=0,link[0]=-1;
        lst=0;
    }
    
    void extend(int c)
    {
        int cur=++sz;
        len[cur]=len[lst]+1;
        
        int p=lst;
        while(p!=-1 && !nxt[p][c])
            nxt[p][c]=cur,p=link[p];
        
        if(p==-1)
            link[cur]=0;
        else
        {
            int q=nxt[p][c];
            if(len[p]+1==len[q])
                link[cur]=q;
            else
            {
                int clone=++sz;
                len[clone]=len[p]+1;
                link[clone]=link[q];
                memcpy(nxt[clone],nxt[q],C*4);
                
                while(p!=-1 && nxt[p][c]==q)
                    nxt[p][c]=clone,p=link[p];
                link[q]=link[cur]=clone;
            }
        }
        lst=cur;
    }
    
    int id[N<<1];
    
    //构建完成后,id顺序为len递增(逆拓扑序)【仅可排一次】 
    void sort()
    {
        static int bucket[N<<1];
        
        memset(bucket,0,sizeof(bucket));
        for(int i=1;i<=sz;i++)
            bucket[len[i]]++;
        for(int i=1;i<=sz;i++)
            bucket[i]+=bucket[i-1];
        for(int i=1;i<=sz;i++)
            id[bucket[len[i]]--]=i;
    }
}sam;

int n,m;
char s[N],t[N];

int ans[N<<1],curlen[N<<1];

int main()
{
    scanf("%s",s+1);
    n=strlen(s+1);
    
    for(int i=1;i<=n;i++)
        sam.extend(s[i]-'a');
    sam.sort();
    
    for(int i=1;i<=sam.sz;i++)
        ans[i]=sam.len[i];
    
    while(~scanf("%s",t+1))
    {
        m=strlen(t+1);
        
        int cur=0,len=0;
        for(int i=1;i<=m;i++)
        {
            int p=cur,c=t[i]-'a';
            while(p!=-1 && !sam.nxt[p][c])
                p=sam.link[p],len=sam.len[p];
            
            if(p!=-1)
                cur=sam.nxt[p][c],len++;
            else
                cur=len=0;
            
            //对于同一cur,可能从多个p转移过来,故不能用len直接赋值cur
            curlen[pos]=max(curlen[pos],len);
        }
        
        for(int i=sam.sz;i>=1;i--)
        {
            int id=sam.id[i];
            curlen[sam.link[id]]=max(curlen[sam.link[id]],curlen[id]);
        }
        for(int i=1;i<=sam.sz;i++)
            ans[i]=min(ans[i],curlen[i]),curlen[i]=0;
    }
    
    int fans=0;
    for(int i=1;i<=sam.sz;i++)
        fans=max(fans,ans[i]);
    printf("%d\n",fans);
    return 0;
}
View Code

 

 

~ 从 SAM 到广义 SAM ~

 

在上文中处理多个字符串的问题时,我们先对于其中一个字符串建立 SAM,对于剩余的字符串依次在这个 SAM 上面跑。有没有办法对于多个字符串直接建立一个自动机,使得任意字符串的任意子串都能被自动机接受?答案是肯定的。

我们考虑这样的实现:对于字符串 $s_1,\cdots,s_n$,我们首先对于 $s_1$ 使用 $extend()$ 函数依次插入建立 SAM,接着我们将 $last$ 移回 $state_0$ 并重新依次插入 $s_2$。这个实现的逻辑比较简单,对于能够“接受任意字符串的任意子串”的正确性也是比较显然的,不过在细节上存在一些问题:对于多次经过同一个状态时反复插入会产生孤立联通块,这些孤立连通块内的点在进行从 $state_0$ 开始的匹配时永远不会被走到(详见 ix35 - 悲惨故事 长文警告 关于广义 SAM 的讨论)。这种现象的存在会导致进行拓扑排序等操作时出错(一个例子为 yy1695651 - ICPC2019 G题广义SAM的基排为什么不行?)。

另一方面,上述方法的复杂度与 $\sum_i |s_i|$ 线性相关,因为相当于依次插入 $n$ 个字符串。但假如题目中的字符串不是直接给出的,而是采用 trie 树的形式,那么 $n$ 个节点的 trie 树可以产生总长为 $n^2$ 级别的字符串,其形状大概如下:

出于适用范围更广的原则,本文介绍基于 trie 树和离线 bfs 的广义 SAM 实现方法。其能够建立正确的广义 SAM,同时时空复杂度仅与 trie 树的节点数相关。

 

 

~ 广义 SAM(GSA)的建立 ~

 

我们从 trie 树的角度来考虑 SAM 的实现。由于模式串仅为单一字符串,所以对于 SAM 而言是基于退化为一条链的 trie 树建立的。

那么对于一个非退化的 trie 树,我们也可以在修改很少几处的情况下建立 GSA。

 

1. 插入顺序

在建立 SAM 的过程中,我们是按照字符串中从前到后的顺序依次扩展 SAM 的,从 trie 树的角度看就是从浅到深。那么在建立 GSA 时,我们为了保证先扩展较浅的节点、再扩展较深的节点,可以采用 bfs 进行实现。

如果不采用从浅到深的顺序,就可能出现由状态 $x$ 转移到状态 $y$、而状态 $y$ 已经先求完 $\mathrm{nxt}$ 与 $\mathrm{link}$ 的情况,此时会产生类似每次重置 $last$ 一样的孤立连通块问题,需要额外的不少讨论才能保证正确性。使用从浅到深的 bfs 则不会出现这种情况。

 

2. 节点标号

在建立 GSA 前,我们首先需要用 $\mathrm{nxt}$ 数组来建立 trie 树。因此每步插入的 $cur$ 都是已经建好的节点,无须开新点;只有 $clone$ 需要动态开点。同时,由于我们按照深度扩展,每次的 $last$ 都是通过 bfs 得知的,无须再额外维护 $last$ 这个变量。

 

3. 信息隔离

在建立 SAM 的过程中,我们在插入 $s[i]$ 时不会知道任何 $s[j]$($j>i$)的信息。但是在建立 GSA 时,我们有可能获得更深节点的信息,具体来说就是 $\mathrm{nxt}_{last}$ 可能指向比 $last$ 更深的节点。因此我们需要仔细检查之前 SAM 的板子中所有用到 $\mathrm{nxt}$ 的地方:

void extend(int c)
{
    int cur=++sz;
    len[cur]=len[lst]+1;
    
    int p=lst;
    //nxt[lst][c]一定存在,且等于cur【会产生问题 #1】 
    while(p!=-1 && !nxt[p][c])
        nxt[p][c]=cur,p=link[p];
    
    if(p==-1)
        link[cur]=0;
    else
    {
        int q=nxt[p][c]; 
        //这里我们认为dep[p]<dep[lst],
        //那么假如nxt[p][c]是某次clone赋值的,说明那次clone在当前扩展操作之前,那么直接使用是ok的; 
        //假如nxt[p][c]是构建trie树时赋值的,则dep[q]=dep[p]+1<=dep[lst]<dep[cur],说明一定已经被扩展到了。
        //综上,【不会产生问题】 
        if(len[p]+1==len[q])
            link[cur]=q;
        else
        {
            int clone=++sz;
            len[clone]=len[p]+1;
            link[clone]=link[q];
            memcpy(nxt[clone],nxt[q],C*4);
            //假如nxt[q][i]是构建trie树时赋值的,dep[nxt[q][i]]最大可以等于dep[cur], 
            //即使nxt[q][i]与cur同深度,在广义SAM的扩展过程中也有先后之分,
            //所以如果扩展cur时未扩展过nxt[q][i],那么nxt[q][i]在此时应该视为0【会产生问题 #2】 
            
            while(p!=-1 && nxt[p][c]==q)
                nxt[p][c]=clone,p=link[p];
            //同q=nxt[p][c]处的分析,对于nxt[p][c]的使用都是安全的【不会产生问题】 
            link[q]=link[cur]=clone;
        }
    }
    lst=cur;
}

我们需要对于会产生问题的两处稍加修改,具体而言:

  (1) 在 #1 处,初始令 $p=\mathrm{link}[lst]$ 即可(因为 $\mathrm{nxt}[lst][c]=cur$ 在构建 trie 树时就被赋值了;而其余所有 $\mathrm{nxt}[p][c]$ 假如为非 clone 点,在 trie 树上的深度一定小于 $dep[cur]$,则可以直接使用 );

  (2) 在 #2 处,我们对于 $q$ 的每一个 $\mathrm{nxt}$ 转移,均检查是否被扩展过(等价于检查 $\mathrm{len}$ 是否为 $0$),若未被扩展过,即使 $\mathrm{nxt}[q][i]\neq 0$ 也必须视为 $0$。

 

综合以上三点,我们可以写出如下的 GSA 的实现:

typedef pair<int,int> pii;

const int N=1000005;
const int C=26; //注意检查字符集大小!

//在结构题外开任何与SAM状态相关的数组,都需要N<<1
struct GeneralSuffixAutomaton
{
    int sz; //状态数上限=2|S|-1 
    
    int len[N<<1],link[N<<1];
    int nxt[N<<1][C];
    //map<char,int> nxt[N<<1];
    //extend(int,char)
    
    GeneralSuffixAutomaton()
    {
        len[0]=0,link[0]=-1;
        sz=0;
    }
    
    //向trie树上插入一个字符串 
    void insert(char *s,int n)
    {
        for(int i=1,cur=0;i<=n;i++)
        {
            int c=s[i]-'a';
            if(!nxt[cur][c])
                nxt[cur][c]=++sz;
            cur=nxt[cur][c];
        }
    }
    
    //扩展GSA,对于lst的c转移进行扩展 
    int extend(int lst,int c)
    {
        int cur=nxt[lst][c];
        len[cur]=len[lst]+1;
        
        int p=link[lst];
        while(p!=-1 && !nxt[p][c])
            nxt[p][c]=cur,p=link[p];
        
        if(p==-1)
            link[cur]=0;
        else
        {
            int q=nxt[p][c];
            if(len[p]+1==len[q])
                link[cur]=q;
            else
            {
                int clone=++sz;
                len[clone]=len[p]+1;
                link[clone]=link[q];
                for(int i=0;i<C;i++)
                    nxt[clone][i]=len[nxt[q][i]]?nxt[q][i]:0;
                
                while(p!=-1 && nxt[p][c]==q)
                    nxt[p][c]=clone,p=link[p];
                link[q]=link[cur]=clone;
            }
        }
        return cur;
    }
    
    //基于trie树建立GSA
    void build()
    {
        queue<pii> Q;
        for(int i=0;i<C;i++)
            if(nxt[0][i])
                Q.push(pii(0,i));
        
        while(!Q.empty())
        {
            pii tmp=Q.front();
            Q.pop();
            
            int cur=extend(tmp.first,tmp.second);
            for(int i=0;i<C;i++)
                if(nxt[cur][i])
                    Q.push(pii(cur,i));
        }
    } 
    
    int id[N<<1];
    
    //构建完成后,id顺序为len递增(逆拓扑序)【仅可排一次】 
    void sort()
    {
        static int bucket[N<<1];
        
        memset(bucket,0,sizeof(bucket));
        for(int i=1;i<=sz;i++)
            bucket[len[i]]++;
        for(int i=1;i<=sz;i++)
            bucket[i]+=bucket[i-1];
        for(int i=1;i<=sz;i++)
            id[bucket[len[i]]--]=i;
    }
}gsa;

 

 

~ 经典应用(GSA) ~

 

1. 本质不同的子串数量及出现次数

本质不同的子串数量还是利用 GSA 的性质即可,等于 $\sum_{state}\mathrm{len}(state)-\mathrm{len}\left(\mathrm{link}(state)\right)$。

例题:洛谷P6139 【模板】广义后缀自动机(广义 SAM)

#include <queue>
#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;

typedef long long ll;
typedef pair<int,int> pii;

const int N=1000005;
const int C=26; //注意检查字符集大小!

//在结构题外开任何与SAM状态相关的数组,都需要N<<1
struct GeneralSuffixAutomaton
{
    int sz; //状态数上限=2|S|-1 
    
    int len[N<<1],link[N<<1];
    int nxt[N<<1][C];
    //map<char,int> nxt[N<<1];
    //extend(int,char)
    
    GeneralSuffixAutomaton()
    {
        len[0]=0,link[0]=-1;
        sz=0;
    }
    
    //向trie树上插入一个字符串 
    void insert(char *s,int n)
    {
        for(int i=1,cur=0;i<=n;i++)
        {
            int c=s[i]-'a';
            if(!nxt[cur][c])
                nxt[cur][c]=++sz;
            cur=nxt[cur][c];
        }
    }
    
    //扩展GSA,对于lst的c转移进行扩展 
    int extend(int lst,int c)
    {
        int cur=nxt[lst][c];
        len[cur]=len[lst]+1;
        
        int p=link[lst];
        while(p!=-1 && !nxt[p][c])
            nxt[p][c]=cur,p=link[p];
        
        if(p==-1)
            link[cur]=0;
        else
        {
            int q=nxt[p][c];
            if(len[p]+1==len[q])
                link[cur]=q;
            else
            {
                int clone=++sz;
                len[clone]=len[p]+1;
                link[clone]=link[q];
                for(int i=0;i<C;i++)
                    nxt[clone][i]=len[nxt[q][i]]?nxt[q][i]:0;
                
                while(p!=-1 && nxt[p][c]==q)
                    nxt[p][c]=clone,p=link[p];
                link[q]=link[cur]=clone;
            }
        }
        return cur;
    }
    
    //基于trie树建立GSA
    void build()
    {
        queue<pii> Q;
        for(int i=0;i<C;i++)
            if(nxt[0][i])
                Q.push(pii(0,i));
        
        while(!Q.empty())
        {
            pii tmp=Q.front();
            Q.pop();
            
            int cur=extend(tmp.first,tmp.second);
            for(int i=0;i<C;i++)
                if(nxt[cur][i])
                    Q.push(pii(cur,i));
        }
    } 
    
    int id[N<<1];
    
    //构建完成后,id顺序为len递增(逆拓扑序)【仅可排一次】 
    void sort()
    {
        static int bucket[N<<1];
        
        memset(bucket,0,sizeof(bucket));
        for(int i=1;i<=sz;i++)
            bucket[len[i]]++;
        for(int i=1;i<=sz;i++)
            bucket[i]+=bucket[i-1];
        for(int i=1;i<=sz;i++)
            id[bucket[len[i]]--]=i;
    }
}gsa;

int n,m;
char s[N];

int main()
{
    scanf("%d",&n);
    for(int i=1;i<=n;i++)
    {
        scanf("%s",s+1);
        m=strlen(s+1);
        gsa.insert(s,m);
    }
    
    gsa.build();
    
    ll ans=0;
    for(int i=1;i<=gsa.sz;i++)
        ans+=gsa.len[i]-gsa.len[gsa.link[i]];
    printf("%lld\n",ans);
    return 0;
}
View Code

出现次数也是在 parent 树上 DP。

例题:洛谷P3181 [HAOI2016]找相同字符

这道题中,具体一个状态产生的贡献是 $\left(\mathrm{len}(state)-\mathrm{minlen}(state)+1\right)\cdot cnt_1[state]\cdot cnt_2[state]$。

#include <queue>
#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;

typedef long long ll;
typedef pair<int,int> pii;

const int N=400005;
const int C=26; //注意检查字符集大小!

//在结构题外开任何与SAM状态相关的数组,都需要N<<1
struct GeneralSuffixAutomaton
{
    int sz; //状态数上限=2|S|-1 
    
    int len[N<<1],link[N<<1];
    int nxt[N<<1][C];
    //map<char,int> nxt[N<<1];
    //extend(int,char)
    
    GeneralSuffixAutomaton()
    {
        len[0]=0,link[0]=-1;
        sz=0;
    }
    
    //向trie树上插入一个字符串 
    void insert(char *s,int n)
    {
        for(int i=1,cur=0;i<=n;i++)
        {
            int c=s[i]-'a';
            if(!nxt[cur][c])
                nxt[cur][c]=++sz;
            cur=nxt[cur][c];
        }
    }
    
    //扩展GSA,对于lst的c转移进行扩展 
    int extend(int lst,int c)
    {
        int cur=nxt[lst][c];
        len[cur]=len[lst]+1;
        
        int p=link[lst];
        while(p!=-1 && !nxt[p][c])
            nxt[p][c]=cur,p=link[p];
        
        if(p==-1)
            link[cur]=0;
        else
        {
            int q=nxt[p][c];
            if(len[p]+1==len[q])
                link[cur]=q;
            else
            {
                int clone=++sz;
                len[clone]=len[p]+1;
                link[clone]=link[q];
                for(int i=0;i<C;i++)
                    nxt[clone][i]=len[nxt[q][i]]?nxt[q][i]:0;
                
                while(p!=-1 && nxt[p][c]==q)
                    nxt[p][c]=clone,p=link[p];
                link[q]=link[cur]=clone;
            }
        }
        return cur;
    }
    
    //基于trie树建立GSA
    void build()
    {
        queue<pii> Q;
        for(int i=0;i<C;i++)
            if(nxt[0][i])
                Q.push(pii(0,i));
        
        while(!Q.empty())
        {
            pii tmp=Q.front();
            Q.pop();
            
            int cur=extend(tmp.first,tmp.second);
            for(int i=0;i<C;i++)
                if(nxt[cur][i])
                    Q.push(pii(cur,i));
        }
    } 
    
    int id[N<<1];
    
    //构建完成后,id顺序为len递增(逆拓扑序)【仅可排一次】 
    void sort()
    {
        static int bucket[N<<1];
        
        memset(bucket,0,sizeof(bucket));
        for(int i=1;i<=sz;i++)
            bucket[len[i]]++;
        for(int i=1;i<=sz;i++)
            bucket[i]+=bucket[i-1];
        for(int i=1;i<=sz;i++)
            id[bucket[len[i]]--]=i;
    }
}gsa;

int n,m;
char s[N],t[N];
int cnt1[N<<1],cnt2[N<<1];

int main()
{
    scanf("%s",s+1);
    n=strlen(s+1);
    scanf("%s",t+1);
    m=strlen(t+1);
    
    gsa.insert(s,n);
    gsa.insert(t,m);
    
    gsa.build();
    gsa.sort();
    
    for(int i=1,cur=0;i<=n;i++)
    {
        cur=gsa.nxt[cur][s[i]-'a'];
        cnt1[cur]++;
    }
    for(int i=1,cur=0;i<=m;i++)
    {
        cur=gsa.nxt[cur][t[i]-'a'];
        cnt2[cur]++;
    }
    
    ll ans=0;
    for(int i=gsa.sz;i>=1;i--)
    {
        int id=gsa.id[i];
        int lst=gsa.link[id];
        ans+=1LL*(gsa.len[id]-gsa.len[lst])*cnt1[id]*cnt2[id];
        cnt1[lst]+=cnt1[id];
        cnt2[lst]+=cnt2[id];
    }
    printf("%lld\n",ans);
    return 0;
}
View Code

 

多个字符串的最长公共子串也许不能用 GSA 做

SPOJ 的那题是比较弱的,因为 $n\leq 10$,所以可以存个长度为 $10$ 的数组搞一搞。以下讨论的是 $n, \sum_i |s_i|\leq 1\times 10^6$ 级别的该问题。

在 SAM 上,疑似对于每个状态暴力跳 fail 是 $\mathrm{O}(n^{1.5})$ 的。那么在 GSA 上,是否能够保证所有串各自占有的状态之和是 $|s_i|^{1.5}$ 呢?这个结论是存疑的...就算成立,时间复杂度也明显比 SAM 要差。

 

 

~ 一些练习题 ~

 

SAM / GSA 的裸题一般都是一眼会,稍微难一点的题目会在此基础上套线段树合并LCT access 维护树链等内容,但此时难点就不再 SAM 上了。做到啥补充啥吧。

 

Kattis-firstofhername First of Her Name(2019 ICPC WF)【GSA 裸题】

题目中给定了一棵 trie 树,所以建立 GSA 是比较自然的想法。考虑回答查询,某个字符串的前缀如果是 $s$,就相当于反转的 $s$ 是从 trie 树的对应位置到 $state_0$ 的路径的前缀。所以我们只需要用类似求出现次数的方法,对于每个字符串在结尾处令 $cnt=1$(在这题中就是将 $cnt[1:n]$ 均置 $1$),然后按照拓扑序 DP 即可。

#include <queue>
#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;

typedef long long ll;
typedef pair<int,int> pii;

const int N=1000005;
const int C=26; //注意检查字符集大小!

ll cnt[N<<1];

//在结构题外开任何与SAM状态相关的数组,都需要N<<1
struct GeneralSuffixAutomaton
{
    int sz; //状态数上限=2|S|-1 
    
    int len[N<<1],link[N<<1];
    int nxt[N<<1][C];
    //map<char,int> nxt[N<<1];
    //extend(int,char)
    
    GeneralSuffixAutomaton()
    {
        len[0]=0,link[0]=-1;
        sz=0;
    }
    
    //扩展GSA,对于lst的c转移进行扩展 
    int extend(int lst,int c)
    {
        int cur=nxt[lst][c];
        len[cur]=len[lst]+1;
        
        int p=link[lst];
        while(p!=-1 && !nxt[p][c])
            nxt[p][c]=cur,p=link[p];
        
        if(p==-1)
            link[cur]=0;
        else
        {
            int q=nxt[p][c];
            if(len[p]+1==len[q])
                link[cur]=q;
            else
            {
                int clone=++sz;
                len[clone]=len[p]+1;
                link[clone]=link[q];
                for(int i=0;i<C;i++)
                    nxt[clone][i]=len[nxt[q][i]]?nxt[q][i]:0;
                
                while(p!=-1 && nxt[p][c]==q)
                    nxt[p][c]=clone,p=link[p];
                link[q]=link[cur]=clone;
            }
        }
        return cur;
    }
    
    //基于trie树建立GSA
    void build()
    {
        queue<pii> Q;
        for(int i=0;i<C;i++)
            if(nxt[0][i])
                Q.push(pii(0,i));
        
        while(!Q.empty())
        {
            pii tmp=Q.front();
            Q.pop();
            
            int cur=extend(tmp.first,tmp.second);
            for(int i=0;i<C;i++)
                if(nxt[cur][i])
                    Q.push(pii(cur,i));
        }
    } 
    
    int id[N<<1];
    
    //构建完成后,id顺序为len递增(逆拓扑序)【仅可排一次】 
    void sort()
    {
        static int bucket[N<<1];
        
        memset(bucket,0,sizeof(bucket));
        for(int i=1;i<=sz;i++)
            bucket[len[i]]++;
        for(int i=1;i<=sz;i++)
            bucket[i]+=bucket[i-1];
        for(int i=1;i<=sz;i++)
            id[bucket[len[i]]--]=i;
    }
}gsa;

int n,q;
char s[N];

int main()
{
    scanf("%d%d",&n,&q);
    for(int i=1;i<=n;i++)
    {
        char ch=getchar();
        while(ch<'A' || ch>'Z')
            ch=getchar();
        
        int lst;
        scanf("%d",&lst);
        gsa.nxt[lst][ch-'A']=i,cnt[i]=1;
    }
    
    gsa.sz=n;
    gsa.build();
    gsa.sort();
    
    for(int i=gsa.sz;i>=1;i--)
    {
        int id=gsa.id[i];
        cnt[gsa.link[id]]+=cnt[id];
    }
    
    while(q--)
    {
        scanf("%s",s+1);
        int n=strlen(s+1);
        
        int cur=0;
        for(int i=n;i>=1;i--)
        {
            int c=s[i]-'A';
            if(!gsa.nxt[cur][c])
            {
                cur=0;
                break;
            }
            cur=gsa.nxt[cur][c];
        }
        printf("%d\n",gsa.len[cur]>=n?cnt[cur]:0);
    }
    return 0;
}
View Code

 

计蒜客42551 Loli, Yen-Jen, and a cool problem(2019 ICPC 徐州)【GSA 基础性质,去重处理】

这题本身挺裸的。对于每次询问 $s[x-L+1:x]$ 的出现次数,我们只需要找到其所在的状态就行了,从 $x$ 沿着 $\mathrm{link}$ 暴力跳、直到当前状态的 $\mathrm{len}$ 区间包含 $L$ 即可。理论上需要将 $x$ 相同的询问同时处理,不过数据貌似没卡。

这题稍微有点恶心的是 trie 树上会有重合的节点。对于重合的节点,我们选择最早出现的一个作为代表,其余节点全部映射到代表。这会导致非代表点变成孤立连通块,在进行拓扑排序的时候需要注意忽略这些 $\mathrm{len}=0$ 的孤立点。

#include <queue>
#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;

typedef long long ll;
typedef pair<int,int> pii;

const int N=300005;
const int C=26; //注意检查字符集大小!

//在结构题外开任何与SAM状态相关的数组,都需要N<<1
struct GeneralSuffixAutomaton
{
    int sz; //状态数上限=2|S|-1 
    
    int len[N<<1],link[N<<1];
    int nxt[N<<1][C];
    //map<char,int> nxt[N<<1];
    //extend(int,char)
    
    GeneralSuffixAutomaton()
    {
        len[0]=0,link[0]=-1;
        sz=0;
    }
    
    //扩展GSA,对于lst的c转移进行扩展 
    int extend(int lst,int c)
    {
        int cur=nxt[lst][c];
        len[cur]=len[lst]+1;
        
        int p=link[lst];
        while(p!=-1 && !nxt[p][c])
            nxt[p][c]=cur,p=link[p];
        
        if(p==-1)
            link[cur]=0;
        else
        {
            int q=nxt[p][c];
            if(len[p]+1==len[q])
                link[cur]=q;
            else
            {
                int clone=++sz;
                len[clone]=len[p]+1;
                link[clone]=link[q];
                for(int i=0;i<C;i++)
                    nxt[clone][i]=len[nxt[q][i]]?nxt[q][i]:0;
                
                while(p!=-1 && nxt[p][c]==q)
                    nxt[p][c]=clone,p=link[p];
                link[q]=link[cur]=clone;
            }
        }
        return cur;
    }
    
    //基于trie树建立GSA
    void build()
    {
        queue<pii> Q;
        for(int i=0;i<C;i++)
            if(nxt[0][i])
                Q.push(pii(0,i));
        
        while(!Q.empty())
        {
            pii tmp=Q.front();
            Q.pop();
            
            int cur=extend(tmp.first,tmp.second);
            for(int i=0;i<C;i++)
                if(nxt[cur][i])
                    Q.push(pii(cur,i));
        }
    } 
    
    int id[N<<1];
    
    //构建完成后,id顺序为len递增(逆拓扑序)【仅可排一次】 
    void sort()
    {
        static int bucket[N<<1];
        
        memset(bucket,0,sizeof(bucket));
        for(int i=1;i<=sz;i++)
            if(len[i])
                bucket[len[i]]++;
        for(int i=1;i<=sz;i++)
            bucket[i]+=bucket[i-1];
        for(int i=1;i<=sz;i++)
            if(len[i])
                id[bucket[len[i]]--]=i;
    }
}gsa;

int n,q;
char s[N];

int to[N];
ll cnt[N<<1];

int main()
{
    scanf("%d%d",&n,&q);
    scanf("%s",s+1);
    
    gsa.nxt[0][s[1]-'A']=1,to[1]=1,cnt[1]++;
    for(int i=2;i<=n;i++)
    {
        int x,c=s[i]-'A';
        scanf("%d",&x);
        
        x=to[x];
        if(!gsa.nxt[x][c])
            gsa.nxt[x][c]=i;
        to[i]=gsa.nxt[x][c];
        cnt[to[i]]++;
    }
    
    gsa.sz=n;
    gsa.build();
    gsa.sort();
    
    for(int i=gsa.sz;i>=1;i--)
    {
        int id=gsa.id[i];
        if(id)
            cnt[gsa.link[id]]+=cnt[id];
    }
    
    while(q--)
    {
        int x,y;
        scanf("%d%d",&x,&y);
        
        x=to[x];
        while(gsa.len[gsa.link[x]]>=y)
            x=gsa.link[x];
        printf("%lld\n",cnt[x]);
    }
    return 0;
}
View Code

 

牛客4010D 卡拉巴什的字符串2020 Wannafly Camp Day2)【题意转化,SAM 基础性质】

我们考虑依次加入字符,那么在加到 $s[i]$ 时,我们需要考虑 $s[1:i-1]$ 的答案与 $s[1:i]$ 的区别。假如加入 $s[i]$ 会使得答案发生改变,那么最长的 LCP 一定是一个 $s[1:i]$ 的后缀,否则该 LCP 也一定在 $s[1:i-1]$ 中出现过了。

通过这样的转化后,我们只需要关心最长的、出现过多次的 $s[1:i]$ 的后缀是什么。在插入完 $s[i]$ 后,SAM 的 $last$ 指针指向 $\mathrm{endpos}=\{i\}$ 的状态,那么根据后缀链接的性质,$\mathrm{link}(last)$ 恰为我们所要的东西(因为其 endpos 集合大于 $1$,即出现多次),所以答案即为 $\mathrm{len}(\mathrm{link}(last))$ 的 MEX。有一个细节是,当字符串为形如 "aaaa" 时答案应该全是 $0$,因为 LCP 的值为 $-,1,2,3$,并没有出现过 $0$。

#include <cstdio>
#include <vector>
#include <cstring>
#include <algorithm>
using namespace std;

const int N=1000005;
const int C=26; //注意检查字符集大小!

//在结构题外开任何与SAM状态相关的数组,都需要N<<1
struct SuffixAutomaton
{
    int sz,lst; //状态数上限=2|S|-1 
    
    int len[N<<1],link[N<<1];
    int nxt[N<<1][C];
    //map<char,int> nxt[N<<1];
    //extend(char),并使用nxt[clone]=nxt[q]替换memcpy
    
    SuffixAutomaton()
    {
        len[0]=0,link[0]=-1;
        lst=sz=0;
    }
    
    void init()
    {
        for(int i=0;i<=sz;i++)
        {
            len[i]=link[i]=0;
            for(int j=0;j<C;j++)
                nxt[i][j]=0;
        }
        link[0]=-1;
        lst=sz=0;
    }
    
    void extend(int c)
    {
        int cur=++sz;
        len[cur]=len[lst]+1;
        
        int p=lst;
        while(p!=-1 && !nxt[p][c])
            nxt[p][c]=cur,p=link[p];
        
        if(p==-1)
            link[cur]=0;
        else
        {
            int q=nxt[p][c];
            if(len[p]+1==len[q])
                link[cur]=q;
            else
            {
                int clone=++sz;
                len[clone]=len[p]+1;
                link[clone]=link[q];
                memcpy(nxt[clone],nxt[q],C*4);
                
                while(p!=-1 && nxt[p][c]==q)
                    nxt[p][c]=clone,p=link[p];
                link[q]=link[cur]=clone;
            }
        }
        lst=cur;
    }
    
    int id[N<<1];
    
    //构建完成后,id顺序为len递增(逆拓扑序)【仅可排一次】 
    void sort()
    {
        static int bucket[N<<1];
        
        memset(bucket,0,sizeof(bucket));
        for(int i=1;i<=sz;i++)
            bucket[len[i]]++;
        for(int i=1;i<=sz;i++)
            bucket[i]+=bucket[i-1];
        for(int i=1;i<=sz;i++)
            id[bucket[len[i]]--]=i;
    }
}sam;

int n;
char s[N];

int main()
{
    int T;
    scanf("%d",&T);
    while(T--)
    {
        scanf("%s",s+1);
        n=strlen(s+1);
        
        sam.init();
        
        printf("0");
        sam.extend(s[1]-'a');
        
        int mx=0,zero=0;
        for(int i=2;i<=n;i++)
        {
            sam.extend(s[i]-'a');
            
            int val=sam.len[sam.link[sam.lst]];
            mx=max(mx,val);
            zero|=(val==0);
            printf(" %d",zero?mx+1:0);
        }
        printf(" \n");
    }
    return 0;
}
View Code

 

牛客38727C Cmostp(2022 牛客暑期多校 加赛)【LCT access 维护 parent 树的链并,离线转化】

这题其实只有模型抽象用到了 SAM,所以题解写在 LCT 那里了。右端点在区间 $[l,r]$ 的状态等价于字符串 $s[1:l]$ 到 $s[1:r]$ 这 $r-l+1$ 个状态在 parent 树上的链并。

 

(待续)

 

SAM 在处理字典序相关问题时比较乏力,在这种情况下从后缀树、后缀数组的角度往往具有更好的性质。

当问题中的询问数是 $|s|$ 级别的时,可以往 SAM 想一想。

posted @ 2022-09-13 01:05  LiuRunky  阅读(368)  评论(0编辑  收藏  举报