【后綴數組真難懂啊啊……就20+行的代碼搞了幾天才理解……不知是不是我太沙茶了】
【1】一些定義:
字符串:廣義的字符串是指“元素類型有序,且元素值有一定范圍的序列”,其元素不一定非要是字符,可以是數字等,因此整數、二進制數等也是字符串;
字符集:字符串的元素值的范圍稱為字符集,其大小記為SZ。
字符串的長度:字符串中元素的個數,一般記為N,長度為N的字符串A第一次提到時一般用A[0..N-1]來表示;
前綴:字符串A[0..N-1]的從A[0]開始的若干個連續的字符組成的字符串稱為A的前綴,以下“前綴i”或者“編號為i的前綴”指的都是A[0..i];
后綴:字符串A[0..N-1]的到A[N-1]終止的若干個連續的字符組成的字符串稱為A的后綴,以下“后綴i”或者“編號為i的后綴”指的都是A[i..N-1];
對于一個長度為N的字符串,將其N個后綴按字典序大小進行排序,得到兩個數組sa[i]和rank[i],sa[i]為排在第i位的后綴的編號(也就是一般說的ord[i]),rank[i]為排在后綴i排在的位置(稱為后綴i的名次)。sa、rank值的范圍均為[0..N-1]。sa和rank互逆,即sa[i]=j等價于rank[j]=i,或者說成sa[rank[i]]=rank[sa[i]]=i。這里,sa稱為
后綴數組,rank稱為
名次數組。
【2】用倍增算法求后綴數組:
在論文里,后綴數組有兩種求法:倍增算法和DC3算法,前者的時間復雜度為O(NlogN),但常數較小,后者的時間復雜度為O(N),但常數較大,在實際應用中,兩者的總時間相差不大,且后者比前者難理解得多(本沙茶理解前者都用了幾天時間……后者就木敢看了)。這里就總結一下倍增算法吧囧……
首先,貼一下本沙茶的用倍增算法求后綴數組的模板:
void suffix_array()
{
int p, v0, v1, v00, v01;
re(i, SZ) S[i] = 0;
re(i, n) rank[i] = A[i];
re(i, n) S[A[i]]++;
re2(i, 1, SZ) S[i] += S[i - 1];
rre(i, n) sa[--S[A[i]]] = i;
for (int j=1; j<n; j<<=1) {
p = 0; re2(i, n-j, n) tmp[p++] = i;
re(i, n) if (sa[i] >= j) tmp[p++] = sa[i] - j;
re(i, SZ) S[i] = 0;
re(i, n) S[rank[i]]++;
re2(i, 1, SZ) S[i] += S[i - 1];
rre(i, n) sa[--S[rank[tmp[i]]]] = tmp[i];
tmp[sa[0]] = p = 0;
re2(i, 1, n) {
v0 = sa[i - 1]; v1 = sa[i];
if (v0 + j < n) v00 = rank[v0 + j]; else v00 = -1;
if (v1 + j < n) v01 = rank[v1 + j]; else v01 = -1;
if (rank[v0] == rank[v1] && v00 == v01) tmp[sa[i]] = p; else tmp[sa[i]] = ++p;
}
re(i, n) rank[i] = tmp[i];
SZ = ++p;
}
}
這里A是待求sa和rank的字符串。
<1>倍增算法的思想:
記R[i][j]為A[i..i+2
j-1](如果越界,則后面用@填充)在A的所有長度為2
j的子串(越界則后面用@填充)中的名次(rank)值。倍增算法就是按階段求出所有R[i][j]的值,直到2
j>N為止。首先,R[i][0]的就是字符A[i]在A[0..N-1]中的名次,是可以直接用計數排序來實現的。然后,若R[0..N-1][j-1]已知,則可以按照以下方法求出R[0..N-1][j]的值:對每個i(0<=i<N),構造一個二元組<X
i, Y
i>,其中X
i=R[i][j-1],Y
i=R[i+2
j][j-1](若i+2
j>=N,則Y
i=-∞),然后對這N個二元組按照第一關鍵字為X,第二關鍵字為Y(若兩者都相等則判定為相等)進行排序(可以用基數排序來實現),排序后,<X
i, Y
i>的名次就是的R[i][j]的值。
<2>一開始,對A中的各個字符進行計數排序:
re(i, SZ) S[i] = 0;
re(i, n) rank[i] = A[i];
re(i, n) S[A[i]]++;
re2(i, 1, SZ) S[i] += S[i - 1];
rre(i, n) sa[--S[A[i]]] = i;
這個木有神馬好說的,在搞懂了基數排序之后可以秒掉。唯一不同的是這里加了一句:rank[i]=A[i],這里的rank[i]是初始的i的名次,MS不符合rank[i]的定義和sa與rank間的互逆性。這里就要解釋一下了囧。因為在求sa的過程中,rank值可能不符合定義,因為長度為2
j的子串可能會有相等的,此時它們的rank值也要相等,而sa值由于有下標的限制所以不可能有相等的。因此,在過程中,rank其實是用來代替A的子串的,這樣rank值只需要表示一個“相對順序”就行了,也就是:rank[i0]>(=, <)rank[i1],當且僅當A[i0..i0+2
j-1]>(=, <)A[i1..i1+2
j-1]。這樣,可以直接將A[i]值作為初始的rank[i]值。
<3>j(代替2
j)的值從1開始不斷倍增,對二元組進行基數排序求出新階段的sa值:
for (int j=1; j<n; j<<=1) {
p = 0; re2(i, n-j, n) tmp[p++] = i;
re(i, n) if (sa[i] >= j) tmp[p++] = sa[i] - j;
re(i, SZ) S[i] = 0;
re(i, n) S[rank[i]]++;
re2(i, 1, SZ) S[i] += S[i - 1];
rre(i, n) sa[--S[rank[tmp[i]]]] = tmp[i];
注意這個基數排序的過程是很特別的。首先,它并不是對A在進行排序,而是對上一階段求出的rank在進行排序。因為前面已經說過,在求sa的過程中,rank就是用來代替A的對應長度的子串的,由于不能直接對子串進行排序(那樣的話時間開銷很恐怖的),所以只能對rank進行排序。另外,這里在對二元組<x, y>的第二關鍵字(y)進行排序的過程中加了優化:這些y其實就是把上一階段的sa整體左移了j,右邊空出的部分全部用@(空串)填充得到的,由于空串的字典序肯定最小,因此將右邊的空串按照下標順序先寫入臨時sa(代碼中用tmp表示的就是臨時sa,也就是對第二關鍵字y排序后的ord結果),然后,上一階段的sa如果左移后還木有消失的(也就是sa值大于等于j的),再按順序寫入臨時sa,就得到了排序結果。剩下的對x的排序結果就是上一階段的sa,唯一不同的是對于x相同的,按照臨時名次遞增的順序。
<4>求出新階段的rank值:
tmp[sa[0]] = p = 0;
re2(i, 1, n) {
v0 = sa[i - 1]; v1 = sa[i];
if (v0 + j < n) v00 = rank[v0 + j]; else v00 = -1;
if (v1 + j < n) v01 = rank[v1 + j]; else v01 = -1;
if (rank[v0] == rank[v1] && v00 == v01) tmp[sa[i]] = p; else tmp[sa[i]] = ++p;
}
re(i, n) rank[i] = tmp[i];
SZ = ++p;
由于下一階段需要使用本階段的rank值,因此在求出了本階段的sa值以后,需要求rank值。(代碼中的tmp起了臨時rank的作用,目的是節省空間)
因為sa值已經求出,因此只要依次掃描sa就可以得到rank值,唯一要做的工作就是找到哪些子串是相等的,它們的rank值應該相等,除此之外,rank值只要依次加1即可。判定相等的方法:只需判定rank[i]和rank[i+j]是否都對應相等即可。若rank[i+j]越界,用-∞(當然任何一個負數都行,代碼中用了-1)來表示。
最后還有一個優化:由于本階段的名次的范圍只有[0..p]這么多,下一階段的“字符集”(其實就是rank集)的大小SZ可以設為p+1,這樣可以省一些時間。
這樣后綴數組sa和名次數組rank就全部求完了。
以后還有一些更重要的東東就是AC自動機、后綴數組等的應用問題,算了,以后再搞吧囧。