扩展KMP与Manacher

KMP已经介绍过了，这次主要介绍的是Manacher

最长回文子串

Manacher算法要解决的问题就是求最长回文子串，用到的思维和扩展KMP实在是像，不过理解起来比扩展KMP简单。

先定义数组v，v[i]表示以第i个字符为中心，到回文串一端的距离，我们以字符串"acabacab"为例，如下表（index是下标）

i是当前要计算的指针，m是上次计算的指针，e是下一个要比较的位置的指针

然后++i，注意这时候，由于以字符b两边是对称的，所以在求v[4]的值的时候，可以先查v[2]的值，是1，所以v[4]也是1

再继续++i，同样的，在求v[5]的值的时候，可以先查v[1]的值，是2，但是，这个长度达到了e指针的位置，即i+2>=e，这时候就更新指针m，并扩展e的位置，即比较str[7]与str[3]，找到以下标5为中心的回文串边界。

然后，v[5]的值就是e-i，再接着，求v[6]的值就查v[4]的

再接着，求v[7]的值就查v[3]的，不过v[3]的值是4，而i+4>=e又满足了，就再次更新指针m，并扩展e

最后得到的v数组最大值乘以2减1就是答案，以上指针在算法里只增不减，所以时间复杂度是 $O(n)$ 。但问题是，这只能求出奇数长度的回文串，所以，这里用到一个技巧，把原来的字符串"acabacab"转成"a#c#a#b#a#c#a#b"，其中#号是任意一个在原字符串中不会出现的字符，这样一来，任意原串的回文串都存在中心字符了。

还有一个细节，如果要计算的串是"a"，那么回文串的长度应该是1，但是，在这个算法里，会产生越界访问，因为要判断字符a后面的'\0'是不是与a前面的字符相同，如果正好都是'\0'，就会导致计算错误。所以，计算的时候还要对原来的字符串在最前面添加第二个在原字符串中不会出现的字符，例如转成"$a#c#a#b#a#c#a#b"

扩展KMP模板和Manacher模板对比

代码高度相似，你看看就知道了

void extkmp_z(const char* str, int z[])
{
    int s_len = strlen(str);
    z[0] = s_len;
    for (int i = 1, p = 1, e = 1; i < s_len; ++i)
    {
        if (i + z[i - p] >= e)
        {
            e = std::max(i, e);
            p = i;

            while (e < s_len && str[e] == str[e - i])
                ++e;

            z[i] = e - i;
        }
        else
            z[i] = z[i - p];
    }
}

void extkmp_ext(const char* str, int ext[], const char* pattern, int z[])
{
    int s_len = strlen(str);
    extkmp_z(pattern, z);
    for (int i = 0, p = 0, e = 0; i < s_len; ++i)
    {
        if (i + z[i - p] >= e)
        {
            e = std::max(i, e);
            p = i;

            while (e < s_len && str[e] == pattern[e - i])
                ++e;

            ext[i] = e - i;
        }
        else
            ext[i] = z[i - p];
    }
}

void init_str(const char* str, char* t)
{
    int len = strlen(str);
    t[0] = '$';
    t[1] = '#';
    for (int i = 0; str[i]; ++i)
    {
        t[(i + 1) << 1] = str[i];
        t[((i + 1) << 1) + 1] = '#';
    }
    t[(len + 1) << 1] = 0;
}

int manacher(const char* s, int v[])
{
    int len = strlen(s);
    int max_len = 0;
    v[0] = 1;
    v[1] = 1;
    for (int i = 1, m = 1, e = 1; i < len; ++i)
    {
        if (i + v[m - (i - m)] >= e)
        {
            e = std::max(i, e);
            m = i;

            while (e < len && s[e] == s[i - (e - i)])
                ++e;

            v[i] = e - i;
        }
        else
            v[i] = v[m - (i - m)];
        max_len = std::max(max_len, v[i]);
    }
    return max_len - 1;
}

以上是第一种写法的模板，再对比写法二

