后缀数组学习

缘起

终于要重新复习SA了. 开心，撒花~ 这次一定要把SA学的清清楚楚明明白白的~

分析

按照罗穗骞大神的09年国家队论文 《后缀 数 组—-处理字符串的有力工具》 的说法. 

1
2
3
4
缀数组是处理字符串的有力工具。 后缀数组是后缀树的一个非常精巧的替代品， 
它比后缀树容易编程实现， 能够实现后缀树的很多功能而时间复杂度也并不逊色，
而且它比后缀树所占用的内存空间小很多。
可以说，在信息学竞赛中后缀数组比后缀树要更为实用。

其实本蒟蒻并不认为后缀树编码复杂（敲顺了也就那样），但是后缀树吃内存是真的. 所以SA还是有必要学习的.

什么是后缀数组

就是后缀的数组. 这好像没说~ 囧 ~ 

1
2
一个字符串的后缀数组指的是将它所有的后缀（包括空串和原串）按字典序升序排成的数组就叫做后缀数组.
其中认为空串比任何字符串的字典序都小.

例如 S=”abracadabra” 的后缀数组如下

即排在sa第i位的是sa[i], 它对应的后缀是S[sa[i],…], 即sa[i] 回答的是字典序排第i小的后缀是原串中哪个后缀.

后缀数组的实现

后缀数组正如上面罗菊苣介绍的那样，可以在一大票的复杂字符串问题中大显神威, 所以能快速计算出后缀数组很重要. 这就是SA的实现.

假设|S|=n, 则最朴素的算法是将n个后缀排序, 而将n个长度O(n)的字符串排序的复杂度是O(n^2*logn).

但是这种复杂度n一旦达到了1e4的规模就是不可接受的. 所以我们要好好利用他们都是S的后缀这一特点寻求更加高效的算法. 一般的，sa的实现有 倍增和dc3两种算法.

倍增

倍增法是由Manber&Myers发明的. 倍增法的思想是首先计算从每个位置开始长度为2的子串的字典序升序排序, 再利用这个结果计算长度为4的子串的字典序升序排序, 接下去再用4的结果计算长度为8的子串的排序. 如此这般, 不断倍增. 直至长度>=n就得到了后缀数组. 一图胜千言, 还是以上面那个字符串”abracadabra”为例来看看倍增算法的过程（下图顶部的单个a是字符串S最后那个a, 因为它后面没有字符了. 所以不论如何倍增, 一直是孤零零的一个’a’在那里）

下面我们用S[i,k]表示从S[i]开始的长度为k的子串. 其中, 当从位置i开始剩余字符数量不足k时, 表示的是从位置i到字符串末尾的子串。

要计算长度为2的子串的顺序. 只需要排序2个字符组成的数对就好了. 现在假设已经求得了长度为k的子串的顺序, 要求长度为2k的子串的顺序. 于是不得不引入一个数组 rank, 其中 rank[i] 表示S[i,k] 排名第几. 即rank和sa是互逆的. 要计算长度为2k的子串的顺序, 就只需要对2个rank组成的数对进行排序就可以了. 即我们通过对 (rank[i], rank[i+k]) 和 (rank[j], rank[j+k])这两个数对进行字典序比较即可. 这就实现了高效的比较字符串大小. 这就是倍增的核心思想. 同样是一图胜千言

注意, 这里rank的排名采用并列制. 例如上图rank一列出现了好多2（即并列第二名）.

下面上倍增（double augment da）的板子

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
#include "stdafx.h"
#include <stdio.h>
#include <string.h>
#include <algorithm>
using namespace std;
#define LOCAL

const int maxn = 1e5+5; // 字符串的最大长度
int n,k, r[maxn], t[maxn], sa[maxn]; // n是字符串的长度, k是每轮比较的后缀的长度（k将倍增）,r是原字符串的rank数组. t用于暂存中间结果, sa用于保存s的后缀数组
char s[maxn]; // 字符串

