字符串hash-RK算法讲解二
算法分析:预处理时间Θ(m),即求h,p,t的时间为,匹配时间在最坏情况下为Θ((n-m-1)m),因为可能出现每次都是可能命中点的情况。如T=a^n,P=a^m,此种情况下验证时间为Θ((n-m-1)m)。当然实际中,可能的命中点一般很少。假设有c个,则算法的期望匹配时间为O(n-m+1 +cm)=O(m+n),当m<<n时,期望匹配时间为O(n).
Rabin-karp算法是朴素字符串匹配算法的一个特例。当字母表∑为d进制数时,即∑={0,1,2,…d-1}。如当d=10时字母表中的每个字符都是一个十进制数。我们在比较两个长度为m的子串时,可以把这两个子串当作整数进行比较,而不用逐个字符比较,从而在某种程度上减少算法时间。
把一个由d进制数字组成的字符串转换成相应的十进制整数,这是大家曾经都写过的东西,一个可能的简单实现如下:
int ConvertToInt(const char* str, int d) { int ans =0; int m = strlen(str); for(int i=0; i<m; i++) { ans = ans*d + str[i]-'0'; } return ans; }
这就是所谓的霍纳法则(Horner’s rule)
对一个长度为m的子串,求其对应的整数值,其复杂度为O(m),如果不做一些特殊的处理,每次把子串和模式P比较时,先求子串对应的整数t,再与模式P对应的整数值p比较,算法的复杂度并没有得到改进。如果用t(s)表示当前位移下子串T[s+1..s+m]的值,我们注意到把模式向右滑动一个窗口之后,下一个m位的子串T[s+2..s+m+1]对应的整数值t(s+1)和t(s)相比,只是去掉了最高位数字T[s+1],增加了一个最低位的数字T[s+m+1]。因而t(s+1)和t(s)之间有如下关系式:
t(s+1)=d*(t(s)-T[s+1]*d^(m-1))+T[s+m+1]
例如对于文本十进制数字组成的文本”2359023141526739921”,d=10,m=5,当s=6时,t(s)=”31415”,T[s+1]=T[7]=3,T[s+m+1]=T[12]=2,所以
t(s+1) =10*(31415-3*10^4)+2=14152
因而我们初始时,只用求p=P[1..m],t0=T[1..m],t(s+1)的求值可以在每次循环迭代中去求,避免了多次函数调用的开销。因而一个可能的实现如下:
void Bad_Rabin_Karp_Matcher(const char* T, const char* P, int d) { int n = strlen(T); int m = strlen(P); int h = pow(static_cast<double>(d),m-1);//h=d^(m-1) int p=0; int t=0; for(int i=0; i<m; i++) { p = p*d + P[i]-'0';//p=P[0,m-1] t = t*d + T[i]-'0';//t=T[0,m-1] } for(int s=0; s<=n-m; s++) { if(p==t) { cout<<"Pattern occurs with shift"<<s<<endl; } if(s<n-m) { t = d*(t-h*(T[s]-'0')) + T[s+m]-'0'; } } }
注:
(1) 在前面讲解过程中数组下标从1开始,而在实际编程中数组下标是从0开始的,所以代码和描述的有点差异,但是原理是一样的。
(2)注意到我们求p,t,h的过程,不管用的是int还是long long型,当m增大时,肯定会产生溢出,比如说一个32位的int型表示的最大整数也就是20多亿。
(3)一般做过一些ACM题目的同学,知道在这种情况我们可以通过对一个数求模来防止溢出,记其为q。
(4)这里又出现了一个问题,当p≠t (mod q)时,的确可以推出p≠q;但是当p=t (mod q)时,推不出p=t。即可能出现所谓的伪命中点,即p=t (mod q)但是p≠q。所以当p=t (mod q)时,我们要进行额外检查,依次比较这m个字符,看其是否真的命中。
(5)通过把q设置为一个较大的素数,可以有效的减少伪命中点数,因而可以减小额外检查的开销。
通过以上分析,改进后一个可能的算法实现如下:
void Rabin_Karp_Matcher(const char* T, const char* P, int d, int q) { int n = strlen(T); int m = strlen(P); int p = 0; int t = 0; int h = d;//当m=1时有问题,此处应该为h=1,下面循环中k初始为1 for(int k=2; k<m; k++) { h = h*d % q; } for(int i=0; i<m; i++) { p = (p*d + P[i]-'0') % q; t = (t*d + T[i]-'0') % q; } for(int s=0; s<=n-m; s++) { // cout<<"p="<<p<<" t="<<t<<endl; if(p==t) { if(strncmp(T+s,P,m)==0) { cout<<"Pattern occurs with shift"<<s<<endl; } } if(s<n-m) { t = (d*(t-h*(T[s]-'0'))+T[s+m]-'0') % q;//此处也有问题,t可能小于0 } } }
但是运行结果和和预期不符,比如说你输入P=”111111”,T=”1”,结果只输出了一个位移值0,输入T=”1234”,P=”34”居然没有有效的位移输出。
分析代码,唯一可能的原因就是求余过程出现到了问题,因此我们在循环中比较p和t的值之前先输出它们的值,便于分析。果不其然,拿P=”111111”,T=”1”的输出作为示例,除了第一个位移p和t相等,后面的位移t居然是负值。也就是说对q求余的结果可能出现[-q+1,-1]范围的负值,为了计算的正确性,我们要保证迭代的过程中t都为正值。当t<0时,只需加上q即可。当然初始计算p和t的值时,p和t是不可能为负值的。修改完代码后,上面举的第二个输入用例运行正确,但是第一个仍然有问题,分析发现,原来当模式P只含一个字符时m=1,h=d^(m-1)=d^0=1,而前面我在求h时,只想到h是m-1个的相乘,把h初始化为了d,疏忽了m可能取1,所以出现了此种情况。到此分析结束,完整正确的函数及测试代码如下:💗
#include <cstdlib> #include <iostream> #include <cstring> #include <cmath> using namespace std; void Rabin_Karp_Matcher(const char* T, const char* P, int d, int q) { int n = strlen(T); int m = strlen(P); int p = 0; int t = 0; int h = 1; for(int k=1; k<m; k++) //预处理 { h = h*d % q; } for(int i=0; i<m; i++) { p = (p*d + P[i]-'0') % q; t = (t*d + T[i]-'0') % q; } for(int s=0; s<=n-m; s++)//匹配 { cout<<"p="<<p<<" t="<<t<<endl; if(p==t) { if(strncmp(T+s,P,m)==0) { cout<<"Pattern occurs with shift"<<s<<endl; } } if(s<n-m) { t = (d*(t-h*(T[s]-'0'))+T[s+m]-'0') % q; if(t<0) { t = t+q; } } } } int main(int argc, char *argv[]) { const int Max_Length = 1000; const int d = 10; const int q = 13; char T[Max_Length]; char P[Max_Length]; while(gets(T)) { gets(P); Rabin_Karp_Matcher(T,P,d,q); cout<<"next case:"<<endl; } system("PAUSE"); return EXIT_SUCCESS; }