后缀数组的倍增算法(Prefix Doubling)

2017-01-08 08:19:29来源:cnblogs.com作者:janmelin人点击

第七城市

后缀数组的倍增算法(Prefix Doubling)


文本内容除特殊注明外,均在知识共享署名-非商业性使用-相同方式共享 3.0协议下提供,附加条款亦可能应用。


最近在自学习BWT算法(Burrows-Wheeler transform),其中涉及到对字符串循环移位求编码。直观的办法就是模拟,使用O(n3)的时间求出BWT编码。经过简单的简化后也要O(n2logn)的时间,显然当字符串长度很大时这种方法的效率很低。
由于循环移位的结果类似后缀(二者有所不同,所以在字符串结尾添加了一个字典序严格小于所有字符的符号,例如'/0',使得循环移位的有效部分等同于后缀),因此可以使用后缀树或后缀数组的方式优化BWT的过程。


关于学习倍增算法,你应该:

  1. 理解朴素的后缀数组生成方法
  2. 理解基数排序(本文使用了基数排序,至于原始的倍增算法是否是使用基数排序本人也不清楚)
  3. 了解KMP算法的原理

先来谈谈KMP算法。它之所以能有效减少比对次数是因为它利用了之前比对的结果——利用前缀的自相似性跳过必然失败的匹配,直接进行有可能成功的尝试。
而倍增算法同样拥有类似的思想,例如cake拥有后缀

cakeakekee

当我们比较了每个后缀第一个字母后(2nd 1st 4th 3rd),实际上我们也知道了每个后缀的第二个字母的比较结果(1st 4th 3rd -)。类似的,后续结果也就知道了。因此,我们可以得到逐步扩展每个后缀的前缀比较结果(2 1 4 3)->(21 14 43 3-)->(214 143 43- 3--)……参考图1(本例和图中所示不同,但思路是一样的)
上述延伸过程是线性增加的。若是再贪心一点,则可以利用上一回的比较结果将该回的前缀比较长度增加一倍,即指数级增长。这也就是倍增算法的核心思路。


图1引用自NOCOW

再来谈谈利用基数排序的算法实现。基数排序分为LSD(Least significant digital)和MSD(Most significant digital)两大类。乍一看后缀数组的比较是从高位开始的(p.s. 为什么不从低位开始呢?删除一个整数而不改变相对大小关系很简单,但添加一个整数而不改变相对大小关系比较麻烦),很适合MSD。但MSD的时间开销随序列复杂度和长度增长很快,仅适用于短序列,所以LSD是个无奈之选。可以说,倍增算法的代码之所以晦涩很大一部分原因就是使用LSD的缘故。


code in C++

#include<stdio.h>#include<string.h>#define rank r_saconst int MAXN=21;char str[MAXN];int sa[MAXN];//suffix arrayint l_sa[MAXN];//low of saint r_sa[MAXN];//reverse mapping of sa, also known as rank arrayint t_r_sa[MAXN];//temperary copy of r_sachar BWT[MAXN];int c[MAXN+128];//数组长度必须大于字符串长度和字符总数的最大值bool sa_cmp(int *r,int sa1,int sa2,int j){        //此处完美地进行了越界判断        return (r[sa1]==r[sa2] && r[sa1+j]==r[sa2+j]);}int prefixdouble(char *s,int l){        int i,j,k,m;                //对后缀的第一位进行基数排序        memset(c,0,sizeof(c));        for(i=0;i<l;i++)                c[ s[i] ]++;        for(m=1;m<MAXN+128;m++)                c[m]+=c[m-1];        for(i=l-1;i>=0;i--)                sa[ --c[s[i]] ]=i;        //r_sa[i]=k 即第i个后缀排名第k        for(i=0;i<l;i++)                r_sa[i]=s[i];//此时仅需反映相对大小顺序        int p;        for(j=1;j<=l;j*=2){                                //由于采用LSD,先对低位进行排序                p=0;                //l_sa[k]=i 即排名第k的是第i个后缀                for(i=l-j;i<l;i++)                    l_sa[p++]=i;//长度小于j的后缀无低位关键字,直接排在最前                for(k=0;k<l;k++)                    if(sa[k]>=j) l_sa[p++]=sa[k]-j;//第i-j个后缀的低位关键字等于第i个后缀的高位关键字,并且高位关键字在之前已有序                //再对高位进行排序                memset(c,0,sizeof(c));                for(k=0;k<l;k++)                        c[ r_sa[ l_sa[k] ] ]++;                for(m=1;m<MAXN+128;m++)                        c[m]+=c[m-1];                for(k=l-1;k>=0;k--)                        sa[ --c[ r_sa[ l_sa[k] ] ] ]=l_sa[k];                //更新r_sa                memcpy(t_r_sa,r_sa,4*MAXN);                r_sa[ sa[0] ]=p=0;                //相邻后缀如果前缀相同,那么其rank也相同                for(k=1;k<l;k++)                        r_sa[sa[k]]=sa_cmp(t_r_sa,sa[k-1],sa[k],j)?p:++p;                if(p==l-1) break;        }/*test        for(k=0;k<l;k++)                printf("%2d:%s/n",k,s+sa[k]);        BWT[0]=s[l-2];        for(i=1;i<l;i++)                BWT[i]=s[sa[i]-1];        BWT[l]='/0';        printf("trans:");        for(i=0;i<l;i++)                printf("%c",BWT[i]);*/}int main(){        printf("The string inputed should short than 20 symbols./n");        scanf("%s",str);        int l=strlen(str);        prefixdouble(str,l+1);        for(int i=0;i<l;i++)               printf("%d ",sa[i]);        return 0;}



第七城市

最新文章

123

最新摄影

微信扫一扫

第七城市微信公众平台