bool c0(int i, int j) // 比较(r[i], r[i+k]) 和 (r[j], r[j+k])
{
	if (r[i]!=r[j])
	{
		return r[i]<r[j];
	}
	int ri = i+k<=n?r[i+k]:-1, rj = j+k<=n?r[j+k]:-1;
	return ri<rj;
}

void da() // double augment 倍增算法求后缀数组
{
	k = 1; // 初始长度为1
	n = 0;
	while(s[n])
	{
		sa[n] = n; // 最初sa是没有任何排序的
		r[n] = s[n]-'a'+1;
		n++;
	} // 最后s[n]=0
	r[n]=0,sa[n] = n; // 空串
	while(k<=n)
	{
		sort(sa, sa+n+1, c0); // 最初从长度为1的后缀去排序长度为2的后缀
		t[sa[0]] = 0; // 先在t中暂存新计算出的r, 最后再存回r中, t(可以理解为r)是sa的互逆
		for (int i = 1;i<=n;i++) // sa[i]的rank肯定没有sa[i-1]的高啊
		{
			t[sa[i]] = t[sa[i-1]]+c0(sa[i-1],sa[i]); // 因为这里c0要用原先的r,所以只能暂存新的r到t中,注意, 因为c0可能返回1或者0,所以新的t两个相邻的可能一样(即并列排名)
		}
		memcpy(r,t,sizeof(int)*(n+1)); // 存回r中
		k<<=1;
	}
}

int main() {
#ifdef LOCAL
	freopen("d:\\data.in", "r", stdin);
#endif
	while(~scanf("%s",s))
	{
		da(); // 倍增算法求s的后缀数组, 结果存在 sa[1,...,n]中, sa的取值在[0,...,n-1]范围内,和下面的dc3是一样的
	}
	return 0;
}

测试数据

1
2
3
aabaaaab
aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa
abracadabra

dc3（Difference Cover mod 3）

dc3算法是细节党. 思想并不难懂, 但是细节比较多. 建议初学者（譬如在下）学习的时候理解其思想, 而并不要纠结于细节. 不然容易自闭.

