FFT 处理字符串匹配问题

介绍一种新的处理字符串匹配问题的思路

前言

KMP 算法做字符串匹配相信大家都不陌生,然后最近 CF 上学到一个新技术:用 FFT 处理字符串匹配。我们都知道 FFT 是快速傅里叶变换算法,通常用于处理多项式乘法问题。那么 FFT 处理字符串又是怎么回事呢?下面就让小编来带大家看一看吧。

模式匹配问题

字模式匹配问题是几乎所有字符串问题的基础。该问题的目标是在主串中寻找所有子串出现的位置,即对字符串 ,寻找所有 满足 。首先回忆一下传统的几个算法,brute force、KMP 和 Boyer-Moore 算法。

朴素的暴力直接枚举每个 ,按位依次比较各个字符串是否相等。暴力算法的最坏复杂度是 。1977年三位算法家 Knuth、Morris 和 Pratt 提出的 KMP 算法,利用每次配对过程中成功配对的字符信息,通过 的预处理得到一个 fail 数组,表示子串的每个位置失配时跳转到哪里继续尝试,实现了最坏复杂度 的模式匹配操作。

虽然在广大高校的数据结构课中学习的都是 KMP 算法,但它绝不是模式匹配问题的最优解,甚至在随机的数据下 KMP 往往还没有暴力跑得快。同样在1977 年提出的 Boyer-Moore 算法(也称 Horspool 算法)则优秀许多,它引入了从右到左匹配的思想,通过分别保存每个字符和子串的每个后缀在子串中出现的最后位置,该算法实现了 的亚线性平均复杂度,在实际应用中效果出色。BM 算法的代价是空间占用较高,而且它在精心设计的特殊串上可以被卡出 的最坏复杂度,当原串不包含循环节时最坏复杂度则为

感兴趣的同学可以参考原论文《A Fast String Searching Algorithm》;相当感兴趣的同学推荐 Dan Gusfield 的一本书《Algorithms on Strings, Trees, and Sequences》,其中对各种传统字符串算法有详尽的介绍。

以上算法本质上都是对暴力匹配算法的改进,而本文介绍的思路有所不同,它起源于如下的构造:

首先定义一种映射,将 字符串转化为多项式。将字符串中的每个字母映射到一个数字,然后将这些数字看做多项式各项的系数,字符串从左往右的字母依次对应多项式中从低次到高次的项。

比如对于 ,如果将 a 转化为 1,b 转化为 2,然后将每个数字看做多项式的系数,则 就转化为了 ,不妨记做 则转化为多项式

接下来构造以下多项式:

为什么要算它呢?注意到,当且仅当原字符串 时,有多项式的对应系数 (因为对应字母各个相等,平方和才会等于 0 )。也就是说,计算出 W 后,看其中哪些系数等于 0,就可以 地找到原字符串中匹配的位置。

那么接下来只要计算 即可。对平方式进行展开,有:

其中的 是可以通过预处理前缀和 地求出的,而 是一个定值。剩下的关键问题就是求 ,这时候就想到了用 FFT 。

FFT 算法可以在对数复杂度内对两个多项式相乘,也就是给定长度为 n 的多项式 ,在 的复杂度内计算

把 (1) 式和 (2) 式凑一下,有

其中 的逆序字符串, 为两个多项式之积。这样就将 的计算在一次多项式乘法中完成了。

最后总结一下,使用 FFT 进行字符串匹配,首先将 S 和 T 转换为多项式,将 T 逆序并补到和 S 同样的长度(若 T 比 S 长则直接返回无配对),然后将两个多项式相乘;同时预处理一下 T 平方的前缀和。然后依次根据 (1) 式和 (3) 式计算每个位置 W 的系数,若为 0,就表示找到一个匹配位置。该算法的复杂度为

Bonus,由于我们只关心系数是否等于 0,因此用 NTT 代替 FFT 也是可以的,NTT 在没有板子的情况下比 FFT 好写一点,而且不用担心奇怪的精度问题,但有可能被构造一个 W 正好等于 998244353 的情况卡 WA。一种解决办法是对把每个字母对应的数值 random_shuffle 一下,或者用两个模数各算一遍。

参考代码如下:

