后缀算法杂烩
前置知识
后缀排序
不考虑暴力了,直接搞上正解。
我们设 \(sa[i],rk[i]\) 分别表示第 \(i\) 名的子串初始点在哪,以及以 \(i\) 开头的子串的排名。
我们考虑倍增的做法。先将长度为 \(1\) 的子串排序求出。
然后每次倍增长度,设长度为 \(w\),然后我们对于每个长度为 \(w\) 的子串以当前位置 \(x\) 的长度为 \(w/2\) 的子串的排名为第一关键字,以 \(x+w/2\) 位置的长度为 \(w/2\) 子串的排名为第二关键字,然后搞一波基数排序即可。
最后 \(w>n\) 时当场终止。
我们发现只是中途在长度刚好大于等于 \(w\) 时终止,也可以求出长度为 \(w\) 的子串的排序。
不难证明,这个复杂度是 \(O(n\log n)\) 的。
给出一个参考代码。
#include<bits/stdc++.h>
#define ll long long
#define db double
#define filein(a) freopen(#a".in","r",stdin)
#define fileot(a) freopen(#a".out","w",stdout)
#define sky fflush(stdout);
#define gc getchar
#define pc putchar
namespace IO{
inline bool blank(char c){
return c==' ' or c=='\n' or c=='\t' or c=='\r' or c==EOF;
}
inline void gs(char *s){
char ch=gc();
while(blank(ch) ) {ch=gc();}
while(!blank(ch) ) {*s++=ch;ch=gc();}
*s=0;
}
inline void gs(std::string &s){
char ch=gc();s+='#';
while(blank(ch) ) {ch=gc();}
while(!blank(ch) ) {s+=ch;ch=gc();}
}
inline void ps(char *s){
while(*s!=0) pc(*s++);
}
inline void ps(const std::string &s){
for(auto it:s)
if(it!='#') pc(it);
}
template<class T>
inline void read(T &s){
s=0;char ch=gc();bool f=0;
while(ch<'0'||'9'<ch) {if(ch=='-') f=1;ch=gc();}
while('0'<=ch&&ch<='9') {s=s*10+(ch^48);ch=gc();}
if(ch=='.'){
db p=0.1;ch=gc();
while('0'<=ch&&ch<='9') {s=s+p*(ch^48);ch=gc();}
}
s=f?-s:s;
}
template<class T,class ...A>
inline void read(T &s,A &...a){
read(s);read(a...);
}
};
using IO::read;
using IO::gs;
using IO::ps;
const int S=1e6+3;
int sa[S],rk[S<<1],lark[S];
int main(){
//filein(a);fileot(a);
char *c=new char[S];
gs(c+1);
int n=strlen(c+1);
int *id=new int[n+3];
int m=std::max(n,300);
int *cnt=new int[S+3];
for(int i=0;i<=m;++i) cnt[i]=0;
for(int i=1;i<=n;++i) ++cnt[rk[i]=c[i] ];
for(int i=1;i<=m;++i) cnt[i]+=cnt[i-1];
for(int i=n;i>=1;--i) sa[cnt[rk[i] ]--]=i;
delete []cnt; cnt=NULL;
for(int w=1;w<n;w<<=1){
cnt=new int[m+3];
for(int i=0;i<=m;++i) cnt[i]=0;
for(int i=1;i<=n;++i) id[i]=sa[i];
for(int i=1;i<=n;++i) ++cnt[rk[id[i]+w] ];
for(int i=1;i<=m;++i) cnt[i]+=cnt[i-1];
for(int i=n;i>=1;--i) sa[cnt[rk[id[i]+w] ]--]=id[i];
/*----cutline----*/
for(int i=0;i<=m;++i) cnt[i]=0;
for(int i=1;i<=n;++i) id[i]=sa[i];
for(int i=1;i<=n;++i) ++cnt[rk[id[i] ] ];
for(int i=1;i<=m;++i) cnt[i]+=cnt[i-1];
for(int i=n;i>=1;--i) sa[cnt[rk[id[i] ] ]--]=id[i];
/*----cutline----*/
for(int i=1;i<=n;++i) lark[i]=rk[i];
int p=0;
for(int i=1;i<=n;++i){
if(
i!=1 and
lark[sa[i] ]==lark[sa[i-1] ] and
lark[sa[i]+w]==lark[sa[i-1]+w]
){
rk[sa[i] ]=p;
}else{
rk[sa[i] ]=++p;
}
}
m=p;
delete []cnt; cnt=NULL;
}
for(int i=1;i<=n;++i) printf("%d ",sa[i]);
pc('\n');
delete []c; c=NULL;
delete []cnt; cnt=NULL;
return 0;
}
我们交一下上面的代码,发现常数巨大,所以我们考虑卡一下常数。
1.对于第二关键字的排序不需要基数排序
我们上一回求的 \(sa\) 可以直接用。发现后面的位置不存在长度为 \(w\) 的子串,直接倒着打入前几名。
然后其他的直接枚举 \(sa[i]\), 这个是枚举的排名为 \(i\) 的第二关键字的位置,然后 \(sa[i]-w\) 就是当前位置,直接接着打入排名即可。
这个是常数大的次要原因。
2.动态改变基数排序关键字枚举的值域
每次给 \(rk\) 上值的时候记录一下有几个不同的,我们需要值域就为那么大。
有稍微的常数优化。
3.使其尽量访问连续内存
我们可以用数组存储一下,防止反复算相同的东西,具体可以看看代码。
4.达到目的直接终止,不做多余操作
如果每个子串排名互不相同就可以退了。
优化后的代码如下。
#include<bits/stdc++.h>
#define ll long long
#define db double
#define filein(a) freopen(#a".in","r",stdin)
#define fileot(a) freopen(#a".out","w",stdout)
#define sky fflush(stdout);
#define Better_IO 1
namespace IO{
inline bool blank(const char &c){
return c==' ' or c=='\n' or c=='\t' or c=='\r' or c==EOF;
}
#if Better_IO==true
char buf[(1<<20)+3],*p1(buf),*p2(buf);
char buf2[(1<<20)+3],*p3(buf2);
const int lim=1<<20;
inline char gc(){
if(p1==p2) p2=(p1=buf)+fread(buf,1,lim,stdin);
return p1==p2?EOF:*p1++;
}
#define pc putchar
#else
#define gc getchar
#define pc putchar
#endif
inline void gs(char *s){
char ch=gc();
while(blank(ch) ) {ch=gc();}
while(!blank(ch) ) {*s++=ch;ch=gc();}
*s=0;
}
inline void gs(std::string &s){
char ch=gc();s+='#';
while(blank(ch) ) {ch=gc();}
while(!blank(ch) ) {s+=ch;ch=gc();}
}
inline void ps(char *s){
while(*s!=0) pc(*s++);
}
inline void ps(const std::string &s){
for(auto it:s)
if(it!='#') pc(it);
}
template<class T>
inline void read(T &s){
s=0;char ch=gc();bool f=0;
while(ch<'0'||'9'<ch) {if(ch=='-') f=1;ch=gc();}
while('0'<=ch&&ch<='9') {s=s*10+(ch^48);ch=gc();}
if(ch=='.'){
db p=0.1;ch=gc();
while('0'<=ch&&ch<='9') {s=s+p*(ch^48);ch=gc();}
}
s=f?-s:s;
}
template<class T,class ...A>
inline void read(T &s,A &...a){
read(s);read(a...);
}
};
using IO::read;
using IO::gs;
using IO::ps;
const int S=1e6+3;
int sa[S],rk[S<<1],kc[S],lark[S];
inline bool cmps(int a,int b,int c){
return
b!=0 and
lark[a]==lark[b] and
lark[a+c]==lark[b+c];
}
int main(){
//filein(a);fileot(a);
char *c=new char[S];
gs(c+1);
int n=strlen(c+1);
int *id=new int[n+3];
int m=std::max(n,128);
int *cnt=new int[S+3];
for(int i=0;i<=m;++i) cnt[i]=0;
for(int i=1;i<=n;++i) ++cnt[rk[i]=c[i] ];
for(int i=1;i<=m;++i) cnt[i]+=cnt[i-1];
for(int i=n;i>=1;--i) sa[cnt[rk[i] ]--]=i;
delete []cnt; cnt=NULL;
for(int w=1;w<n;w<<=1){
cnt=new int[m+3];
int p=0;
for(int i=n;i>n-w;--i) id[++p]=i;
for(int i=1;i<=n;++i){
if(sa[i]>w) id[++p]=sa[i]-w;
}
/*----cutline----*/
for(int i=0;i<=m;++i) cnt[i]=0;
for(int i=1;i<=n;++i) ++cnt[kc[i]=rk[id[i] ] ];
for(int i=1;i<=m;++i) cnt[i]+=cnt[i-1];
for(int i=n;i>=1;--i) sa[cnt[kc[i] ]--]=id[i];
/*----cutline----*/
for(int i=1;i<=n;++i) lark[i]=rk[i];
m=0;
for(int i=1;i<=n;++i){
rk[sa[i] ]=cmps(sa[i],sa[i-1],w)?m:++m;
}
delete []cnt; cnt=NULL;
if(m==n) break;
}
for(int i=1;i<=n;++i) printf("%d ",sa[i]);
pc('\n');
delete []c; c=NULL;
delete []cnt; cnt=NULL;
return 0;
}
码风更新一下:
#include<bits/stdc++.h>
#define ll long long
#define db double
#define file(a) freopen(#a".in","r",stdin),freopen(#a".out","w",stdout)
#define sky fflush(stdout)
#define gc getchar
#define pc putchar
namespace IO{
template<class T>
inline void read(T &s){
s=0;char ch=gc();bool f=0;
while(ch<'0'||'9'<ch) {if(ch=='-') f=1;ch=gc();}
while('0'<=ch&&ch<='9') {s=s*10+(ch^48);ch=gc();}
if(ch=='.'){
T p=0.1;ch=gc();
while('0'<=ch&&ch<='9') {s=s+p*(ch^48);p/=10;ch=gc();}
}
s=f?-s:s;
}
template<class T,class ...A>
inline void read(T &s,A &...a){
read(s);read(a...);
}
};
using IO::read;
const int N=1e6+3;
char c[N];
int n,m;
int rk[N],id[N],cnt[N],sa[N],lark[N];
inline bool cmps(int x,int y,int c){
return (
y!=0 and
lark[x]==lark[y] and
lark[x+c]==lark[y+c]
);
}
int main(){
file(a);
scanf("%s",(c+1) );
n=strlen(c+1);
m=std::max(128,n);
for(int i=0;i<=m;++i) cnt[i]=0;
for(int i=1;i<=n;++i) ++cnt[rk[i]=c[i] ];
for(int i=1;i<=m;++i) cnt[i]+=cnt[i-1];
for(int i=n;i>=1;--i) sa[cnt[rk[i] ]--]=i;
for(int w=1;w<n;w<<=1){
int p=0;
for(int i=n;i>n-w;--i) id[++p]=i;
for(int i=1;i<=n;++i)
if(sa[i]>w) id[++p]=sa[i]-w;
for(int i=0;i<=m;++i) cnt[i]=0;
for(int i=1;i<=n;++i) ++cnt[rk[id[i] ] ];
for(int i=1;i<=m;++i) cnt[i]+=cnt[i-1];
for(int i=n;i>=1;--i) sa[cnt[rk[id[i] ] ]--]=id[i];
for(int i=1;i<=n;++i) lark[i]=rk[i];
m=0;
for(int i=1;i<=n;++i){
rk[sa[i] ]=cmps(sa[i],sa[i-1],w)?m:++m;
}
if(m==n) break;
}
for(int i=1;i<=n;++i){
printf("%d ",sa[i]);
}pc('\n');
return 0;
}
这下清爽多啦。
H数组
定义 \(high[i]\) 为 \(sa[i]\) 和 \(sa[i-1]\) 的 LCP
的长度。
\(h[i]\) 表示 \(i\) 开头的后缀与排名在其前一位的后缀的 LCP
于是有 \(h[i]=high[rk[i] ]\) 和 \(high[i]=h[sa[i] ]\)
性质1
.\(high[rk[i]]\ge high[rk[i-1] ]-1\)
通过这个可以直接暴力求解 \(high\) 数组。
inline void getHigh(char *s){
int t=0;
for(int i=1;i<=n;++i){
if(t) --t;
while(s[i+t]==s[sa[rk[i]-1]+t])
++t;
high[rk[i] ]=t;
}
}
应用
两后缀LCP
区间 \(high\) 数组取 min
即可。
求子串出现次数
对于子串 \(s[l,r]\) ,求出前缀为 \(s\) 的字典序最大最小的后缀的排名 \(x,y\),答案 \(y-x+1\)。
等价求有多少个后缀与 \(s[l,n]\) 的 \(LCP\ge r-l+1\)。
比较子串大小
对于 \(A=[a,b],B=[c,d]\),若 \(LCP(a,c)\ge \min\{|A|,|B|\}\) ,那么比较 \(|A|,|B|\) 大小,否则比较 \(rk[a],rk[c]\)。
本质不同子串个数
\(\frac{n(n+1)}{2}-\sum_{i=2}^nhigh[i]\)
计算所有新增的 LCP
即可。
最小表示法
求一个串循环同构中字典序最小的那个。
倍长原串,其同构转为后缀后字典序相对大小不改变,于是可得。
为什么相对大小不变?因为发现多出来的部分是个前缀,你原来比它字典序小,那之后仍然是小的。特殊情况是两个相同的循环同构,但是由于两者相同求出来的结果不变。
最长公共子串
等价于两(或多个)个串拼起来求最大的后缀 LCP
。两个串间添加不同的不在字符集中的元素作为分隔符。
我们考虑二分答案后,检验是否可达。具体就是找到一段连续段 \([l,r]\) 满足 \(x\in[l,r],high[x]\ge ans\) ,且对于 \(x\in[l-1,r]\) ,使得每个串都有至少一个元素在里面。
后缀自动机(SAM)
等价类:这类中最长子串 \(s\) 的后缀在原串中出现位置相同。
一个字符串出现且仅出现在一个等价类。
包含类:子串 \(s\) 的后缀 \(s'\) 在原串的出现次数比 \(s\) 多,那么称 \(s'\) 是 \(s\) 的包含类。
SAM
的节点的父节点总是其包含类。
SAM
的 parent树
,parent指针
指向其出现次数多于它的最长后缀(字符串长度最长的包含类)。
性质1
一个等价类里的串长度连续。
因此可用二元组 \((s,l)\) 表示一个等价类。(最长字符串和最短字符串长度)
为什么呢?你可以发现当两个串出现位置(结尾出现位置,之后称之为 endpos
)相同,当且仅当其中一个是另一个的后缀且只以其后缀的形式出现。
性质2
两个子串的 endpos
要么没有交集,要么满足包含关系(当其中一个是另一个的后缀)。
性质3
\(len_{min}(u)=len_{max}(prt(u) )+1\)
这样的话我们就只需要在节点储存最长的字符长度即可了,最短可由父亲推出。
如果我们定义初始状态为 \(t0\) , 其为空字符并规定其 endpos
为 \(\{0,\dots,|S|\}\)。然后可以发现 parent指针
构成一棵以 \(t0\) 为根的树。(虽然方向是反的)
然后我们可以把加字符的操作看作分裂一个点。
可以证明,点数和边数都是 \(O(n)\) 的。(具体来说点数最多 \(2n-1\),边数至多 \(3n-4\))
SAM的构造
我们对于每个点保存其转移、最长字符串长度及其 parent指针
。
考虑增量法构造,也就是把字符一个一个加进来动态维护 SAM
的形态。
我们记当前新建的状态为 \(cu\) ,上一个状态是 \(la\) 。首先 \(cu\) 是 \(la\) 后面加了一个字符,可以更新 \(len_{max}(cu)=len_{max}(la)+1\)。
然后我们考虑如何更新 prt指针
。
对于新加的字符 \(c\), \(cu=\{S+c,|S|+1\},p=\{S,|S|\}\)。
容易发现,\(p\) 及其 \(prt\) 树祖先必定出现位置全覆盖 \(cu\) 的出现位置。(毕竟 \(cu\) 是 \(p\) 加了一个字符 \(c\))
那么如果有的话,我们在 \(p\) 及其祖先中查找一个尽可能长的 \(q=\delta(p',c)\not=\emptyset\) 是一定可以找到 \(prt(cu)\) 的。
首先 \(q\) 一定是 \(cu\) 的包含类。首先 \(q\) 删掉 \(c\) 一定是 \(p\) 的后缀,然后 \(q\) 就是 \(cu\) 的后缀了。然后就肯定是包含类了。
我们插入 \(c\) 后,\(q\) 中属于 \(cu\) 的部分就会出现次数 \(+1\) 。(剩下的是包含 \(cu\) 的串)
相当于要把这个节点给拆开成两部分。
分类讨论:
1
. \(p=t0\)
直接没有了,\(prt(cu)=t0\)
2
. \(p'+c=q\)
整个直接全部属于 \(cu\) , \(prt(cu)=q\)
3
. \(p'+c\subseteq q\)
那么接下来,就要分裂了。
考虑把 \(q=\{s,l\}\) 拆成 \(q1=\{s,|p'|+2\},q2=\{s[|s|-|p'|,|s|],l\}\)。(即 \(q1\) 是包含 \(cu\) 的部分,\(q2\) 是被 \(cu\) 包含的部分,这里说的包含是指字符串上的)
然后令 \(prt(q1)=prt(cu)=q2\) 即可。
再考虑更新 \(\delta\) 转移。
由 性质2
可得对于向新字符 \(c\) 的转移,只有其 \(prt\) 树上的祖先可以通过转移到达。
我们需要更新的只有 \(cu\) 和 \(q\) 。首先是 \(p\) 到 \(p'\) 路径上的需要更新(在 SAM
中转移是转向添加字符后最短的以其为后缀的子串,但本质是在那个串在那里面),然后是分裂了的话,上面的转移到 \(q\) 的需要接到 \(q2\) 。(因为这部分出现次数一起 \(+1\) 了)
实现的时候 \(q2\) 的 \(len_{max}\) 赋值为 \(len_{max}(p')+1\) ,因为其最大的串都可以转移到 \(q2\)。当然从上面的式子 \(q2=\{s[|s|-|p'|,|s|],l\}\) 也不难看出。
代码实现
struct Suf_Auto_Machine{
Suf_Auto_Machine(){
tot=1;
}
#define fa(x) t[x].fa
int tot;
struct node{
int len,fa,ch[S];
}t[N<<1];
inline void ins(int x){
static int la=1;
int cu=++tot,p=la;
t[cu].len=t[la].len+1;la=cu;
while(p and !t[p].ch[x]){
t[p].ch[x]=cu;
p=fa(p);
}
if(!p) fa(cu)=1;
else{
int q1=t[p].ch[x];
if(t[q1].len==t[p].len+1) fa(cu)=q1;
else{
int q2=++tot;
t[q2].len=t[p].len+1;fa(q2)=fa(q1);
memcpy(t[q2].ch,t[q1].ch,sizeof(t[q1].ch) );
fa(q1)=fa(cu)=q2;
while(p and t[p].ch[x]==q1){
t[p].ch[x]=q2;
p=fa(p);
}
}
}
}
#undef fa
}s;
广义SAM
支持多模式串问题!性质什么的和普通 \(SAM\) 基本没区别。本质是把很多个串楞压缩在一起。(但是要求 \(endpos\) 的话 \(sz\) 要分开记录)(没听懂我在说什么?可以先看下下面的应用再来看这个)
优势是不管其他串是否有这个状态,总之是统一状态的游走。
怎么构造?
偷懒做法
直接一个一个串插入, 每个串插入前 \(la=1\) 。好写,但是复杂度在串很多的情况下不一定保证。
正确做法
在线:
直接一个一个串插入, 每个串插入前 \(la=1\) ,但是要加上特判。要知道,之前加入是没有 \(ch[la][x]\) 非空的情况的,我们只需特判这个即可。
如果,\(len(ch[la][x])=len(la)+1\),那么直接就是(同一个东西)啊。
否则就直接拆开。
挺好啊。
离线:
建 \(trie\) 然后按 \(BFS\) 顺序加入。使得 \(len\) 单调不降。
很好啊,又慢又难写,爪巴!
这里只实现在线版:
struct EXSuf_Auto_Machine{
#define fa(x) t[x].fa
#define len(x) t[x].len
EXSuf_Auto_Machine(){
tot=1;
}
int tot;
struct node{
int len,fa,ch[26];
}t[N<<1];
inline void ins(int x,int &la){
int p=la;
if(t[p].ch[x]){
int q1=t[p].ch[x];
if(len(q1)==len(p)+1) {la=q1;return;}
int q2=++tot;
fa(q2)=fa(q1);len(q2)=len(p)+1;
memcpy(t[q2].ch,t[q1].ch,sizeof(t[q1].ch) );
fa(q1)=q2;
while(p and t[p].ch[x]==q1){
t[p].ch[x]=q2;
p=fa(p);
}
la=q2;return;
}
int cu=++tot;
len(cu)=len(la)+1;
while(p and !t[p].ch[x]){
t[p].ch[x]=cu;
p=fa(p);
}
if(!p) fa(cu)=1;
else{
int q1=t[p].ch[x];
if(len(q1)==len(p)+1) fa(cu)=q1;
else{
int q2=++tot;
fa(q2)=fa(q1);len(q2)=len(p)+1;
memcpy(t[q2].ch,t[q1].ch,sizeof(t[q1].ch) );
fa(q1)=fa(cu)=q2;
while(p and t[p].ch[x]==q1){
t[p].ch[x]=q2;
p=fa(p);
}
}
}
la=cu;return;
}
#undef fa
#undef len
}s;
应用
等价类出现次数
计算 \(prt\) 树上的子树(算上当前这个点)中的前缀的个数即可。需要注意的是分裂出来的 \(q2\) 是不算的,不然会算重。因为它一定是其子树中串的后缀,而且我们知道枚举所以前缀看后缀是否匹配就是出现次数。但是一个等价类被拆开,其中只有 \(q1\) 含有那个作为前缀的部分,于是分裂出来的 \(q2\) 不算。
本质不同子串个数
由于每个不同子串在 SAM
上对应一条从 \(t0\) 出发的路径,于是 DP
计数即可。
另一种解法是直接计算:
算这个的话甚至可以动态计算。动态计算的话不要算 \(q2\) ,它之前已经算过了。
第 \(k\) 大子串
即找到第 \(k\) 大路径。计算每个节点出发的路径数后可以方便处理。
最小表示法
倍长字符串,找到最小的长度为 \(|S|\) 的路径。直接从 \(t0\) 一直往小的走即可。
出现位置
输出 \(endpos\) 集合,只需要知道子树中的前缀的位置即可。位置是哪里?\(len\)。
最长公共子串
先考虑两个串 \(S\) 和 \(T\) ,用 \(T\) 的每个前缀去匹配 \(S\) 的后缀,求最大的匹配即可。
具体是 \(T\) 一位一位加入进来,在 \(S\) 的 \(SAM\) 上匹配,具体来说就是有这位的转移就直接跳,否则向上跳 \(prt\) 直到找到一个有这位转移的。
LCS - Longest Common Substring
再考虑多个串,我们建立广义后缀自动机。
然后转移看是否每个串都在这个点的子树有前缀即可。记得特判里面也要打 \(sz\) 标记。
LCS2 - Longest Common Substring II
后缀树
所有后缀建一个 trie
叫 后缀树
。
度数不为 \(2\) 的点和结尾节点(接受状态点)为关键点建虚树,这个叫 简化后缀树
。
你发现后缀树是把一个关键点到其上方的关键点的点省略掉,且这一堆点的在原串的出现次数相同,这些省去的部分是其前缀,于是可以发现如下结论。
性质1
.SAM
的 parent
树是其反串的简化后缀树。
那么要建后缀树只要对其反串建 SAM
即可。
等价类:后缀树上一个关键点及其到上方关键点间的点(不包括上方关键点)。
包含类:后缀树上的祖先。
那么其实分裂就是加字符进来的时候增添了一个关键点 \(LCA\),于是产生分裂的效果。
在后缀树上,只有开始位置相同的两个位置会被压缩。
大概就这样吧,之后应用自己练题便是。