suffix-array

發表於 更新於 Disqus：

用途

考慮”ababcad”
那麼就可以分成7個這樣的字串：

1
2
3
4
5
6
7
0 - "ababcad"
1 - "babcad"
2 - "abcad"
3 - "bcad"
4 - "cad"
5 - "ad"
6 - "d"

然後再把這7個字串按照ASCII碼的順序排序，變成這樣：

1
2
3
4
5
6
7
0 - "ababcad"
2 - "abcad"
5 - "ad"
1 - "babcad"
3 - "bcad"
4 - "cad"
6 - "d"

所以順序就是0251346

算法

考慮倍增
先把字串在後面補個極小字元(假設$)，更後面的就亂擺(方便起見就讓它變成環吧)，變成

1
2
3
4
5
6
7
8
0 - "ababcad$"
1 - "babcad$a"
2 - "abcad$ab"
3 - "bcad$aba"
4 - "cad$abab"
5 - "ad$ababc"
6 - "d$ababca"
7 - "$ababcad"

可以發現大小不會變，因為有一個極小字元。
然後先取第一個字元

1
2
3
4
5
6
7
8
0 - a
1 - b
2 - a
3 - b
4 - c
5 - a
6 - d
7 - $

直接用std::sort $\mathcal{O}(NlogN)$後再編號進rk陣列(可以想成離散化)。

1
2
3
4
5
6
7
8
7 - $ rk[7]=0
0 - a rk[0]=1
2 - a rk[2]=1
5 - a rk[5]=1
1 - b rk[1]=2
3 - b rk[3]=2
4 - c rk[4]=3
6 - d rk[5]=4

之後要怎麼從長度$2^k$轉移成$2^{k+1}$倍呢?
把他們給變成一個pair(rk[$i$], rk[$i + 2^k$])。

1
2
3
4
5
6
7
8
7 - $a rk[7]=0 (0, 1) // (rk[7], rk[0])
0 - ab rk[0]=1 (1, 2)
2 - ab rk[2]=1 (1, 2)
5 - ab rk[5]=1 (1, 3)
1 - ba rk[1]=2 (2, 1)
3 - bc rk[3]=2 (2, 3)
4 - ca rk[4]=3 (3, 1)
6 - d$ rk[6]=3 (4, 0)

然後使用std::sort排序他們。把排序的結果寫進rk

1
2
3
4
5
6
7
8
7 - $a rk[7]=0 (0, 1)
0 - ab rk[0]=1 (1, 2)
2 - ab rk[2]=1 (1, 2)
5 - ad rk[5]=2 (1, 3)
1 - ba rk[1]=3 (2, 1)
3 - bc rk[3]=4 (2, 3)
4 - ca rk[4]=5 (3, 1)
6 - d$ rk[6]=6 (4, 0)

再蓋pair

1
2
3
4
5
6
7
8
7 - $aba rk[7]=0 (0, 3) //(rk[7], rk[1])
0 - abab rk[0]=1 (1, 1)
2 - abca rk[2]=1 (1, 5)
5 - ad$a rk[5]=2 (2, 0)
1 - babc rk[1]=3 (3, 4)
3 - bcad rk[3]=4 (4, 2)
4 - cad$ rk[4]=5 (5, 6)
6 - d$ab rk[6]=6 (6, 1)

然後排序… 做logn次
最後得出的順序70251346，扣掉最左邊新加的，0251346即是答案。
複雜度:logn次中，每次排序都$\mathcal{O}(nlogn)$，所以是$\mathcal{O}(nlog^2n)$。

優化

  1. 因為只是要對pair來sort，且pair的值域只有$[0,n]$，因此可以先排序右邊再排序左邊(radix sort)，複雜度$O(nlogn)$。
  2. 考慮到pair(rk[$i$], rk[$i+2^k$])中，rk[$i$]已經排序過了，所以可以改成pair(rk[$i-2^k$], rk[$i$])，這樣只需要排序左邊(counting sort)，常數優化。

Code

回傳的tmp沒拔掉最左邊新加的，然後我直接拿TIOJ-1497的code來用

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
#include <bits/stdc++.h>
using namespace std;
#define EmiliaMyWife ios::sync_with_stdio(0); cin.tie(NULL);

/*-----------------------------------------------------------------------------------------------------*/
void csort(vector<int> &tmp, vector<int> &rk) {//counting sort
	int n = tmp.size();
	vector<vector<int>> has(n);
	for(int a : tmp)
		has[rk[a]].push_back(a);
	int it = 0;
	for(int i = 0; i < n; i++)
		for(int &a : has[i])
			tmp[it++] = a;
}

vector<int> solve(string &s) {
	s.push_back(0);//smallest character
	int n = s.size(), k = __lg(n) + 1;
	vector<int> rk(n), tmp(n), nrk(n);
	
	for(int i = 0; i < n; i++)
		tmp[i] = i;
	sort(tmp.begin(), tmp.end(), [&](auto a, auto b) {
		return s[a] < s[b];
	});
	
	for(int i = 0, cnt = 0; i < n; i++) {
		rk[tmp[i]] = cnt;
		if(i + 1 < n && s[tmp[i]] != s[tmp[i + 1]])
			cnt++;
	}
	for(int j = 0; j < k; j++) {
		for(int i = 0; i < n; i++)
			tmp[i] = (tmp[i] - (1 << j) % n + n) % n;
		csort(tmp, rk);
		for(int i = 0, cnt = 0; i < n; i++) {
			nrk[tmp[i]] = cnt;
			auto cur = make_pair(rk[tmp[i]], rk[(tmp[i] + (1 << j)) % n]);
			auto nxt = make_pair(rk[tmp[i + 1]], rk[(tmp[i + 1] + (1 << j)) % n]);
			if(i + 1 < n && cur != nxt)
				cnt++;
		}
		rk = nrk;
	}
	return tmp;
}
/*-----------------------------------------------------------------------------------------------------*/

signed main() {
	EmiliaMyWife

	string s;
	getline(cin, s);
	auto w = solve(s);
	for(int i = 1; i < s.size(); i++)
		cout << w[i] << '\n';

	return 0;
}

習題