后缀数组 倍增算法代码详解

关于后缀数组，在刘汝佳的蓝书上已经有过介绍。关于算法的流程是比较详细的。但不过关于求sa数组的代码看起来就有一定的难度。

我看了很多博客，又在B站看了几遍视频，终于有一个大致的框图了，至于后缀数组是拿来干什么的就不要问我了。

首先关于算法流程这张经典的图

相信大家看过很多遍了吧

现在来看看代码

首先讲几个数组;

sa数组：在完成后缀数组前，表示关键字的排名的位置，下标表示名次，数组的值表示关键字的首字符的位置（0到n-1)，遇到值相同的情况，名次由他们在原串中相对位置的先后决定。在完成后缀数组后，下标就表示名次对应的后缀的首字符的位置了。注意这个数组并不是上面的算法图中rank值对应的。

x数组：表示一个rank，下标表示关键字的位置，数组的值，表示关键字的大小，相同的值有相同的rank。初始化为字符串的每个字符大小，每次迭代后，根据sa的值重新赋值。

y数组：表示上一次排序后组合关键字的第二关键字的数组，下标表示第二关键字名次，数组的值表示组合关键字的首字符的位置。

        for (int i = 0; i<m; i++) c[i] = 0;
	for (int i = 0; i<n; i++) 
		c[x[i] = s[i] - 'a' + 1]++;
	for (int i = 1; i<m; i++) c[i] += c[i - 1];
	for (int i = n - 1; i >= 0; i--) 
		sa[--c[x[i]]] = i;

这一段就是一个简单的基数排序，或者说计数排序吧，其中c数组表示桶吧，这段代码不难的，需要注意的就是先根据x的值排序，x值相等时根据出现的先后顺序排序的。这样第一次的sa数组就完成了。

        for (int k = 1; k <= n; k <<= 1)

开始进入迭代，根据倍增思想，k每一次翻倍。

	int p = 0;
	for (int i = n - k; i<n; i++) y[p++] = i;
	for (int i = 0; i<n; i++) 
	if (sa[i] >= k) 
		y[p++] = sa[i] - k;

这段代码开始就有一点难懂了，首先第一个for循环的意思，因为每次迭代，总会出现有的后缀第二关键字不存在的，这些后缀的第二关键字，我们都认为是最小的，因此在根据y数组的定义，他就是让这些后缀出现在前面。

第二个for循环，其他的后缀的第二关键字是可以用其他后缀的第一关键字表示的，因为sa数组表示的地址，因此他表示的组合关键字他当前的sa数组的值减去k，如果他的sa数组的值小于k，则他不能表示为其他后缀的第二关键字。

	for (int i = 0; i<m; i++) c[i] = 0;
	for (int i = 0; i<n; i++) c[x[y[i]]]++;
	for (int i = 0; i<m; i++) c[i] += c[i - 1];
	for (int i = n - 1; i >= 0; i--) 
		sa[--c[x[y[i]]]] = y[i];

这段代码也不好理解，看起来和最初的计数排序差不多，关键就是弄懂那个x[y[i]]的所表达的意思；

我们可以猜到代码就是排序，猜嘛。肯定就是合并第一第二关键字排序了。让我们回忆一下基数排序的特点，就是先对x的值排序，x的值相等时根据出现的先后顺序进行排序，x里面存储了上一次关键字的排序，然后就是对第二关键字进行排序了，这里x[y[i]]，根据y的值改变了第一关建字的排序，所以对x[y[i]]进行排序就是对第一二关键字，进行排序。

其实就是y数组先对第二关键字进行了排序，然后在对第一关键字进行排序。

        swap(x, y);
	p = 1; x[sa[0]] = 0;
	for (int i = 1; i<n; i++)
		x[sa[i]] = y[sa[i - 1]] == y[sa[i]] && y[sa[i - 1] + k] == y[sa[i] + k] ? p - 1 : p++;
	if (p >= n)
		break;
	m = p;

每次更新rank数组x，因为有了关键字排序数组sa,要得到rank数组也很容易，其中那个三目表达式的前面部分，是为了判断两个合并后的后缀是否完全相同。没有别的意思。最后当p>=n的时候表明已经排好序了。后面的下一次桶的大小为p。

完整代码：

#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
using namespace std;
const int N = 100010;
char s[N];
int sa[N], c[N], x[N], y[N], n, m, tmp[N];

