后綴排序
Task
Description
給定一個字符串,要求按字典序升序輸出它的所有后綴子串的第一個字符所在位置。
Requirements & Limitations
字符集大小為常數,要求時間復雜度 \(O(n \log n)\),其中 \(n\) 為字符串長度
Algorithm
這就是大(ren)名(lei)鼎(zhi)鼎(hui)的后綴排序了。這里記錄倍增法。
定義后綴 \(i\) 為首字符所在原字符串位置為 \(i\) 的后綴子串。
定義 \(sa_i\) 為排名為 \(i\) 的后綴的首字符所在位置,\(rnk_i\) 為后綴 \(i\) 的排名。也即對 \(sa_i\) 做排序相當於對原字符串做排序。
注意到 \(sa\) 和 \(rnk\) 是可以在 \(O(n)\) 時間內互相推出的,因為 排名為 (后綴 \(i\) 的排名) 的后綴為 \(i\),形式化的說,\(sa_{rnk_i} = i\),
考慮進行多次排序,每次排序固定一個長度,對所有后綴子串的該長度的前綴排序(長度不足補 \(0\))。我們倍增這個長度。
首先,我們對所有后綴子串的第一位進行排序,當字符相同時,我們要求首字符所在位置小的在前面,例如對於字符串 \(ababa\) 的所有后綴子串
1 ababa
2 baba
3 aba
4 ba
5 a
進行一次排序的結果為
1 ababa
3 aba
5 a
2 baba
4 ba
即我們 \(sa\) 數組當前的值為 \(\{1,~3,~5,~2,~4\}\),然后計算 \(rnk\) 數組,但是需要注意的是,兩字符串第一位相同時,他們的 \(rnk\) 也應該相同。例如,\(rnk_1 = rnk_3 = 1\)。類似的得到 \(rnk\) 數組當前的值為 \(\{1,~2,~1,~2,~1\}\)。
然后進行下一次排序。
由於第一次排序的長度為 \(1\),所以本次排序的長度為 \(1 \times 2 = 2\)
對於上次的排序結果,我們它們的前兩位進行排序。但是注意到一個性質:我們已經排好了第一位,現在對第二位進行大小比較,我們發現后綴 \(i\) 的第 \(2\) 位正是后綴 \((i + 1)\) 的第 \(1\) 位,因此第二位的字符的相對大小關系我們也是知道的,具體的,比較 \((sa_i,~sa_j)\),若 \(rnk_{sa_i} \neq rnk_{sa_j}\),則他們的大小關系為 \(rnk_{sa_i}\) 和 \(rnk_{sa_j}\) 的大小關系。否則則說明它們的第一位是相同的,那么比較第二位,則它們的大小關系是 \(rnk_{sa_i + 1}\) 與 \(rnk_{sa_j + 1}\) 的大小關系。
排序以后的結果為
5 a
1 ababa
3 aba
4 ba
2 baba
即 \(sa\) 數組當前的值為 \(\{5,~1,~3,~4,~2\}\),而 \(rnk\) 當前的值為 \(\{2,~3,~2,~3,~1\}\)
然后進行下一次排序。
由於上次的排序長度為 \(2\),所以本次排序的長度為 \(2 \times 2 = 4\)
類似於上次的排序,我們已經知道了前兩位的相對大小關系,現在要求前四位的相對大小關系,而后綴 \(i\) 的后兩位正是后綴 \((i + 2)\) 的前兩位,在上次排序已經被排好了。於是他們后兩位字符的相對大小關系也是知道的。
類似的,得到排序的結果為
5 a
3 aba
1 ababa
4 ba
2 baba
即 \(sa = \{5,~3,~1,~4,~2\}\),\(rnk = \{3,~5,~2,~4,~1\}\)
至此, \(rnk\) 數組已經是一個 \(1 \sim n\) 的排列,這意味着所有的后綴已經都排好了序,輸出即可。
形式化的,我們倍增所排序的長度,設當前需要排序的長度為 \(2len\),已經對 \(len\) 長度排好了序,且 \(sa\) 和 \(rnk\) 存儲的是長度為 \(len\) 時的排序信息,那么對於兩個后綴 \((i, j)\),他們的前半部分大小關系為 \((rnk_i, rnk_j)\),而注意到他們后半部分即為后綴 \((i + len)\) (或后綴 \((j + len)\))的前半部分,因此他們后半部分的大小關系為 \((rnk_{i + len}, rnk_{j + len})\)。根據字符串比較規則,如果前半部分大小不同,則大小關系為前半部分大小關系,否則為后半部分大小關系。也就相當於對 \((rnk_i, rnk_{i + len})\) 這個二元組進行排序。一直到所倍增的長度大於 \(n\) 時,則所有后綴的大小已被排好。
Algorithm \(1\)
Solution
我會 std::sort
!
我們倍增長度,對二元組 \((rnk_i, rnk_{i + len})\) 排序時,直接使用 std::sort
。注意到我們倍增了 \(O(\log n)\) 次,每次需要進行一次 \(O(n \log n)\) 的排序,因此總時間復雜度 \(O(n \log^2 n)\)。但是依然能跑過你谷的模板題。
這個做法非常易於理解,因為這只是翻譯了上面推導過程,代碼也好記好寫。
Code
#include <cstdio>
#include <algorithm>
const int maxn = 2000005;
int n;
char S[maxn];
int rnk[maxn], sa[maxn];
std::pair<int, int> MU[maxn];
int ReadStr(char *p);
bool cmp(const int &a, const int &b);
int main() {
freopen("1.in", "r", stdin);
n = ReadStr(S);
for (int i = 1; i <= n; ++i) {
rnk[i] = 1;
sa[i] = i;
MU[i].second = S[i];
}
for (int len = 1; len <= n; len <<= 1) {
std::sort(sa + 1, sa + 1 + n, cmp);
for (int i = 1; i <= n; ++i) {
if (MU[sa[i]] == MU[sa[i - 1]]) {
rnk[sa[i]] = rnk[sa[i - 1]];
} else {
rnk[sa[i]] = i;
}
}
for (int i = 1; i <= n; ++i) {
MU[i].first = rnk[i];
MU[i].second = rnk[i + len];
}
}
for (int i = 1; i <= n; ++i) {
qw(sa[i], ' ', true);
}
putchar('\n');
return 0;
}
int ReadStr(char *p) {
auto beg = p;
do *(++p) = IPT::GetChar(); while (((*p >= 'a') && (*p <= 'z')) || ((*p >= '0') && (*p <= '9')) || ((*p >= 'A') && (*p <= 'z')));
*p = 0;
return p - beg - 1;
}
inline bool cmp(const int &a, const int &b) {
if (MU[a] != MU[b]) {
return MU[a] < MU[b];
} else {
return a < b;
}
}
Algorithm \(2\)
考慮是對雙關鍵字進行排序,因此我們可以先對第二關鍵字做桶排序,然后對第一關鍵字做桶排序,在排序時,記錄每個桶所壓入的序號序列,取出是按照序列的順序取出每個下標,作為排好序的序列。由於桶排序是穩定的排序,這樣在第一關鍵字相同的時候先取出的是第二關鍵字較小的下標。這樣我們得到了一個 \(O(n)\) 的基數排序,但是由於 std::vector
常數太大,這樣的寫法甚至比 \(O(n \log^2 n)\) 的做法還要慢一倍,在你谷的板板題上會T三個點。
這個做法只是用這個大常數基數排序替換上面的 std::sort
,雖然很慢但是他的復雜度已經正確了。
Code
其余部分相同,只是將 std::sort
換成 RadixSort
,排序部分代碼如下
std::vector<int>bk[maxn];
void RadixSort(int *const beg, int *const ed) {
for (auto it = beg; it != ed; ++it) {
bk[MU[*it].second].push_back(*it);
}
auto p = beg;
for (int i = 0; i <= n; ++i) {
for (auto u : bk[i]) {
*(p++) = u;
}
bk[i].clear();
}
for (auto it = beg; it != ed; ++it) {
bk[MU[*it].first].push_back(*it);
}
p = beg;
for (int i = 0; i <= n; ++i) {
for (auto u : bk[i]) {
*(p++) = u;
}
bk[i].clear();
}
}
Algorithm \(3\)
考慮使用人類智慧優化我們的基數排序部分。
我們再定義一個數組 \(tp\),\(tp_i\) 代表 第二關鍵字 為 \(i\) 時字符串的位置。
(注意我們下方的假設為已經對 \(len\) 做好了排序,現在要排序 \(2len\))
例如,\(len = 2,~rnk_3 = 2\),則有 \(tp_2 = 1\),因為后綴 \(1\) 的第二關鍵字為 \(rnk_{1 + 2 = 3} = 2\),因此 \(tp_2 = 1\)。
形式化的,當 \(i + len \leq n\) 時,\(tp_{rnk_{i + len}} = i\)。
而當 \(i + len > n\) 時,我們發現這一部分的 \(i\) 對應的第二關鍵字都是 \(0\)(因為后半部分是空串),比不是零的都要小,考慮到要求字典序相同時序號小的在前,我們只需要令 \(tp_{i + len - n} = i\) 即可。
當然,賦值時我們考慮枚舉 \(sa\) 的值而不是 \(rnk\) 的(因為初始 \(rnk\) 的賦值會是字符集,需要一些額外處理),於是有
\(\forall i \leq len,~tp_i = n - len + 1\),這一部分是對應第二關鍵字為 \(0\) 的部分。
\(\forall sa_i > len,~~tp_{++pos} = sa_i - len\)。其中 \(pos\) 為計數器,初值為 \(len\)。注意我們是從小到大枚舉 \(i\) 而不是枚舉 \(sa_i\)。
現在用另一種方法解釋第二個式子:
我們從小到大枚舉了 \(sa_i\) ,即我們按照排名從小到大枚舉了每個字符串,那么他們作為第二關鍵字的值當然是從小到大的(因為第二關鍵也是排名),所以我們按順序記錄這些位置,而這些位置作為第二關鍵字,對應第一關鍵字的位置是 \(sa_i - len\),因為當前需要排序的長度為 \(2len\),所以第一關鍵字與第二關鍵字的位置相差 \(len\)。
上面兩句求 \(tp\) 的方法,是人類智慧基數排序的關鍵前提條件。
現在考慮進行基數排序。
基數排序的代碼如下
void RadixSort() {
for (int i = 0; i <= m; ++i) tax[i] = 0;
for (int i = 1; i <= n; ++i) ++tax[rnk[i]];
for (int i = 1; i <= m; ++i) tax[i] += tax[i - 1];
for (int i = n; i; --i) sa[ tax[rnk[tp[i]]]-- ] = tp[i];
}
\(tax\) 代表桶,首先將桶清零,然后枚舉每個字符串,將字符串的排名放入桶內。
然后對桶做一個前綴和。我們發現我們現在放入桶內的是第一關鍵字,這意味着從這個桶 \(tax_i\) 里出來的字符串所在的排名定大於 \(\sum_{x = 1}^{i - 1} tax_x\),因為這部分的第一關鍵字小於 \(i\),他們的排名一定小於第一關鍵字為 \(i\) 的字符串。同理,他們的排名一定不大於 \(\sum_{x = 1}^i tax_x\)。
然后進入第四行,我們倒序枚舉第二關鍵字的值,用 \(tp\) 來確定這個第二關鍵字是哪個字符串的,然后用 \(rnk\) 找到這個字符串對應哪個桶。這個桶即為 \(tax_{rnk_{tp_i}}\)。然后我們注意到從這個桶里出來的元素下標上界為當前的前綴和值,而我們是倒序枚舉的第二關鍵字,所以它的排名應該是桶里剩下字符串中第二關鍵字最大也即排名最大的,因此它的排名就是 \(tax_{rnk_{tp_i}}\)。完成后,將這個桶的值自減 \(1\),代表將它的上界減一。這么做下去就可以完成基數排序。
好人類智慧啊QAQ
Code
#include <ctime>
#include <cstdio>
#include <cstring>
#include <vector>
#include <algorithm>
int n, m;
char S[maxn];
int rnk[maxn], sa[maxn], tp[maxn], tax[maxn];
void RadixSort();
int ReadStr(char *p);
bool cmp(const int &a, const int &b);
int main() {
freopen("1.in", "r", stdin);
n = ReadStr(S);
for (int i = 1; i <= n; ++i) {
rnk[i] = S[i];
tp[i] = i;
}
m = 1000;
RadixSort();
for (int len = 1, p = 0; p < n; len <<= 1, m = p) {
p = 0;
for (int i = 1; i <= len; ++i) tp[++p] = n - len + i;
for (int i = 1; i <= n; ++i) if (sa[i] > len) tp[++p] = sa[i] - len;
RadixSort();
std::swap(tp, rnk);
rnk[sa[1]] = p = 1;
for (int i = 2; i <= n; ++i) {
rnk[sa[i]] = ((tp[sa[i - 1]] == tp[sa[i]]) && (tp[sa[i - 1] + len] == tp[sa[i] + len])) ? p : ++p;
}
}
for (int i = 1; i <= n; ++i) {
qw(sa[i], ' ', true);
}
putchar('\n');
return 0;
}
int ReadStr(char *p) {
auto beg = p;
do *(++p) = IPT::GetChar(); while (((*p >= 'a') && (*p <= 'z')) || ((*p >= '0') && (*p <= '9')) || ((*p >= 'A') && (*p <= 'z')));
*p = 0;
return p - beg - 1;
}
void RadixSort() {
for (int i = 0; i <= m; ++i) tax[i] = 0;
for (int i = 1; i <= n; ++i) ++tax[rnk[i]];
for (int i = 1; i <= m; ++i) tax[i] += tax[i - 1];
for (int i = n; i; --i) sa[ tax[rnk[tp[i]]]-- ] = tp[i];
}