扩展KMP算法 Extend KMP

From http://noalgo.info/340.html 

前面介绍过经典的KMP算法，本文继续介绍KMP算法的扩展，即扩展KMP算法：给定一个较长的主字符串S和一个较短的模式串T，用m、n分别表示S和T的长度（即m=|S|，n=|T|），定义extend[i]=S[i..n]与T的最长公共前缀的长度，求extend[i]。
为什么说这是KMP算法的扩展呢？显然，如果在S的某个位置x有extend[x]==n，则可知在S中可以找到T，并且首位置是x，这正是KMP算法要解决的问题。而且，扩展KMP算法能找到S中所有T的匹配，更一般地，可以知道S中以每个字符开始的后缀与T的最大的匹配长度（最长公共前缀长度）。

以下所有例子及代码实现中，字符串下标均从0开始！

一 基本思想

扩展KMP算法可以在线性时间内求出extend数组的值。

令S=aaaaabbb，T=aaaaac，如下

S:aaaaabbbT:aaaaac  

首先，通过6次比较可以知道,extend[0]=5。然后要计算extend[1]，如下

S:aaaaabbbT: aaaaac 

可以看到，extend[1]=4。但是这里是否要经过5次比较才可以得到4这个结果呢？其实，正如KMP算法一样，这里也可以充分利用在前面匹配中已得到的信息来简化当前的匹配。

事实上，由extend[0]=5我们可以知道，S[0..4]=T[0..4]

因为现在extend[1]考虑的是S中第1个字符开始的后缀，则取出S的下标从1开始的部分，有S[1..4]=T[1..4]

计算extend[1]，实际上用S中S[1]开始的后缀去匹配T，根据上式可以知道，S中前面4个字符的匹配就是T中T[1]开始的后缀的匹配，即用T来跟T自己匹配！

假设我们初始化了一个数组next，next[i]表示T[i..n]跟T的最长公共前缀的长度。对于这里，有next[1]=4，即T[1..4]=T[0..3]

综合以上的两个等式，有S[1..4]=T[1..4]=T[0..3]

即S的第1到4个字符跟T的前4个字符是匹配的，我们可以直接从S的第5个字符开始比较。

S:aaaaabbbT: aaaaac 

这里发现，S[5]=b，T[4]=a，于是S[5]!=T[4]，于是知道extend[1]=4。

二 具体算法之extend数组

以上只是算法一个简单直观的描述，下面介绍算法的具体实现方法。

扩展KMP算法顺序遍历S中的每个字符，遍历过程中依次计算每个位置的extend值。

假设当前遍历到位置i，即extend[0..i-1]这i个位置的值已经计算得到。算法在比遍历过程中记录了匹配成功的字符的最远位置p及这次匹配的起始位置a，事实上，匹配位置x时匹配成功的最远位置是x+extend[x]-1，p就是x=0…i-1时得到的最大值，a就是取到最大值时对应的那个x。

S串： 0 1 2 … a … i-1 i … p …

对于a这个位置，我们有 S[a..p] = T[0..p-a]

取i开始的那一段，则有 S[i..p] = T[i-a..p-a]

要求extend[i]，即S[i...]匹配T[0...]，我们可以先用T[i-a...]来匹配T[0...]，这正是next[i-a]表示的含义：T中第位置i-a开始的后缀跟T的最长公共前缀长度。令L=next[i-a]，分以下两种情况讨论：

1. i + L < p，如下所示
   S:012…a…i-1ii+1…i+L……P…T:       01…LL+1………

   S跟T在绿色部分一定匹配，且根据next数组的含义，红色部分一定不匹配，因为如红色部分匹配的话，S[i..p]可以看成是T[i-a..p-a]，则T[i-a...]与T[0...]的匹配会大于L，与最长公共前缀矛盾。
   于是，无需任何比较就可以得到，extend[i]=L，同时a和p保持不变。

2. i + L >= p，如下所示
   S:012…ai-1i…pp+1………T:      0…p-ip-i+1…L…

   这时候绿色部分是已经匹配成功的，但因为之前匹配到最远的位置只有p，所以红色部分p+1及其以后位置的匹配结果是未知的。所以接下来要从S[p+1]和T[p-i+1]开始逐个向右匹配，直到匹配失败为止。匹配完之后，可以得到extend[i]的值和p的值，同时需要更新a的值为当前的i。

这就是扩展KMP算法的主要思想，下面是其C++实现的代码。注意程序中为了实现的方便，p表示的是匹配成功的字符的最远位置的下一个位置，即下次匹配直接从S[p]开始（上面的是从S[p+1]开始）。同时使用了变量j跟踪T中跟S中p位置对应的位置，注意遍历S时i在递增的同时，j需要递减。

可以看到，扩展KMP算法主要依赖指针p的值计算extend数组，而p只会扫描主串S一遍，于是算法的复杂度是线性的，O(m+n)。

//扩展KMP算法求extend数组，s为文本串，t为模式串void getExtend (const char *s, const char *t, int *extend){int ls = strlen(s), lt = strlen(t);int next[mx]； getExtendNext(t, next);//计算next数组。mx是表示最大长度的常数。for (int i = 0, j = -1, a, p; i < ls; i++, j--)if (j < 0 || i + next[i - a] >= p){if (j < 0)j = 0, p = i;while (p < ls && j < lt && s[p] == t[j])j++, p++;extend[i] = j, a = i;}elseextend[i] = next[i - a];}

三 具体算法之next数组

但还有一个问题没有解决，next数组怎么得到呢？

根据定义，next[i]表示T[i..n]和T[0..n]的最长公共前缀的长度，相当于以上算法中把S换成T的版本，即T自己和自己的一个匹配。不过这里next数组是从1开始的，next[0]匹配的肯定是整个串的长度n。

