KMP算法

任何优秀的算法都是简约而美丽的。KMP更是如此。

下面这些定义是十分重要的，功欲善其事，必先利其器。

一、基本定义（Basic Definitions）

设AA为一个字符集，并且x=x0…xk−1x=x_0…x_{k-1}，kk为自然数，xx是长度为 kk 的在AA上的一个字符串。

xx 的前缀(prefix)为一个字串uu,其中：

u=x0x1…xb−1,b∈{0,1,…,k}
u = x_0x_1…x_{b-1}, b\in\{0,1, …,k\}

xx 的后缀(suffix)为一个字串uu,其中：

u=xk−bxk−b+1…xk−1,b∈{0,1,…,k}
u = x_{k-b}x_{k-b+1}…x_{k-1}, b\in\{0,1, …,k\}

如果u≠x(即b<k)u \neq x (即b < k), 则 xx 的一个前缀或者后缀 uu 被称为特征前缀(proper prefix)或特征后缀(proper suffix)。

设 xx 的一个边界框(border)(这个概念相当重要，往后读读就会明白)是一个子串 rr

其中：

r=x0x1…xb−1且r=xk−bxk−b+1…xk−1,b∈{0,1,…,k−1},k=|x|
r = x_0x_1…x_{b-1} \color{red}且 r=x_{k-b}x_{k-b+1}…x_{k-1}, b∈\{0,1, …,k-1\},k=|x|

边界框是 xx 的一个真子串，这个真子串既是特征前缀也是特征后缀，且特征前缀和特征后缀相等。我们称长度bb为边界框的宽度(width)（说是厚度也许更确切，更能表达意思），如果某个边的宽度是所有边中最大的，则称该边为最宽边界框(widest border)。

可以看到边界框具有一个很好的性质：即左边界=右边界。

例：设 x=abacabx = abacab.则 xx 的特征前缀分别为 ε,a,ab,aba,abac,abaca,\varepsilon, a, ab, aba, abac, abaca, 特征后缀分别为 ε,b,ab,cab,acab,bacab，x\varepsilon, b, ab, cab, acab, bacab，x 具有两个边界框 ε和ab\varepsilon 和 ab。

其中边界框 ε\varepsilon 表示宽度为0的串，边界框 abab 的宽度为2，最宽边界框很明显为abab.

把这个最宽边界框加红一下：x=abacabx =\color {red}{ab}ac\color {red} {ab}，可以看到，这个红色部分非常像一个边框，这也是 border 的含义。

对于任意字符串 x∈Ax\in A（AA是一个字符集），其中空串 ε\varepsilon 总为 xx 的一个边界框。空串ε\varepsilon本身没有边界框。

KMP算法中的位移量(shift distance)将会用到字符串中边界框的概念。

例：

i	0	1	2	3	4	5	6	7	8	9
t	aa	bb	cc	aa	bb	cc	aa	bb	dd
p	a\color{blue}a	b\color{blue}b	c\color{green}c	a\color{blue}a	b\color{blue}b	d\color{red}d
p				a\color{blue}a	b\color{blue}b	cc	a\color{blue}a	b\color{blue}b	dd

位置 0,…,40,…,4 的字符已经完全匹配，但是在位置5，cc和dd不匹配。于是模式串向右移动三个位置，接着在位置5继续进行比较。

其中位移量取决于已匹配的字符串（如上图中 abcababcab 部分）的最宽边界框(上图的蓝色部分)，在这个例子中，已经匹配的 abcababcab 的串长 j=5j=5 ，其中最宽边界框 abab 的宽度为 b=2b=2 ，于是位移量为 j−b=5−2=3j-b = 5-2=3.

这样做的原因是边界框有左右两个部分，如上面的 abcab\color{red}{ab}c\color{red}{ab}，阅读前面的部分可以知道左右两个部分是完全一样的。如果p0p1...pj−1p_0p_1...p_{j-1} 部分已经完全匹配，但是在 jj 处失配(pj≠tip_j \neq t_i)，只需要把左边界的部分移动到原先右边界的位置，然后接着看是否可以匹配。

在预处理阶段，应该先求出模式串的每个前缀的最宽域的宽度。然后在搜索阶段，位移量可以根据已经匹配上的前缀计算得出。

二、预处理(Preprocessing)

