CF528D. Fuzzy Search [FFT]

本文介绍了解决CF528D.FuzzySearch问题的方法,通过预处理和使用bitset进行暴力匹配来提高效率。文章详细解释了如何对DNA序列进行匹配,并利用卷积技巧简化计算过程。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >

题意:DNA序列,在母串s中匹配模式串t,对于s中每个位置i,只要s[i-k]到s[i+k]中有c就认为匹配了c。求有多少个位置匹配了t


预处理\(f[i][j]\)表示位置i可以匹配字符j
分别考虑每一个字符c,对s的每个位置i求出用\(s[i,i+m-1]\)匹配t,这个字符匹配了几次
\(a_i=[s的位置i匹配c],\ b_i=[t_i=c]\)
那么c的匹配次数就是\(c_j=\sum\limits_{i=0}^{m-1}a_{j+i}b_i\),位置i匹配了t当且仅当四种字符的匹配次数和等于t的长度m


这时候就可以考虑bitset暴力过了

一个常用技巧是,反转模式串(或母串),然后就成了卷积的形式:
\[ c_j=\sum\limits_{i=0}^{m-1}a_{j+i}b_{m-1-i}=D_{m+j-1} \]
这样计算是没有问题的,因为b只有\([0,m-1]\)有值其他地方为0

注意处理每个字符前memset a和b!!!!!

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
typedef long long ll;
const int N=(1<<20)+5, INF=1e9;
const double PI=acos(-1);
inline int read(){
    char c=getchar();int x=0,f=1;
    while(c<'0'||c>'9'){if(c=='-')f=-1;c=getchar();}
    while(c>='0'&&c<='9'){x=x*10+c-'0';c=getchar();}
    return x*f;
}

struct meow{
    double x, y;
    meow(double a=0, double b=0):x(a), y(b){}
};
meow operator +(meow a, meow b) {return meow(a.x+b.x, a.y+b.y);}
meow operator -(meow a, meow b) {return meow(a.x-b.x, a.y-b.y);}
meow operator *(meow a, meow b) {return meow(a.x*b.x-a.y*b.y, a.x*b.y+a.y*b.x);}
meow conj(meow a) {return meow(a.x, -a.y);}
typedef meow cd;

namespace FFT{
    int n, rev[N];
    void ini(int lim) {
        n=1; int k=0;
        while(n<lim) n<<=1, k++;
        for(int i=0; i<n; i++) rev[i] = (rev[i>>1]>>1) | ((i&1)<<(k-1));
    }
    void dft(cd *a, int flag) {
        for(int i=0; i<n; i++) if(i<rev[i]) swap(a[i], a[rev[i]]);
        for(int l=2; l<=n; l<<=1) {
            int m=l>>1; 
            cd wn = meow(cos(2*PI/l), flag*sin(2*PI/l));
            for(cd *p=a; p!=a+n; p+=l) {
                cd w(1, 0);
                for(int k=0; k<m; k++) {
                    cd t = w*p[k+m];
                    p[k+m] = p[k] - t;
                    p[k] = p[k] + t;
                    w=w*wn;
                }
            }
        }
        if(flag==-1) for(int i=0; i<n; i++) a[i].x/=n;
    }
    void mul(cd *a, cd *b) {
        dft(a, 1); dft(b, 1);
        for(int i=0; i<n; i++) a[i]=a[i]*b[i];
        dft(a, -1);
    }
}using FFT::mul; using FFT::ini;

int n, m, k, lim, f[N][5], cnt[5], id[300];
cd a[N], b[N], c[N];
char s[N], t[N];
int ans[N];
void solve(int now) { 
    memset(a, 0, sizeof(a)); memset(b, 0, sizeof(b));
    for(int i=0; i<n; i++) a[i].x = f[i][now];
    for(int i=0; i<m; i++) b[m-1-i].x = id[(int)t[i]]==now; 
    mul(a, b);
    for(int i=0; i<n; i++) ans[i] += int(a[m-1+i].x+0.5);
}
int main() {
    freopen("in","r",stdin);
    n=read(); m=read(); k=read(); 
    lim=n+m-1; ini(lim);
    scanf("%s%s",s,t);
    id['A']=0; id['T']=1; id['C']=2; id['G']=3;
    int l=0, r=0; cnt[ id[(int)s[0]] ]++;
    for(int i=0; i<n; i++) {
        while(l<i-k) cnt[ id[(int)s[l++]] ]--;
        while(r<n-1 && r<i+k) cnt[ id[(int)s[++r]] ]++;
        for(int j=0; j<4; j++) if(cnt[j]) f[i][j]=1;
    }
    for(int i=0; i<4; i++) solve(i);
    int sum=0;
    for(int i=0; i<n; i++) if(ans[i]==m) sum++;
    printf("%d",sum);
}
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值