多项式套路入门 (1) - 单模式串匹配

普通单模式串匹配

可以直接 \(\text{KMP}\),但是考虑为下文作铺垫。

定义字符串匹配函数 \(C(x,y)=A(x)-B(y)\),表示 \(A\) 串的第 \(x\) 个字符和 \(B\) 串的第 \(y\) 个字符匹配。然后定义完全匹配函数 \(P(x)=\sum\limits_{i=0}^{m-1}[C(i,x-m+i+1)]^2\) ,若 \(P(x)\) 表示 \(0\),则 \(B\) 串以第 \(x\) 位结束的连续 \(m\) 位与 \(A\) 串完全匹配。翻转 \(A\) 串后记为 \(A'\) 串,则可以暴力拆开 \(P(x)\),易得到:

\[\qquad P_{x}=\sum\limits_{i=0}^{m-1}A'(m-i-1)+\sum\limits_{i=0}^{m-1}B(x-m+1+i)-2\sum\limits_{i=0}^{m-1}A'(m-i-1)B(x-m+1+i) \qquad \]

发现前两项都是暴力可做,第三项显然是个卷积,直接 \(\text{FFT}\) 即可做。

带通配符的单模式串匹配

\(A\) 串和 \(B\) 串中可能带通配符的问题。

考虑存在通配符可匹配,重新定义匹配函数 \(C(x,y)=[A(x)-B(y)]^2A(x)B(y)\)。则 \(P(x)=\sum\limits_{i=0}^{m-1}[C(i,x-m+i+1)]^2A(x)B(x-m+1+i)\)

暴力翻转和展开一下完全匹配函数,得到:

\[\qquad P(x)=\sum\limits_{i=0}^{m-1}A'(m-i-1)^3B(x-m+1+i)+\sum\limits_{i=0}^{m-1}B(x-m+1+i)^3A'(m-i-1)-2\sum\limits_{i=0}^{m-1}A'(m-i-1)^2B(x-m+1+i)^2 \qquad \]

显然三项皆可卷积,直接暴力卷积做即可。\(6\)\(\text{DFT}\)\(1\)\(\text{IDFT}\)

//残缺的字符串,模板题。
#include <cstring>
#include <cmath>
#include <iostream>
#include <algorithm>
#include <queue>
#include <set>
#include <vector>
#include <stack>
#include <map>
#include <bitset>
#define ri register
#define inf 0x7fffffff
#define E (1)
#define mk make_pair
#define int long long
//#define double long double
using namespace std; const int N=1200010; const double pi=acos(-1.0);
inline int read()
{
    int s=0, w=1; ri char ch=getchar();
    while(ch<'0'||ch>'9') {if(ch=='-') w=-1; ch=getchar(); }
    while(ch>='0'&&ch<='9') s=(s<<3)+(s<<1)+(ch-'0'), ch=getchar();
    return s*w;
}
void print(int x) { if(x<0) x=-x, putchar('-'); if(x>9) print(x/10); putchar(x%10+'0'); }
int m,n,ss[N],tt[N],rev[N],T,cnt,pos[N]; char s[N],t[N];
struct Complex
{
    double x,y;
    inline Complex ( double pp=0, double qq=0 ) { x=pp, y=qq; }
}a[N],b[N],f[N];
inline Complex operator + (Complex a,Complex b) { return Complex(a.x+b.x,a.y+b.y); }
inline Complex operator - (Complex a,Complex b) { return Complex(a.x-b.x,a.y-b.y); }
inline Complex operator * (Complex a,Complex b) { return Complex(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x); }
inline void Get_Rev() { for(ri int i=0;i<T;i++) rev[i]=(rev[i>>1]>>1)|((i&1ll)?(T>>1):0); }
inline void FFT(Complex *S,int type)
{
    for(ri int i=0;i<T;i++) if(i<rev[i]) swap(S[i],S[rev[i]]);
    for(ri int mid=1;mid<T;mid<<=1)
    {
        Complex wn(cos(pi/mid),type*sin(pi/mid));
        for(ri int r=mid<<1, j=0;j<T;j+=r)
        {
            Complex w(1,0);
            for(ri int k=0;k<mid;k++, w=w*wn)
            {
                Complex x=S[j+k], y=w*S[j+mid+k];
                S[j+k]=x+y, S[j+mid+k]=x-y;
            }
        }
    }
}
signed main()
{
    m=read(), n=read();
    scanf("%s%s",s,t);
    for(ri int i=m-1;~i;i--) ss[i]=(s[m-i-1]=='*')?(0):(s[m-i-1]-'a'+1);
    for(ri int i=0;i<n;i++) tt[i]=(t[i]=='*')?(0):(t[i]-'a'+1);
    T=1; while(T<=n) T<<=1; Get_Rev();
    for(ri int i=0;i<m;i++) a[i]=Complex(ss[i]*ss[i]*ss[i],0);
    for(ri int i=0;i<n;i++) b[i]=Complex(tt[i],0);
    for(ri int i=m;i<T;i++) a[i]=Complex(0,0);
    for(ri int i=n;i<T;i++) b[i]=Complex(0,0);
    FFT(a,1), FFT(b,1);
    for(ri int i=0;i<T;i++) f[i]=f[i]+a[i]*b[i];
    for(ri int i=0;i<m;i++) a[i]=Complex(ss[i],0);
    for(ri int i=0;i<n;i++) b[i]=Complex(tt[i]*tt[i]*tt[i],0);
    for(ri int i=m;i<T;i++) a[i]=Complex(0,0);
    for(ri int i=n;i<T;i++) b[i]=Complex(0,0);
    FFT(a,1), FFT(b,1);
    for(ri int i=0;i<T;i++) f[i]=f[i]+a[i]*b[i];
    for(ri int i=0;i<m;i++) a[i]=Complex(ss[i]*ss[i],0);
    for(ri int i=0;i<n;i++) b[i]=Complex(tt[i]*tt[i],0);
    for(ri int i=m;i<T;i++) a[i]=Complex(0,0);
    for(ri int i=n;i<T;i++) b[i]=Complex(0,0);
    FFT(a,1), FFT(b,1);
    for(ri int i=0;i<T;i++) f[i]=f[i]-Complex(2,0)*a[i]*b[i];
    FFT(f,-1);
    for(ri int i=m-1;i<n;i++) if(fabs(f[i].x/T)<0.5) pos[++cnt]=i-m+2;
    printf("%lld\n",cnt);
    for(ri int i=1;i<=cnt;i++) printf("%lld ",pos[i]);
    puts("");
    return 0;
}

考虑用 \(f_{i,c}\) 表示在第 \(i\) 个位置能否放 \(c\) 这个字符。

定义一个匹配函数 \(P(x)\) 表示 \(S\) 串中以第 \(x\) 位结束的长度为 \(m\) 的字符串与 \(T\) 有多少字符匹配。则我们可以考虑把匹配的字符分开考虑,最后统计一下总和是否等于 \(m\) 即可。

考虑构建一个函数 \(b_{i}\) 表示 \(t_{i}\) 这个位置是否等于当前匹配字符 \(c\)。再定义一个匹配函数 \(F(x,c)\),含义为 \(P(x)=F(x,A)+F(x,C)+F(x,G)+F(x,T)\),则 \(F(x,c)=\sum\limits_{i=0}^{m-1}b_{i}\times f_{x-m+1+i,c}\) 。考虑翻转 \(T\) 串,则 \(F(x,c)=\sum\limits_{i=0}^{m-1}b_{m-i-1}\times f_{x-m+1+i,c}\)。显然卷积形式。

#include <cstring>
#include <cmath>
#include <iostream>
#include <algorithm>
#include <queue>
#include <set>
#include <vector>
#include <stack>
#include <map>
#include <bitset>
#define ri register
#define inf 0x7fffffff
#define E (1)
#define mk make_pair
#define int long long
//#define double long double
using namespace std; const int N=1000010; const double pi=acos(-1.0);
inline int read()
{
    int s=0, w=1; ri char ch=getchar();
    while(ch<'0'||ch>'9') {if(ch=='-') w=-1; ch=getchar(); }
    while(ch>='0'&&ch<='9') s=(s<<3)+(s<<1)+(ch-'0'), ch=getchar();
    return s*w;
}
void print(int x) { if(x<0) x=-x, putchar('-'); if(x>9) print(x/10); putchar(x%10+'0'); }
int n,m,k,rev[N],T,f[N]; char s[N],t[N];
struct Complex
{
    double x,y;
    inline Complex ( double pp=0, double qq=0) { x=pp, y=qq; }
}a[N],b[N];
inline Complex operator + (Complex a,Complex b) { return Complex(a.x+b.x,a.y+b.y); }
inline Complex operator - (Complex a,Complex b) { return Complex(a.x-b.x,a.y-b.y); }
inline Complex operator * (Complex a,Complex b) { return Complex(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x); }
inline void Get_Rev() { for(ri int i=0;i<T;i++) rev[i]=(rev[i>>1]>>1)|((i&1)?(T>>1):0); }
inline void FFT(Complex *s,int type)
{
    for(ri int i=0;i<T;i++) if(i<rev[i]) swap(s[i],s[rev[i]]);
    for(ri int mid=1;mid<T;mid<<=1)
    {
        Complex wn(cos(pi/mid),type*sin(pi/mid));
        for(ri int r=mid<<1, j=0;j<T;j+=r)
        {
            Complex w(1,0);
            for(ri int k=0;k<mid;k++, w=w*wn)
            {
                Complex x=s[j+k], y=w*s[j+mid+k];
                s[j+k]=x+y, s[j+mid+k]=x-y;
            }
        }
    }
}
inline void S(char g)
{
    for(ri int i=0;i<T;i++) a[i]=Complex(0,0), b[i]=Complex(0,0);
    int pre=-1e9;
    for(ri int i=0;i<n;i++)
    {
        if(s[i]==g) pre=i;
        if(i-pre<=k) a[i]=Complex(1,0);
    }
    int nxt=1e9;
    for(ri int i=n-1;~i;i--)
    {
        if(s[i]==g) nxt=i;
        if(nxt-i<=k) a[i]=Complex(1,0);
    }
    for(ri int i=m-1;~i;i--) if(t[i]==g) b[m-i-1]=Complex(1,0);
    FFT(a,1), FFT(b,1);
    for(ri int i=0;i<T;i++) a[i]=a[i]*b[i];
    FFT(a,-1);
    for(ri int i=0;i<n;i++) f[i]+=(int)(a[i].x/T+0.5);
}
signed main()
{
    n=read(), m=read(), k=read();
    scanf("%s%s",s,t);
    T=1; while(T<=n) T<<=1;
    Get_Rev(); int res=0;
    S('A'), S('C'), S('T'), S('G');
    for(ri int i=0;i<n;i++) if(f[i]>=m) res++;
    printf("%lld\n",res);
    return 0;
}

例题:[luoguP4199]万径人踪灭

考虑容斥解决问题。答案为:所有满足字母和坐标都关于一条对称轴对称的字符串子序列数量 \(-\) 所有满足字母和坐标都关于一条对称轴对称的字符串子串数量。

后者显然可以用 \(\text{manacher}\) 算法在 \(O(n)\) 的时间内求出;对于前者,我们考虑记 \(F_{i}\) 表示关于第 \(i\) 个位置(包括字符位和字符与字符之间的位置)对答案的贡献,则存在:

\[\qquad F_{i}=\sum\limits_{j=0}^{i}[s_{j}=s_{2i-j}] \qquad \]

则显然第 \(i\) 个位置的答案为 \(2^{f_{i}}-1\),考虑转化一下上式。记 \(G_{i}=\sum\limits_{j=0}^{i}[s_{j}=s_{i-j}]\)。考虑到字符集大小为 \(2\),可以直接记 \(f_{i}=[s_{i}=\) a \(]\)\(g_{i}=[s_{i}=\) b \(]\),则得到 \(G=f*f+g*g\),则求出 \(G\)

现在考虑 \(F\)\(G\) 的关系。考虑一对满足 \(s_{x}=s_{y}\)\((x,y)\)\(x\not=y\),则它在对 \(G_{x+y}\) 的计算中会贡献两次。则当 \(2|G_{i}\) 时,\(F_{i}=\frac{G_{i}}{2}\)。否则会有一对满足 \(s_{x}=s_{y}\)\((x,y)\) 且满足 \(x=y\),这种点对在答案中应该被加入 \(F_{i}\),故当 \(2\not| G_{i}\) 时,\(F_{i}=\frac{G_{i}+1}{2}\)。故 \(F_{i}=\lfloor\frac{G_{i}+1}{2}\rfloor\)

然后直接 \(\text{FFT}\)\(\text{NTT}\) 解决即可。

#include <cstring>
#include <cmath>
#include <iostream>
#include <algorithm>
#include <queue>
#include <set>
#include <vector>
#include <stack>
#include <map>
#include <bitset>
#define ri register
#define inf 0x7fffffff
#define E (1)
#define mk make_pair
//#define int long long
//#define double long double
using namespace std; const int N=400010, Mod=1e9+7; const double pi=acos(-1.0);
inline int read()
{
    int s=0, w=1; ri char ch=getchar();
    while(ch<'0'||ch>'9') {if(ch=='-') w=-1; ch=getchar(); }
    while(ch>='0'&&ch<='9') s=(s<<3)+(s<<1)+(ch-'0'), ch=getchar();
    return s*w;
}
void print(int x) { if(x<0) x=-x, putchar('-'); if(x>9) print(x/10); putchar(x%10+'0'); }
int rev[N],T,n,res,f[N],g[N]; char s[N],t[N];
struct Complex
{
    double x,y;
    inline Complex (double pp=0, double qq=0) { x=pp, y=qq; }
}a[N];
inline Complex operator + (Complex a,Complex b) { return Complex(a.x+b.x,a.y+b.y); }
inline Complex operator - (Complex a,Complex b) { return Complex(a.x-b.x,a.y-b.y); }
inline Complex operator * (Complex a,Complex b) { return Complex(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x); }
inline void Get_Rev() { for(ri int i=0;i<T;i++) rev[i]=(rev[i>>1]>>1)|((i&1)?(T>>1):0); }
inline void FFT(Complex *s,int type)
{
    for(ri int i=0;i<T;i++) if(i<rev[i]) swap(s[i],s[rev[i]]);
    for(ri int mid=1;mid<T;mid<<=1)
    {
        Complex wn(cos(pi/mid),type*sin(pi/mid));
        for(ri int j=0, r=mid<<1;j<T;j+=r)
        {
            Complex w(1,0);
            for(ri int k=0;k<mid;k++, w=w*wn)
            {
                Complex x=s[j+k], y=w*s[j+mid+k];
                s[j+k]=x+y, s[j+mid+k]=x-y;
            }
        }
    }
}
inline int ksc(int x,int p) { int res=1; for(;p;p>>=1, x=1ll*x*x%Mod) if(p&1) res=1ll*res*x%Mod; return res; }
signed main()
{
    scanf("%s",s); n=strlen(s); int cnt=0;
    t[0]=t[1]='#';
    for(ri int i=0;i<n;i++) t[i*2+2]=s[i],t[i*2+3]='#';
    cnt=(n+1)*2;
    int mid,mr=0;
    for(ri int i=1;i<cnt;i++)
    {
        if(i<mr) f[i]=min(f[mid*2-i],f[mid]+mid-i);
        else f[i]=1;
        while(t[i-f[i]]==t[i+f[i]]) f[i]++;
        if(mr<i+f[i]) mr=i+f[i], mid=i;
    }
    for(ri int i=0;i<cnt;i++) res=(res-f[i]/2+Mod)%Mod;
    T=1; while(T<=cnt) T<<=1; Get_Rev();
    for(ri int i=0;i<n;i++) a[i]=Complex(s[i]=='a',0);
    FFT(a,1);
    for(ri int i=0;i<T;i++) a[i]=a[i]*a[i];
    FFT(a,-1);
    for(ri int i=0;i<T;i++) g[i]+=(int)(a[i].x/T+0.5), g[i]%=Mod;
    for(ri int i=0;i<T;i++) a[i]=Complex(0,0);
    for(ri int i=0;i<n;i++) a[i]=Complex(s[i]=='b',0);
    FFT(a,1);
    for(ri int i=0;i<T;i++) a[i]=a[i]*a[i];
    FFT(a,-1);
    for(ri int i=0;i<T;i++) g[i]+=(int)(a[i].x/T+0.5), g[i]%=Mod;
    for(ri int i=0;i<T;i++)
    {
        int p=(g[i]+1)/2;
        res=(res+ksc(2,p)-1)%Mod;
    }
    printf("%d\n",res);
    return 0;
}
posted @ 2020-08-03 15:23  zkdxl  阅读(213)  评论(1编辑  收藏  举报