就像【1】中说的那样, dc3算法分成三步

  1. 将所有的后缀分成两部分，第一部分是模3不等于0的，比如Suffix[1],Suffix[2],Suffix[4],Suffix[7]等等，第二部分是模3等于0的后缀，Suffix[0],Suffix[3]等等。这一步计算的是第一部分的后缀排名（就当第二部分后缀不存在）. 计算的方法是3轮基数排序（, 之所以是3轮, 是因为是dc3算法, 是模掉了3, 而每次对关键字的排序使用计数排序, 即代码中的sort就是计数排序, 之所以不是dc4、dc5、dc6, 是因为那基数排序的轮数就太多了,为什么不是dc2算法? 因为如果是dc2的话, 即计算第一部分后缀的排序的话, 依据的就是前2个字符, 这样不容易区分出来第一部分后缀, 所以经常陷入递归而实际跑的时候性能并不如dc3算法.）

    例如S=”aabaaaaba”, 它的第一部分后缀是S[1,…], S[2,…], S[4,…],S[5,…], S[7,…,],S[8,…] 共计6个后缀, 则我们首先根据这6个后缀的前三个字符作为三个关键字进行3轮计数排序（也就是根据这三个关键字进行基数排序）.首先, 这6个后缀一定不会出现两个字典序相等的状况（毕竟长度都不一样），但是仅仅根据前三个字符进行字典序排序之后可能并不能区分某些个后缀. 碰到这种情况, 就要将继续比较. 具体做法是, 将这6个后缀中的每个后缀的前三个字符”压缩”成1个. 则这6个后缀变成了长度为6的新的字符串（注意, 长度由9变成了6，即2/3衰减），要再将这新的字符串的后缀排序，也就是要对这个新的字符串进行递归的处理. 但是注意, 该新的字符串的长度缩小了(缩小比例恒定为2/3). 所以递归必定会在有限步结束.

    而且据此, 我们就能看出dc3的复杂度是 

    1
    f(n)=O(n)+f(2/3*n)

    所以dc3是O(n)的优秀算法. 这就是为什么dc3的复杂度优于da. 因为da没有用递归. 它仅仅是二分. 所以da的复杂度是O(nlogn).

  2. 利用上一步的结果, 对第二部分后缀进行排序. 为什么能用上一步的结果? （注意，其实dc3和da很类似，都是要能利用前面的结果, 这样算法才能高效）。好了，不卖关子了，因为第二部分的后缀，可以看做一个字符A加上某个第一部分的一个后缀S（确切讲, 是模掉3余1的位置起始的后缀），而将A和S分别视作是基数排序的两个关键字. 但是S这个关键字在第一步是已经排好序了的，不需要再排序了，所以这一步只需要一轮计数排序（即代码中的sort）就可以搞定.

  3. 既然第一部分和第二部分都已经有序了，那么要整体有序，只要学过归并排序的童鞋都知道这最后一步是干什么的———自然是归并第一部分的有序后缀和第二部分的有序后缀咯~ 当然, 整体思想是这样, 但是具体细节依然比较考究——不卖关子了——就是怎么能高效处理归并过程中对来自2部分的后缀之间的大小比较? 如果是遍历的话, 算法的性能还是不行的。事实上，我们可以这样

    分两种情况。第一种是比较Suffix[3i] 和Suffix[3j+1]，我们把它们看做：

    Suffix[3i]=S[3i]+Suffix[3i+1]

    Suffix[3j+1]=S[3j+1]+Suffix[3j+2]

    Suffix[3i+1]和Suffix[3j+2]可以直接比较，因为它们都属于第一部分，而S[3i]和S[3j+1]仅仅是两个字符之间的比较自然O(1)时间；

    第二种情况是Suffix[3i] 和Suffix[3j+2]，把它们看做是

    Suffix[3i]=S[3i]+Suffix[3i+1]

    Suffix[3j+2]=S[3j+2]+Suffix[3(j+1)]

    而S[3i]和S[3j+2]都是字符, 可以直接比较, 从而将问题转换为了Suffix[3i+1]和 Suffix[3(j+1)] 两个来自第一部分和第二部分的后缀的比较. 显然, 因为参与比较的后缀的长度变短了, 所以问题规模缩小——即我们可以用递归来高效完成归并过程中的遇到了来自两部分的后缀之间的比较， 即代码中的c12函数.

好了，说了这么多, 上dc3的板子咯~

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
#include "stdafx.h"
#include <stdio.h>
#define LOCAL
#define F(x) ((x)/3+((x)%3==1?0:tb))
#define G(x) ((x)<tb?(x)*3+1:((x)-tb)*3+2)
// 宏F(x)的作用是计算原后缀S[x,...]在缩减为2/3的"新字符串"(在基数排序前)中的索引, 宏G的作用与F相反. 这个F宏有讲究, 啥讲究呢? tb是啥? 就是模掉3余1的后缀, 而F宏的目的就是将模掉3余1的后缀全部排在模掉3余2的后缀的前面,为什么要这样呢? 因为下面收集第二部分后缀的时候(即模掉3余0的后缀), 采用了if (san[i]<tb){wb[ta++] = san[i]*3;}, 而模掉3余0的后缀其实就是一个字符拼上模掉3余1的后缀, 而这个代码的作用就是对模掉3余0的后缀进行第二关键字的排序,而这个第二关键字不就在对第一部分后缀排序之后的<tb的位置上吗? 所以F宏是这样来设置的
const int maxn = 1e5+5; // 字符串最大长度
char s[maxn]; // 字符串
int n,sa[maxn*3], r[maxn*3], wa[maxn*3], wb[maxn*3], wv[maxn*3], ws[maxn*3]; // 字符串长度为n,sa用于存储最后的后缀数组, r用于存储原字符串-'a'+1的结果(例如s中的'a'变成r中的1),wa和wb用于基数排序倒来倒去(为啥要三倍空间? 因为可能发生递归,而因为每次是对新字符串进行递归, 而新字符串的长度缩减为2/3, 是等比数列, 所以开3maxn使得递归不需要再分配空间了), ws用于计数排序, wv用于第三步的归并排序.