定理（Theorem）：设串 r,sr, s 都是串 xx 的边界框，其中 |r|<|s||r| < |s|，则串 rr 是串 ss 的边界框。

证明：图1中的串 xx 包含了两个边界框 rr 和 ss.由于 rr 是 xx 的前缀，同时它也是 ss 的特征前缀，而且 rr也是 xx 的后缀，rr也是 ss 的特征后缀，又因为 rr 的长度小于 ss，因此 rr 是 ss 的边界框。



图1 边界框中的边界框

例： 串abacabaabacaba有两个非空边界框，一个是 aa ,标红后为 abacaba\color{red}abacab\color{red}a，还有一个是 abaaba，标记成绿 abacaba\color{green}{aba}c\color{green}{aba}，很明显可以看到 aa 也是 aba\color{red}ab\color{red}a 的边界框。

这个性质非常有用，后面可以看到。

如果串 ss 是 xx 的最宽边界框，则 xx 的下一最宽边界框(next-widest border)为 ss 的最宽边界框 rr.

说的直白点，下一个最宽边界框 = (xx的最宽边界框)的最宽边界框

定义：设 xx 是一个字符串，且 aa 是 AA 上的一个字符，如果 rara 是 xaxa 的一个边界框，则称 xx 的边界框 rr 可通过字符 aa 延拓(extend)。



图2 边界框延拓

图2中，如果 xj=ax_j=a，则 xx 的宽为 jj 的边界框可通过字符 aa 延拓为 rara，rara是 xaxa的边界框。

例：s=cbobcbos=cbobcbo, ss的一个子串 x=cbobcbx=cbobcb 可以通过字符 oo 延拓为 ss, 已知 xx 的边界框是 cbcb，边界延拓后 cbocbo 是 ss 的边界框。

在预处理阶段，构造一个长度为 m+1m+1 的数组 bb，数组的每一项 b[i]b[i] 为模式串 p=p0p1⋯pm−1p=p_0p_1\cdots p_{m-1} 的长度为 ii 的前缀(p0p1⋯pi−1)(p_0p_1\cdots p_{i-1})的最宽边界框的宽度 (i=0,…,m)(i=0,…,m).由于长度为 0 的前缀 ε\varepsilon 没有边界框，我们规定b[0]=−1b[0] = -1。



图3 模式串的长度为i的前缀的宽度为b[i]的域

假设 b[0]…b[i]b[0]…b[i] 的值已知，则 b[i+1]b[i+1] 的值可以通过检测串 p0…pi−1p_0…p_{i-1} 是否可以通过字符 pip_i 延拓来计算。在图3中，就是判断是否有 pi=pb[i]p_i = p_{b[i]} ,注意灰色部分是 p0…pi−1p_0…p_{i-1} 的边界框，其宽度就是 b[i]b[i]。如果 pi=pb[i]p_i = p_{b[i]} ，可以得到 b[i+1]=b[i]+1b[i+1]=b[i]+1. 利用 b[i],b[b[i]],b[b[b[i]]]…b[i],b[b[i]],b[b[b[i]]]… 的降序可以获得一个可供检测的边界框集(borders)。注意到这个降序序列，读成一句话就是“边界框的宽度，边界框的边界框的宽度，边界框的边界框的宽度……”

预处理算法包括了一个含变量 jj 的循环，用来遍历这些值，即b[i],b[b[i]],b[b[b[i]]]…b[i],b[b[i]],b[b[b[i]]]…

如果 pj=pip_j=p_i，则宽为 jj 的边界框可以通过字符 pip_i延拓；否则，通过把 jj 设为 b[j]b[j]，即去查找下一最宽边界框，看是否可以延拓，可以那就可以去设定 b[i+1]b[i+1] 的值了，如果不可以，继续查找，直到再也找不到下一最宽边界框为止，即 j=−1j = -1的时候。

每次出现 jj++ 后，jj 的值就是p0…pip_0…p_i的最宽边界框的宽度，因为找到了一个字符 pj=pip_j=p_i 可以把边界框 p0p1…pj−1p_0p_1…p_{j-1} 延拓为新的边界框 p0p1…pjp_0p_1…p_{j}，因此 p0…pip_0…p_i 的最宽边界框的宽度就是 p0p1…pj−1p_0p_1…p_{j-1} 这个边界框的宽度加1。然后把 b[i+1]b[i+1] 的值设置为 jj (也就是 ii++ 后设置 b[i]b[i] 的值)。下面是预处理处算法代码：