void Suffix()
{
	for (int i = 0; i<m; i++) c[i] = 0;
	for (int i = 0; i<n; i++)
		c[x[i] = s[i] - 'a' + 1]++;
	for (int i = 1; i<m; i++) c[i] += c[i - 1];
	for (int i = n - 1; i >= 0; i--)
		sa[--c[x[i]]] = i;
	for (int k = 1; k <= n; k <<= 1)
	{
	    for (int i = 0; i<n; i++) cout << sa[i] << " ";
        cout << endl;
		int p = 0;
		for (int i = n - k; i<n; i++) y[p++] = i;
		for (int i = 0; i<n; i++)
		if (sa[i] >= k)
			y[p++] = sa[i] - k;

		for (int i = 0; i<m; i++) c[i] = 0;
		for (int i = 0; i<n; i++) c[x[y[i]]]++;
		for (int i = 0; i<m; i++) c[i] += c[i - 1];
		for (int i = n - 1; i >= 0; i--)
			sa[--c[x[y[i]]]] = y[i];

		swap(x, y);
		p = 1; x[sa[0]] = 0;
		for (int i = 1; i<n; i++)
			x[sa[i]] = y[sa[i - 1]] == y[sa[i]] && y[sa[i - 1] + k] == y[sa[i] + k] ? p - 1 : p++;
		if (p >= n)
			break;
		m = p;
	}
	for (int i = 0; i<n; i++) cout << sa[i] << " ";
	cout << endl;
}

int main()
{
	cin >> s;
	n = strlen(s);
	s[n++] = 'a' - 1;
	m = 30;
	Suffix();
	return 0;
}

今天贴一个封装好的吧

struct SA {
    int sa[maxn], ra[maxn], height[maxn];
    int t1[maxn], t2[maxn], c[maxn];
    void build(char *str, int n, int m) {
        str[n] = 0;
        n++;
        int i, j, p, *x = t1, *y = t2;
        for (i = 0; i < m; i++) c[i] = 0;
        for (i = 0; i < n; i++) c[x[i] = str[i]]++;
        for (i = 1; i < m; i++) c[i] += c[i - 1];
        for (i = n - 1; i >= 0; i--) sa[--c[x[i]]] = i;
        for (j = 1; j <= n; j <<= 1) {
            p = 0;
            for (i = n - j; i < n; i++) y[p++] = i;
            for (i = 0; i < n; i++) if (sa[i] >= j) y[p++] = sa[i] - j;
            for (i = 0; i < m; i++) c[i] = 0;
            for (i = 0; i < n; i++) c[x[y[i]]]++;
            for (i = 1; i < m; i++) c[i] += c[i - 1];
            for (i = n - 1; i >= 0; i--) sa[--c[x[y[i]]]] = y[i];
            swap(x, y);
            p = 1;
            x[sa[0]] = 0;
            for (i = 1; i < n; i++)
                x[sa[i]] = (y[sa[i - 1]] == y[sa[i]] && y[sa[i - 1] + j] == y[sa[i] + j]) ? p - 1 : p++;
            if (p >= n) break;
            m = p;
        }
        n--;
        for (int i = 0; i <= n; i++) ra[sa[i]] = i;
        for (int i = 0, j = 0, k = 0; i < n; i++) {
            if (k) k--;
            j = sa[ra[i] - 1];
            while (str[i + k] == str[j + k]) k++;
            height[ra[i]] = k;
        }
        st_init(height, n);
    }
    int lg[maxn], table[20][maxn];
    void st_init(int *arr, int n) {
        if (!lg[0]) {
            lg[0] = -1;
            for (int i = 1; i < maxn; i++)
                lg[i] = lg[i / 2] + 1;
        }
        for (int i = 1; i <= n; ++i)
            table[0][i] = arr[i];
        for (int i = 1; i <= lg[n]; ++i)
            for (int j = 1; j <= n; ++j)
                if (j + (1 << i) - 1 <= n)
                    table[i][j] = min(table[i - 1][j], table[i - 1][j + (1 << (i - 1))]);
    }
    int lcp(int l, int r) {
        l = ra[l], r = ra[r];
        if (l > r) swap(l, r);
        ++l;
        int t = lg[r - l + 1];
        return min(table[t][l], table[t][r - (1 << t) + 1]);
    }
} sa;