void sort(int *r, int *a, int *b, int n, int m) // 对a进行计数排序,a的长度是n, 字符集大小是m, 计数排序的结果保存在b中, 但是a仅仅保存的是索引, 真正的值保存在r中
{
	for (int i = 0;i<n; i++)
	{
		wv[i] = r[a[i]];
	} // wv保存参与计数排序的数值, 因为a中保存的仅仅是r中的索引
	for (int i = 0;i<m;i++)
	{
		ws[i] = 0;
	} // 初始化计数排序
	for (int i =0; i<n; i++)
	{
		ws[wv[i]]++;
	}
	for (int i = 1;i<m;i++)
	{
		ws[i]+=ws[i-1];
	}
	for (int i = n-1; ~i; i--)
	{
		b[--ws[wv[i]]]=a[i];
	} // 计数排序
}

int c0(int *r, int a, int b) // 比较两个后缀的前三个字符
{
	return r[a]==r[b]&&r[a+1]==r[b+1]&&r[a+2]==r[b+2];
}

bool c12(int k, int *r, int a, int b) // k是b模掉3的余数, 这里就是比较r[a,...]和r[b,...]这两个后缀,就像文章中说的那样, 用归并
{
	if (k==2) // 如果是 模掉3余0的后缀和模掉3余2的后缀比较, 则转换为模3余0和模3余1的后缀的比较, 即递归
	{
		return r[a]<r[b]||r[a]==r[b]&&c12(1,r,a+1,b+1);
	}
	else // 如果是模掉3余0的后缀和模掉3余1的后缀比较, 则可以转换为模掉3余1和模掉3余2的后缀的比较
	{
		return r[a]<r[b]||r[a]==r[b] && wv[a+1]<wv[b+1]; // 因为wv中已经存了模掉3不为0的后缀的排名
	}
}