#include "bits/stdc++.h"
using namespace std;
typedef complex<double> c;
void fft(vector<c>& f, int on); // fft实现略
vector<int> mul(vector<int>& a, vector<int>& b, int len){
    vector<c> pa(len+len), pb(len+len);
    for(int i=0; i<len; i++){
        pa[i] = {a[i], 0};
        pb[i] = {b[i], 0};
    }
    fft(pa, 1);
    fft(pb, 1);
    for(int i=0; i<pa.size(); i++){
        pa[i] *= pb[i];
    }
    fft(pa, -1);
    vector<int> ans(pa.size());
    for(int i=0; i<pa.size(); i++){
        ans[i] = (pa[i].real()/pa.size() + 0.5);
    }
    return ans;
}
int main(){
    string a, b;
    cin >> a >> b;
    int len = a.length() + b.length();
    vector<int> S(len), T(len), sumS(len+1);
    for(int i=0; i<a.length(); i++){
        S[i] = a[i] - 'a' + 1;
        sumS[i+1] = sumS[i] + S[i] * S[i];
    }
    for(int i=a.length(); i<len; i++){
        sumS[i+1] = sumS[i];
    }
    int sumT = 0;
    for(int i=0; i<b.length(); i++){
        T[i] = b[b.length() - 1 - i] - 'a' + 1;
        sumT += T[i] * T[i];
    }
    vector<int> C = mul(S, T, len);
    for(int i=0; i<a.length(); i++){
        int res = sumS[i+b.length()] - sumS[i] + sumT - 2 * C[i+b.length()-1];
        printf("%d", res==0);
    }
    return 0;
}

那么有人会问了,这个算法搞了半天,复杂度还没 KMP 优秀,到底有什么用呢?下面列举两个应用。

带通配符的模式匹配

假如出一个题: 分别是两个字符串,其中包含一些字母和一些”.“表示通配符。通配符可以匹配任意字母,通配符也可以匹配通配符。

比如说 ="ab.abc",="a.c" 那么 在位置 0 和位置 3 处匹配。

带通配符的情况下,传统的暴力算法依然生效,KMP 等算法就很难生效了。此时 FFT 方法就可以发挥作用,只要把多项式 W 改写为

其中通配符对应的数字为 0。此时多项式展开后化为

三项都使用多项式乘法预处理一下即可。

感兴趣的同学可以去做一下洛谷 - 残缺的字符串,注意用 FFT 的话得开 double 算,会爆 int;也可以 NTT。

参考代码

typedef complex<double> c;
void fft(vector<c>& f, int on);
vector<double> mul(vector<double>& a, vector<double>& b, int len){
	vector<c> pa(len+len), pb(len+len);
	for(int i=0; i<len; i++){
        pa[i] = {a[i], 0.};pb[i] = {b[i], 0.};
	}
	fft(pa, 1);
	fft(pb, 1);
	for(int i=0; i<pa.size(); i++){
        pa[i] *= pb[i];
    }
	fft(pa, -1);
	vector<double> ans(pa.size());
	for(int i=0; i<pa.size(); i++){
        ans[i] = pa[i].real()/pa.size();
	}
	return ans;
}
int main(){
	int n, m;
	cin >> m >> n;
	string a, b;
	cin >> b >> a;
	if(n<m) return puts("0"), 0;
	int len = n;
	vector<double> S(len, 0), T(len, 0), S2(len, 0), T2(len, 0), S3(len, 0), T3(len, 0);
	for(int i=0; i<n; i++)if(a[i]!='*'){
        S[i] = a[i] - 'a' + 1;
        S2[i] = S[i] * S[i];
        S3[i] = S[i] * S2[i];
	}
	for(int i=0; i<m; i++)if(b[m-1-i]!='*'){
        T[i] = b[m-1-i] - 'a' + 1;
        T2[i] = T[i] * T[i];
        T3[i] = T[i] * T2[i];
	}
	vector<double> C1 = mul(S3, T, len), C2 = mul(S, T3, len), C3 = mul(S2, T2, len);
	vector<int> ans;
	for(int i=0; i<=n-m; i++){
        int res = C1[i+m-1] + C2[i+m-1] - 2 * C3[i+m-1];
        if(abs(res) < 0.5)ans.push_back(i);
	}
	printf("%d\n", ans.size());
	for(auto i:ans){
        printf("%d ", i+1);
    }
	return 0;
}

Codeforces 1334G

学完这个技术顺便秒一道 education round 的 G 题。

这题数据范围比较水,用 bitset 压位可以做(建 26 个 bitset 分别代表每个字母是否出现,然后或运算统计贡献,复杂度 ),甚至优化的暴力也能莽过去。

本题的多项式为

暴力化简一波即可,最终化为一个 T 四次方的前缀和、一个常数加上三个 FFT。

本题的官方题解还给了另一种构造 FFT 的方式,只要一次 FFT 即可,挺巧妙的。但我没搞懂是怎么想到的,而且似乎也不具有通用性,所以就没研究了。