void kmpPreprocess()
{
//i: 当前指针。j: 当前边界框的宽度。
int i = 0, j = -1;
b[i] = j; //初始化b[0]
while (i < m)
{
//查找下一最宽边界框，直到可以延拓
while (j >= 0 && p[i] != p[j]) j = b[j];//边界框宽度大于等于0且无法延拓
// 进行延拓。
i++; j++;
b[i] = j;
}
}

例：对于模式串 p=ababaap=ababaa，数组 bb 中保存的边界框集宽度分别如下。例如，长度为55的前缀 ababaababa 有一个宽度为3的域，因此b[5]=3b[5]=3.

j	0	1	2	3	4	5	6
p[j]	a	b	a	b	a	a
b[j]	-1	0	0	1	2	3	1

三、搜索算法(Searching algorithm)

假设我们把上面算法的模式串 pp（长度为mm）改成 ptpt (ptpt 是模式串 pp 和目标串 tt 的连接，见图4)，上面的预处理算法理所当然也适用于计算串 ptpt 的边界框的宽度，如果把串 pp 看成 ptpt 的某个长度为 ii 的前缀 x(=pt0pt1⋯pti−1)x(=pt_0pt_1\cdots pt_{i-1}) 的一个边界框，只要找到了这样的前缀 xx，那就等于找到了匹配串的位置，即 i−mi-m。

串 ptpt 的某个前缀 x(=pt0pt1⋯pti−1)x(=pt_0pt_1\cdots pt_{i-1}) 正好有一个宽度为 mm 的域，那就说明搜索成功了，这个时候匹配位置为 i−mi-m，接着继续匹配下一个位置（如图4）。



图4 ptpt 的一个前缀 xx 的宽度为 mm 的边界框

搜索算法如下：

void kmpSearch()
{
//i: 当前指针 j: 当前边界框的宽度
int i = 0, j = 0;
while (i < n)
{
// 当前位置无法延拓就继续搜索下一最宽边界框，直到可以延拓。
while (j >= 0 && t[i] != p[j]) j = b[j];
i++; j++;
//如果边界框的宽度正好为m，说明匹配到模式串了
if (j == m)
{
//上报匹配结果
report(i – j);

//将j降到下一最宽边界框，然后继续执行下一次匹配。实际上，这行代码也可略去不写。
//因为一下次执行到while(j>=0&&t[i]!=p[j])j=b[j];也会因为无法延拓而执行循环体。
//此时的p[j]肯定等于字符'\0'
j = b[j];
}
}
}

当内循环在 jj 处无法延拓时，则检查下一最宽边界框是否可以延拓，即检查模式串的长度为 jj 的最宽边界框是否可以延拓（如图5）。如果仍无法延拓，则继续检查下一最宽边界框，直至下一最宽边界框为空（即j=−1j = -1时），或者可以延拓为止。



图5 在 jj处失配（无法延拓）后模式串的移动，即看是否可以延拓下一最宽边界框

如果所有的mm个字符都可以匹配上，这时候j=mj=m，函数reportreport的作用就是上报匹配位置。接下来，把最宽边界框降序，继续匹配下一个地方。直到完成外层整个循环。

下面给出一个完整的代码：

#include "string.h"

char t[] = "ababbababacabacababacacbacababacababaa";
char p[] = "ababac";
const int n = strlen(t);
const int m = strlen(p);
int b[7] = {0};

void KmpPreprocess()
{
int i = 0, j = -1;
b[i] = j;
while (i < m)
{
while (j >= 0 && p[i] != p[j]) j = b[j];
i++; j++;
b[i] = j;
}
}
void report(int nIndex)
{
printf("%d ", nIndex);
}
void KmpSearch()
{
int i = 0, j = 0;
while (i < n)
{
while (j >= 0 && t[i] != p[j]) j = b[j];
i++; j++;
if (j == m)
{
report(i - j);
j = b[j];
}
}
}

int _tmain(int argc, _TCHAR* argv[])
{
KmpPreprocess();
KmpSearch();
getchar();
return 0;
}

参考文献：http://www.inf.fh-flensburg.de/lang/algorithmen/pattern/kmpen.htm