void dc3(int *r, int *sa, int n, int m)
{
	int *rn = r+n, *san = sa+n, ta = 0, tb = (n+1)/3, tbc=0, p; // rn用于存储文章中说的"新的字符串", san是新字符串的后缀数组, ta是S[0,...,原字符串长度,即n-1,因为这里考虑空串]中模掉3余0的后缀的个数, tb是模掉3余1的后缀的个数, tbc是模掉3余1或者2的后缀的个数
	r[n] = r[n+1] =0; // 原因在dc3调用的时候说了
	for (int i=0; i<n; i++) // 注意, 传入的n其实是原字符串长度+1(即考虑空串)
	{
		if (i%3!=0)
		{
			wa[tbc++] = i;
		}
	} // 收集模掉3余数不为0的后缀
	sort(r+2, wa, wb, tbc, m); // 根据模掉3不为0的后缀的第3个字符(第三关键字)计数排序,为什么是第三个字符? 因为你将r+2带入到sort函数中去不就是wv[i] = r[a[i]+2]吗?这不就是后缀的第三个字符吗?
	sort(r+1, wb,wa, tbc, m); // 根据模掉3不为0的后缀的第2个字符(第二关键字)计数排序
	sort(r, wa,wb,tbc, m);// 根据模掉3不为0的后缀的第1个字符(关键字)计数排序
	// 上述三轮计数排序构成一次基数排序
	rn[F(wb[0])]=0; p=1; // 开始存储新的字符串rn(里面的值是在新的字符串rn中的索引,从0开始, 其实就是新的字符串的rank数组), wb[0]是字典序最小的模掉3不为0的后缀, wb们是模掉3不为0的后缀在原字符串中的索引. 经过F宏的映射得到的是在新的字符串中的索引, 而因为wb已经排好序, 所以得到的rn就是rank数组
	for (int i = 1; i<tbc; i++)
	{
		rn[F(wb[i])] = c0(r, wb[i-1], wb[i])?p-1:p++; // 如果wb[i]和wb[i-1]相等的话, 即与前面的相同, 就是文章中说的 "仅仅根据前三个字符进行字典序排序之后可能并不能区分某些个后缀" ,即出现并列排名的情况
	} // 计算剩下的新的字符串rn中的索引, 即完成新的字符串的rank数组rn的计算, 注意, 因为c0仅仅比较了每个模掉3不为0的后缀的前三个字符, 所以可能出现并列排名, 这种情况, 我们需要递归
	if (p<tbc)
	{ // 注意, 下面一行递归的时候, rn做为新的r传入, 其实就是在r后面移动,所以r伊始要开三倍空间,同理sa也要开三倍空间
		dc3(rn, san, tbc, p); // 如果前三个字符并不能完全区别这些后缀的话, 即出现并列排名（亦即rn中有相同的值），则就要对rn（注意, 长度已经衰减为原字符串s的2/3, 所以不必担心会超出伊始开的3N空间, 1+2/3+(2/3)^2+...=3, 这就是为什么伊始要开三倍空间的缘故）
	} // 惊为天人的递归, 最后递归出口之后层层返回就会回到最初的rn=r+n, 届时, rn中已经是彼此不同的数字了（即模掉3不为0的后缀们已经完全被区别开来了）
	else // 如果模掉3不为0的后缀们仅仅通过前三个字符已经可以区分开来了, 则可以得到第一部分后缀（即模掉3不为0的后缀）的sa数组——san了
	{
		for (int i =0 ; i<tbc; i++)
		{
			san[rn[i]] = i; // 很好理解——sa和rank互逆
		}
	}
	/************************************************************************/
	/* 至此, 第一部分后缀（即模掉3余1）的后缀数组已经计算完毕               */
	/************************************************************************/
	for (int i = 0; i<tbc; i++)
	{
		if (san[i]<tb) // tb是模掉3余1的后缀的个数,每个模掉3余1的后缀前面加上一个字符就是一个模掉3余0的后缀呀, 所以这里收集tb个, 然后再乘以3得到的就是
		{
			wb[ta++] = san[i]*3;
		}
	} // 开始收集第二部分后缀, 注意, 这里的收集方式是按照第二关键字排好序的——即模掉3余0的后缀的第一个字符后面的模掉3余1的后缀(它作为关键字, 这样才能实现高效判断)
	if (n%3==1) // 则因为上面收集的是<n的模掉3不为0后缀（而没有S[n,...],即空串）, 而这里n-1的确是模掉3为0的后缀,但是没有对应的模掉3余1的后缀跟在后面,所以要补上,其实这里补上的S[n,...]就是空串,因为传入的n其实是字符串长度+1,补上的S[n-1,...]这个模掉3余0的后缀自然是空串啦
	{
		wb[ta++] = n-1; // 则 s[n-1,...]是最后一个模掉3余0的后缀(注意, 因为传入的n其实是数组长度+1,, 所以这里的n-1其实是s的长度, 则S[n-1,...]其实代表空串)
	} // 注意, 可能有读者会问了——作为第二关键字（即排在模掉3余0的后缀后面的那个模掉3余1的后缀）排序的话,将空串排在n-1的位置是不是有点靠后了? 它可是按字典序最前面的呀! 不要紧的, 因为没有任何后缀比空串的第一关键字0来的大(相等的都没有).所以尽管第二关键字的排序对于空串不公,但是在马上下一行进行的第一关键字排序它绝壁能用自己的第一关键字作为实力重新抢回第一名的位置
	sort(r,wb,wa,ta,m); // 将wb中收集的模掉3余0的后缀按照第一关键字排序,则就完成了第一部分后缀的排序,排序结果放在wa中.
	/************************************************************************/
	/* 至此, 第二部分（即模掉3余0）的后缀的后缀数组已经计算完毕             */
	/************************************************************************/
	for (int i =0;i<tbc; i++)
	{
		wv[wb[i]=G(san[i])] = i; // 记录wv的原因是下面的c12要用
	} // 第二部分后缀不是排在wb中了吗? 所以第三步是要归并第一部分后缀和第二部分后缀, 所以现在要将第一部分后缀放进wa中去,而第一部分的后缀的信息已经完全变成san了.而san中保存的是在新的字符串中的索引, 所以要用G宏还原回去
	int i=0,j=0;
	for (p=0; i<ta && j<tbc; p++) // 归并第一部分后缀 和 第二部分后缀
	{
		sa[p] = c12(wb[j]%3, r, wa[i], wb[j])?wa[i++]:wb[j++]; // 按照归并的套路, 自然要比较一番当前准备出列的第一部分后缀和第二部分后缀啦
	}
	while (i<ta)
	{
		sa[p++]=wa[i++];
	}
	while(j<tbc)
	{
		sa[p++]=wb[j++];
	} // 这是归并的老套路了,有木有
	/************************************************************************/
	/* 至此, 第三步归并也已经完成,至此,dc3算法结束                          */
	/************************************************************************/
}