具体C++代码如下：

//计算next数组，保存到next参数中。void getExtendNext (const char *t, int *next){int lt = strlen(t);for (int i = 1, j = -1, a, p; i < lt; i++, j--)if (j < 0 || i + next[i - a] >= p){if (j < 0)j = 0, p = i;while (p < lt && t[j] == t[p])j++, p++;next[i] = j, a = i;}elsenext[i] = next[i - a];}

可以看到，KMP算法和扩展KMP算法再具体求解next数组时都是利用自身和自身匹配的一个思路，于是求next数组的代码跟主体匹配的代码是非常相似的，在编写代码时，可以利用这个相似性快速编码，并减少出错概率。

KMP算法和扩展KMP算法都是充分利用前面匹配的结果，在匹配失败时大幅度移动指针，直接跳过不可能的状态，减少算法的复杂度。

My Summary:

s[a,p-1] = t[0, p-a-1]   -> s[i,p-1]=t[i-a, p-a-1]

t[i-a,next[i-a]+i-a-1] = t[0, next[i-a]-1]

So if next[i-a]+i-a-1 >= p-a-1, extend j,p, extend[i]=j

else extend[i]=next[i-a]

My Code:

  void getNextArr(const char *t, vector<int>& nextarr) {    int lt = strlen(t);    nextarr.resize(lt,lt);    for (int i = 1, j = -1, a = 0, p = 0; i < lt; ++i, --j) {      if (j < 0 || i + nextarr[i - a] >= p) {        if (j < 0)  j = 0, p = i;        while (p < lt && t[p] == t[j]) ++p, ++j;        nextarr[i] = j, a = i;      }      else {        nextarr[i] = nextarr[i - a];      }    }  }  void getExtendArr(const char *s, const char *t, vector<int>& extend) {    int slen = strlen(s), tlen = strlen(t);    extend.resize(tlen);    vector<int> nextarr;    getNextArr(t, nextarr);    for (int i = 0, j = -1, a = 0, p = 0; i < tlen; ++i, --j) {      if (j < 0 || i + nextarr[i - a] >= p) {        if (j < 0) j = 0, p = i;        while (p < slen && j < tlen && s[p] == t[j]) ++p, ++j;        extend[i] = j, a = i;      }      else {        extend[i] = nextarr[i - a];      }    }  }

Debug Code:

#include <stdio.h>#include <iostream>#include <vector>using namespace std;class Solution {public:  void getExtendNext(const char *t, int *next) {    int lt = strlen(t);    for (int i = 1, j = -1, a, p; i < lt; i++, j--) {      if (j < 0 || i + next[i - a] >= p) {        if (j < 0)j = 0, p = i;        while (p < lt && t[j] == t[p])j++, p++;        next[i] = j, a = i;      }      else {        next[i] = next[i - a];      }    }  }  void getExtend(const char *s, const char *t, int *extend) {    int ls = strlen(s), lt = strlen(t);    int next[100]; getExtendNext(t, next);    for (int i = 0, j = -1, a, p; i < ls; i++, j--)    if (j < 0 || i + next[i - a] >= p) {      if (j < 0)j = 0, p = i;      while (p < ls && j < lt && s[p] == t[j])j++, p++;      extend[i] = j, a = i;    }    elseextend[i] = next[i - a];  }  void getNextArr(const char *t, vector<int>& nextarr) {    int lt = strlen(t);    nextarr.resize(lt,lt);    for (int i = 1, j = -1, a = 0, p = 0; i < lt; ++i, --j) {      if (j < 0 || i + nextarr[i - a] >= p) {        if (j < 0)  j = 0, p = i;        while (p < lt && t[p] == t[j]) ++p, ++j;        nextarr[i] = j, a = i;      }      else {        nextarr[i] = nextarr[i - a];      }    }  }  void getExtendArr(const char *s, const char *t, vector<int>& extend) {    int slen = strlen(s), tlen = strlen(t);    extend.resize(tlen);    vector<int> nextarr;    getNextArr(t, nextarr);    for (int i = 0, j = -1, a = 0, p = 0; i < tlen; ++i, --j) {      if (j < 0 || i + nextarr[i - a] >= p) {        if (j < 0) j = 0, p = i;        while (p < slen && j < tlen && s[p] == t[j]) ++p, ++j;        extend[i] = j, a = i;      }      else {        extend[i] = nextarr[i - a];      }    }  }  char *strStr(char *haystack, char *needle) {    // Note: The Solution object is instantiated only once and is reused by each test case.    if (needle == NULL || haystack == NULL)      return NULL;    int i = 0, j = -1, needlelen = strlen(needle), haystacklen = strlen(haystack);    if (needlelen > haystacklen)      return NULL;    if (needlelen == 0)      return haystack;    vector<int> next(needlelen + 1, -1);    while (needle[i]) {      if (j == -1 || needle[i] == needle[j]) {        ++i;        ++j;        next[i] = j;      }      else        j = next[j];    }    i = 0;    j = 0;    while (haystack[i]) {      if (j == -1 || haystack[i] == needle[j]) {        ++i;        ++j;        if (needle[j] == '\0')          return haystack + i - j;      }      else        j = next[j];    }    return NULL;  }};int main() {  Solution solution;  const char *ch = "aaaaabbb";  const char *s = "aaaaabbb";  const char *t = "aaaaac";  int nextarr[8];  vector<int> next, extend;  //solution.getNextArr(ch, next);  //solution.getExtendNext(ch, nextarr);  solution.getExtend(s, t, nextarr);  solution.getExtendArr(s, t, extend);    return 0;}