原文地址:http://blog.csdn.net/xymscau/article/details/8798046
后缀数组号称字符串处理神器,不过发现好多人都只会用模板,其实这不是我们学算法的本质,我们学习算法的本质应该理解其实现原理,并加以实现,特别是算法,更讲究的是一种思想。一年前的我也是只会用别人的模板,最近却静下心来,研究了一下后缀数组,自己写了一份自己的模板。
我基本上是跟着连教的ppt来学习的,当然也少不了百度,先讲一下基本概念。(这里大量引用了连教的ppt)
基本定义:
子串 注:串!=字符串
字符串 S 的子串r[i..j] , i ≤ j ,表示r 串中从 i 到 j 这一段,就是顺次排列r[i],r[i+1],...,r[j] 形成的子串。
后缀
后缀是指从某个位置 i 开始到整个串末尾结束的一个特殊子串。字符串r 的从 第 i 个字 符 开 始 的 后 缀 表 示 为 Suffix(i) ,也 就 是Suffix(i)=r[i..len(r)] 。
后缀数组(SA[i]存放排名第i大的子串首字符下标)
后缀数组 SA 是一个一维数组,它保存1..n 的某个排列 SA[1] ,SA[2] , ……, SA[n] ,并且保证Suffix(SA[i])<Suffix(SA[i+1]), 1 ≤ i<n 。也就是将 S 的 n 个后缀从小到大进行排序之后把排好序的后缀的开头位置顺次放入SA 中。
名次数组(rank[i]存放suffix(i)的优先级)
名次数组 Rank[i] 保存的是 Suffix(i) 在所有后缀中从小到大排列的 “ 名次 ” 。
注:这个是排序的关键字~(这句话是我们排序的重点)
算法目标:
求得串的sa数组和rank数组
易知sa和rank互为逆操作,即sa[rank[i]] = i;
Rank[sa[i]] = i;(所以我们只要求得其一,就能O(n)算出另一个)
注:这个结论只在最后完成排序的时候符合。
但sa和rank的定义一直都是适用的。
原因是最后的时候不会存在相同(rank相等)的两个子串。
算法基本流程
1 //cnt是计数排序的辅助数组,k是第一关键字,id是第二关键字下标数组,r是以下标为第二关键字的新构数组,w存放的是字符串信息。sa保存的是排第i的是谁 2 int *k = rk,*id = height,*r = res, *cnt = wa;//计数排序 3 rep(i,up) cnt[i] = 0; 4 rep(i,len) cnt[k[i] = w[i]]++; 5 rep(i,up) cnt[i+1] += cnt[i]; 6 for(int i = len - 1; i >= 0; i--) { 7 sa[--cnt[k[i]]] = i; 8 }
求第二关键字(想想为什么构造w数组的时候末尾要加个0)
1 //cnt是计数排序的辅助数组,k是第一关键字,id是第二关键字下标数组,r是以下标为第二关键字的新构数组,w存放的是字符串信息,<span style="font-family: Arial, Helvetica, sans-serif;">sa保存的是排第i的是谁</span> 2 for(int i = len - d; i < len; i++) id[p++] = i; 3 rep(i,len) if(sa[i] >= d) id[p++] = sa[i] - d;//id保存了按后h/2排序的的序列,即排第i的后h/2的是原数组中的那个 4 rep(i,len) r[i] = k[id[i]];//构造新的排序数组
对新数组排序
1 //cnt是计数排序的辅助数组,k是第一关键字,id是第二关键字下标数组,r是以下标为第二关键字的新构数组,w存放的是字符串信息,sa保存的是排第i的是谁 2 rep(i,up) cnt[i] = 0; 3 rep(i,len) cnt[r[i]]++; 4 rep(i,up) cnt[i+1] += cnt[i]; 5 for(int i = len - 1; i >= 0; i--) { 6 sa[--cnt[r[i]]] = id[i]; 7 }
得到新的关键字(即按H长度排序后的离散序列)
1 //cnt是计数排序的辅助数组,k是第一关键字,id是第二关键字下标数组,r是以下标为第二关键字的新构数组,w存放的是字符串信息,sa保存的是排第i的是谁 2 swap(k,r); 3 p = 0; 4 k[sa[0]] = p++; 5 rep(i,len-1) { 6 if(sa[i]+d < len && sa[i+1]+d <len &&r[sa[i]] == r[sa[i+1]]&& r[sa[i]+d] == r[sa[i+1]+d]) 7 k[sa[i+1]] = p - 1; 8 else k[sa[i+1]] = p++; 9 }
1 #define rep(i,n) for(int i = 0;i < n; i++) 2 using namespace std; 3 const int size = 200005,INF = 1<<30; 4 int rk[size],sa[size],height[size],w[size],wa[size],res[size]; 5 void getSa (int len,int up) { 6 int *k = rk,*id = height,*r = res, *cnt = wa; 7 rep(i,up) cnt[i] = 0; 8 rep(i,len) cnt[k[i] = w[i]]++; 9 rep(i,up) cnt[i+1] += cnt[i]; 10 for(int i = len - 1; i >= 0; i--) { 11 sa[--cnt[k[i]]] = i; 12 } 13 int d = 1,p = 0; 14 while(p < len){ 15 for(int i = len - d; i < len; i++) id[p++] = i; 16 rep(i,len) if(sa[i] >= d) id[p++] = sa[i] - d; 17 rep(i,len) r[i] = k[id[i]]; 18 rep(i,up) cnt[i] = 0; 19 rep(i,len) cnt[r[i]]++; 20 rep(i,up) cnt[i+1] += cnt[i]; 21 for(int i = len - 1; i >= 0; i--) { 22 sa[--cnt[r[i]]] = id[i]; 23 } 24 swap(k,r); 25 p = 0; 26 k[sa[0]] = p++; 27 rep(i,len-1) { 28 if(sa[i]+d < len && sa[i+1]+d <len &&r[sa[i]] == r[sa[i+1]]&& r[sa[i]+d] == r[sa[i+1]+d]) 29 k[sa[i+1]] = p - 1; 30 else k[sa[i+1]] = p++; 31 } 32 if(p >= len) return ; 33 d *= 2,up = p, p = 0; 34 } 35 } 36 void getHeight(int len) { 37 rep(i,len) rk[sa[i]] = i; 38 height[0] = 0; 39 for(int i = 0,p = 0; i < len - 1; i++) { 40 int j = sa[rk[i]-1]; 41 while(i+p < len&& j+p < len&& w[i+p] == w[j+p]) { 42 p++; 43 } 44 height[rk[i]] = p; 45 p = max(0,p - 1); 46 } 47 } 48 int getSuffix(char s[]) { 49 int len = strlen(s),up = 0; 50 for(int i = 0; i < len; i++) { 51 w[i] = s[i]; 52 up = max(up,w[i]); 53 } 54 w[len++] = 0; 55 getSa(len,up+1); 56 getHeight(len); 57 return len; 58 }