int main() {
#ifdef LOCAL
	freopen("d:\\data.in", "r", stdin);
#endif
	while (scanf("%s", s))
	{
		n = 0;
		while(s[n])
		{
			r[n] = s[n]-'a'+1;
			n++;
		} // r其实就是s的rank数组
        r[n]=0; // 空串
		// 调用的结果是sa[1,...,n]中保存的就是sa, 其中sa的每个值的范围是 [0,...,n-1], 注意, sa[0]一定是n, 即空串, 就不考虑了
		dc3(r,sa,n+1, 27); // 27是26个小写英文字母+0（空字符）, 注意, 其实空字符可以换成是任何不在a-z中出现过的, 但是<a的字符（注意，一定要比a小, 不然例如aabaaaaba最后那个后缀a怎么排第一呢?）注意, 为什么要传入n+1? 因为正如上面说的那样, 我们要将空串这个后缀也考虑进去,所以r[n]、r[n+1]、r[n+2]都要是0, 传入n+1, 然后在dc3中令r[n]=r[n+1]就满足这一点了
	}
	return 0;
}

测试数据

1
2
3
aabaaaab
aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa
abracadabra

倍增和dc3之间的比较

所以dc3的性能理论上比倍增更好（但一般题目不会卡倍增, 但是有些变态题目如 poj 2406 Power Strings会卡）. 但是倍增的编码复杂度低于dc3, 所以自己酌情选择算法.

lcp

提及sa，就不能不提及lcp（longest common prefix）,因为sa+lcp将成为一个更加有力的工具. 以至于在解决实际问题方面能真正的和后缀树媲美

所谓lcp数组指的是后缀数组中相邻2个后缀之间的最长公共前缀组成的数组. 有的文章将lcp称之为”高度数组”.

即

1
2
lcp[i] = S[sa[i],...]和S[sa[i+1],...]的最长公共前缀. 0<=i<n-1,即每个后缀和它后面一名的后缀之间
的LCP， 一共n个后缀,即后缀数组一共n个, 所以lcp数组有n-1个元素

我们可以在O(n)时间内使用类似于尺取的技法高效求出lcp.

因为lcp的定义, 所以lcp的代码显然要引入sa的结果. 而且鉴于da和dc3返回的sa都是存储在[1,…,n], 而且取值在[0,…,n-1]范围内. 所以随便用da或者dc3都行的.

那么O(n)解决lcp问题的契机在哪里呢? 在于如下观察

从S[i,….], i=0 开始计算lcp. 即计算 lcp[r[i]-1] = S[i,…] 和 S[sa[r[i]-1],….]的LCP, for r[i]>=1（排名>=1的后缀）

假设我们已经计算出了S[i,…]和S[sa[r[i]-1],…]的LCP了. 那么考虑S[i+1,…]和排在它前面一名的后缀之间的LCP

令 S[k,…]=S[sa[r[i]-1],…], 即我们已经知道了 S[i,…]和S[k,…]之间的LCP———— 假设是h. 则显然S[i+1]和S[k+1,…]之间的LCP>=h-1（因为去掉了第一个字符嘛），虽然S[k+1]未必是S[i+1]在sa中前一名的后缀. 但是S[i+1]前一名的后缀和S[i+1]的LCP只会比S[k+1…]和S[i+1…]之间的LCP大, 而不会短. 毕竟S[i+1,…]前面一名和S[i+1]距离最近啊~

