CF528D Fuzzy Search（NTT/FFT）

CF528D Fuzzy Search

前言

这个方法不对劲，但是是对的qwq。

我一直在想为什么是 A,C,G,T 四个字符。现在我想通了。这场 CF 大概是 ATcoder 的 AGC 出题人出的（大雾）

解法

首先我们要考虑 \(k=0\) 的时候怎么做。直接 KMP 即可。

但是这样我们没法做 \(k>0\) 的情况。我们考虑用 FFT/NTT。

容易发现题中只有 \(4\) 中字符。我们考虑一个一个处理。假设对于样例，我们先处理 A 字符。

那么我们有：

\[S=A**AA***A*\\T=A*A* \]

其中 * 是通配符，可以为任意字符。

那么怎么处理 \(k\) 呢？

我们容易想到，可以把一个原来是 A 的位置往两旁延伸 \(k\) 个位置，把它们都置为 A。这样 \(k\) 就没有影响啦。

举个例子，对于样例我们就有：

\[S'=AAAAAA*AAA\\T=A*A* \]

为了方便计算，我们可以认为：

\[S'=1111110111\\T=1010 \]

把 \(S'\) 和 \(T\) 转化为多项式 \(f,g\)，那么什么情况下在 \(i\) 的位置珂以成功匹配呢？

如果我们有：

\[h_i=\sum\limits_{j=1}^{|T|} g_j\cdot(f_{i+j}-g_{j})^2 \]

当 \(h_i\) 为 \(0\) 时，在 \(i\) 上可以匹配。

怎样快速求出 \(h_i\) 呢？展开得：

\[h_i=\sum\limits_{j=1}^{|T|} g_j\cdot {f_{i+j}}^2+{g_j}^3-2\cdot {g_j}^2\cdot f_{i+j} \]

我的第一反应是：完了，处理不了。

仔细一想，其实是对的。

由于 \(g,f\) 每一项的取值为 \(0\) 或 \(1\)，所以每个平方或者三次方项都可以消成一次项，不影响计算，因此转化为:

\[h_i=\sum\limits_{j=1}^{|T|} g_j\cdot {f_{i+j}}+{g_j}-2\cdot {g_j}\cdot f_{i+j}\\h_i=\sum\limits_{j=1}^{|T|} {g_j}-{g_j}\cdot f_{i+j} \]

将 \(g\) 翻转就可以求卷积了。

最后四个字符分别处理对应的 \(h_i\)，都为 \(0\) 则这个位置可以匹配。

代码

//Don't act like a loser.
//This code is written by huayucaiji
//You can only use the code for studying or finding mistakes
//Or,you'll be punished by Sakyamuni!!!
#include<bits/stdc++.h>
#define int long long
using namespace std;

int read() {
	char ch=getchar();
	int f=1,x=0;
	while(ch<'0'||ch>'9') {
		if(ch=='-')
			f=-1;
		ch=getchar();
	}
	while(ch>='0'&&ch<='9') {
		x=x*10+ch-'0';
		ch=getchar();
	}
	return f*x;
}

const int MAXN=4e6+10;
const int G=3;
const int MOD=998244353;

int n,m,k;
int d=1,bit=0;
short index[30];
int s[MAXN],t[MAXN],f[MAXN],g[MAXN],r[MAXN],ans[MAXN];

int qpow(int x,int y) {
	int ret=1;
	while(y) {
		if(y&1) {
			ret=x*ret%MOD;
		}
		y>>=1;
		x=x*x%MOD;
	}
	return ret;
}

void NTT(int A[],int d,int inv) {
	for(int i=0;i<d;i++) {
		if(i<r[i]) {
			swap(A[i],A[r[i]]);
		}
	}
	
	for(int len=1;len<d;len<<=1) {
		int gn=qpow(G,(MOD-1)/(len*2));
		for(int i=0;i<d;i+=len<<1) {
			int g=1;
			for(int k=0;k<len;k++,g=g*gn%MOD) {
				int x=A[i+k];
				int y=g*A[i+k+len]%MOD;
				A[i+k]=(x+y)%MOD;
				A[i+k+len]=(x-y+MOD)%MOD;
			}
		}
	}
	
	if(inv==-1) {
		reverse(A+1,A+d);
		int g=qpow(d,MOD-2);
		for(int i=0;i<d;i++) {
			A[i]=A[i]*g%MOD;
		}
	}
}

void FFT_KMP() {
	int cst=0;
	for(int i=0;i<m;i++) {
		cst+=g[i];
	}
	reverse(g,g+m);
	NTT(f,d,1);
	NTT(g,d,1);
	for(int i=0;i<d;i++) {
		f[i]=f[i]*g[i]%MOD;
	}
	NTT(f,d,-1);
	for(int i=0;i<n;i++) {
		ans[i]+=f[m+i-1]+cst-2*f[m+i-1];
	}
}

signed main() {
	index['A'-'A']=1;index['C'-'A']=2;index['G'-'A']=3;index['T'-'A']=4;
	n=read();m=read();k=read();
	
	for(int i=0;i<n;i++) {
		char c;
		cin>>c;
		s[i]=index[c-'A'];
	}
	for(int i=0;i<m;i++) {
		char c;
		cin>>c;
		t[i]=index[c-'A'];
	}
	
	while(d<=n+m) {
		d<<=1;
		bit++;
	}
	for(int i=0;i<d;i++) {
		r[i]=(r[i>>1]>>1)|((i&1)<<(bit-1));
	}
	
	for(int j=1;j<=4;j++) {
		fill(f,f+d,0);
		fill(g,g+d,0);
		int p=0;
		for(int i=0;i<n;i++) {
			if(s[i]==j) {
				p=k+1;
			}
			if(p) {
				f[i]=1;
				p--;
			}
		}
		for(int i=n-1;i>=0;i--) {
			if(s[i]==j) {
				p=k+1;
			}
			if(p) {
				f[i]=1;
				p--;
			}
		}
		for(int i=0;i<m;i++) {
			if(t[i]==j) {
				g[i]=1;
			}
		}
		FFT_KMP();
	}
	
	int tot=0;
	for(int i=0;i<n;i++) {
		tot+=(ans[i]==0? 1:0);
	}
	cout<<tot<<endl;
	
	return 0;
}