void extkmp_ext(const char* str, int ext[], const char* pattern, int z[])
{
    int s_len = strlen(str);
    extkmp_z(pattern, z);
    for (int i = 0, l = 0, r = 0; i < s_len; ++i)
    {
        if (i < r)
        {
            ext[i] = std::min(r - i, z[i - l]);
        }
        else
        {
            ext[i] = 0;
        }
        while (i + ext[i] < s_len && str[i + ext[i]] == pattern[ext[i]])
        {
            ++ext[i];
        }
        if (i + ext[i] > r)
        {
            r = i + ext[i];
            l = i;
        }
    }
}

int manacher(const char* s, int v[])
{
    int len = strlen(s);
    int max_len = 0;
    for (int i = 1, m = 1, r = 1; i < len; ++i)
    {
        if (i < r)
        {
            v[i] = std::min(v[m - (i - m)], r - i);
        }
        else
        {
            v[i] = 1;
        }
        while (s[i + v[i]] == s[i - v[i]])
        {
            ++v[i];
        }
        if (i + v[i] > r)
        {
            r = i + v[i];
            m = i;
        }
        max_len = std::max(max_len, v[i]);
    }
    return max_len - 1;
}

其它解法

用后缀数组来解回文串也是可以的，原串是s的话，把s反向，得到s1，然后构造S = s + '#' + s1，对S求后缀数组后，要注意的是并不是只找height的最大值即可，网上很多文章在这里的算法是错误的，反例是abcdefba，单纯枚举相邻的sa是不正确的。

正确做法是枚举i从1到len(s)/2，分别针对奇数长度和偶数长度要分别计算LCP（最长公共前缀），奇数长度时它的对称位置就是j=len-i+1，偶数长度时对称位置是j=len-i或j=len-i+2，然后对height数组求从rank[i]+1到rank[j]的最小值。最小值维护用ST表，由于ST表的构造是 $O(nlogn)$ ，而查询是 $O(1)$ ，所以生成后缀数组后的计算的总时间复杂度是 $O(nlogn)$

示例代码URAL 1297

int st[12][2048];

void init_st(int a[], int n)
{
    for (int i = 0; i < n; i++)
        st[0][i] = a[i];
    for (int k = 1; (1 << k) < n; k++)
        for (int i = 0, m = 1 << (k - 1); i < n; i++)
            if (i + m < n) st[k][i] = min(st[k - 1][i], st[k - 1][i + m]);
}

int st_query(int l, int r) // [l, r)
{
    if (l > r)swap(l, r);
    int k = (int)(log(r - l - 0.1) / log(2));
    return min(st[k][l], st[k][r - (1 << k)]);
}

int main()
{
    int n, t;
    char s[2100], s1[1010];
    scanf("%s", s1);
    strcpy(s, s1);
    strcat(s, "#");
    strcat(s, strrev(s1));
    SA_2_sort sa;
    sa.init(s);
    int len = sa.size();
    int maxh = 0, spos = 0;
    init_st(&*sa.ht.begin(), len + 1);
    for(int i = 1; i <= len / 2; ++i)
    {
        int j = len - i + 1;
        int lcp = st_query(sa.rk[i] + 1, sa.rk[j] + 1);
        if (lcp * 2 - 1 > maxh)
        {
            spos = i - lcp;
            maxh = lcp * 2 - 1;
        }
        if (i > 1)
        {
            j = len - i + 2;
            lcp = st_query(sa.rk[i] + 1, sa.rk[j] + 1);
            if (lcp * 2 > maxh)
            {
                spos = i - lcp - 1;
                maxh = lcp * 2;
            }
        }
    }
    strncpy(s1, s + spos, maxh);
    s1[maxh] = 0;
    puts(s1);
    return 0;
}

注：后缀数组模板请从之前的文章里复制；函数st_query是左闭右开区间，即求的是[l, r-1]的最小值，采用半开半闭区间目的是为了能直接swap