一图胜千言

从左往右第一个abra… 就是 S[i,…], 第二个 abra… 就是在sa中排在 S[i,…]前面的那个。从左往右第一个bra… 就是 S[i+1,…], 第二个bra… 就是S[k+1,…], 而如果S[k+1,…]排在S[i+1,…]前面的话, 则自然不如排在S[i+1,…]前一名和S[i+1,…]的LCP长——毕竟没有S[i+1,…]前一名和S[i+1,…]距离来的近乎. 如果S[k+1,…]排在S[i+1,…]的后面, 则口胡…. 因此每次考察只需要从h-1开始考察下一组相邻后缀的LCP就行了. 而LCP至多为n, 所以复杂度是O(n)的（因为每次计算LCP就是把h推高，但是最多推到n, 所以复杂度为O（n））.

说了这么多, 下面上美味的lcp的代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
#include "stdafx.h"
#include <stdio.h>
#include <string.h>
#include <algorithm>
using namespace std;
#define LOCAL

const int maxn = 1e5+5;
int n,k, r[maxn], t[maxn], sa[maxn],lcp[maxn];
char s[maxn];

bool c0(int i, int j)
{
	if (r[i]!=r[j])
	{
		return r[i]<r[j];
	}
	int ri = i+k<=n?r[i+k]:-1, rj = j+k<=n?r[j+k]:-1;
	return ri<rj;
}

void da()
{
	k = 1;
	n = 0;
	while(s[n])
	{
		sa[n] = n;
		r[n] = s[n]-'a'+1;
		n++;
	} 
	r[n]=0,sa[n] = n; 
	while(k<=n)
	{
		sort(sa, sa+n+1, c0); 
		t[sa[0]] = 0;
		for (int i = 1;i<=n;i++) 
		{
			t[sa[i]] = t[sa[i-1]]+c0(sa[i-1],sa[i]); 
		}
		memcpy(r,t,sizeof(int)*(n+1));
		k<<=1;
	}
}

void height() // 计算字符串s的lcp, 存储在lcp[1,...,n-1]中,lcp[i]表示S[sa[i],...]和S[sa[i+1],...]之间的LCP（即排名第i和排名第i+1的后缀之间的LCP）, 取值范围是[0,n-1]
{
	int h =0;
	for (int i = 0,j; i<n; i++) // 这里计算的策略是找到每个后缀前一名, 然后计算他俩的lcp长度, 进而为前一名的lcp数组值赋值
	{
		if (r[i]==1) // 第一名没有前一名
		{
			continue;
		}
		j = sa[r[i]-1]; // 排在S[i,...]前一名的后缀是S[j,...]
		if (h)
		{
			h--; // 从h-1开始检查,而不需要每次都从0开始检查
		}
		while(i+h<n && j+h<n)
		{
			if (s[i+h]!=s[j+h])
			{
				break;
			}
			h++;
		}
		lcp[r[i]-1] = h;
	}
}

int main() {
#ifdef LOCAL
	freopen("d:\\data.in", "r", stdin);
#endif
	while(~scanf("%s",s))
	{
		da(); // 计算sa完成之后, r[0,...,n-1]中保存了rank数组, 即r[i]是S[i,...]在后缀数组中的排名, r[i]的范围是[1,...,n]
		height(); // 计算lcp
	}
	return 0;
}

测试数据

1
2
3
aabaaaab
aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa
abracadabra

任意两个后缀之间的LCP

假设rank[i]<rank[j], 则S[i,…]和S[j,…] 的LCP是

lcp[rank[i]], lcp[rank[i]+1] ,…., lcp[rank[j]-1]

之间的RMQ（最小）问题.

参考

【1】https://www.cnblogs.com/jianglangcaijin/p/6035937.html

####