• 【xsy1154】 DNA配对 FFT


    题目大意:给你一个字符串$s$和字符串$w$,字符集为${A,T,C,G}$,你要在字符串$s$中选出一个与$w$长度相同的子串,使得这两个串的差异度最小。

    两个字符$c1$,$c2$的差异度为给定的$c[c1][c2]$。

    字符串长度$≤2*10^5$。

    $FFT$套路题。

    我们将串$w$翻转。

    设$p[i]$为$s$中子串$s[i-|w|+1.......i]$与$w$的差异度。

    显然$p[i]=sum_{j=0}^{i} c[s[j]][w[i-j]]$。(此处的$w$是翻转后的)

    显然的卷积形式。

    五次$FFT$即可。

     1 #include<bits/stdc++.h>
     2 #define M (1<<19)
     3 #define PI acos(-1)
     4 #define INF 19890604
     5 using namespace std;
     6 
     7 struct cp{
     8     double i,r; cp(){i=r=0;}
     9     cp(double rr,double ii){i=ii;r=rr;}
    10     friend cp operator +(cp a,cp b){return cp(a.r+b.r,a.i+b.i);}
    11     friend cp operator -(cp a,cp b){return cp(a.r-b.r,a.i-b.i);}
    12     friend cp operator *(cp a,cp b){return cp(a.r*b.r-a.i*b.i,a.r*b.i+a.i*b.r);}
    13     friend cp operator /(cp a,double b){return cp(a.r/b,a.i/b);}
    14     int num(){return (int)(i+0.499)/2;}
    15 }a[M],b[M],c[M],d[M],ans[M];
    16 
    17 void change(cp a[],int len){
    18     for(int i=0,j=0;i<len-1;i++){
    19         if(i<j) swap(a[i],a[j]);
    20         int k=len>>1;
    21         while(j>=k) j-=k,k>>=1;
    22         j+=k;
    23     }
    24 }
    25 void FFT(cp a[],int len,int on){
    26     change(a,len);
    27     for(int h=2;h<=len;h<<=1){
    28         cp wn=cp(cos(2*on*PI/h),sin(2*on*PI/h));
    29         for(int j=0;j<len;j+=h){
    30             cp w=cp(1,0);
    31             for(int k=j;k<(j+(h>>1));k++){
    32                 cp u=a[k],t=w*a[k+(h>>1)];
    33                 a[k]=u+t; a[k+(h>>1)]=u-t;
    34                 w=w*wn;
    35             }
    36         }
    37     }
    38     if(on==-1) for(int i=0;i<len;i++) a[i]=a[i]/len;
    39 }
    40 
    41 int D[4][4]={0};
    42 int get(char c){
    43     if(c=='A') return 0;
    44     if(c=='T') return 1;
    45     if(c=='C') return 2;
    46     if(c=='G') return 3;
    47 }
    48 
    49 char s[M]={0},w[M]={0};
    50 int n,m;
    51 int main(){
    52     scanf("%s%s",s,w); n=strlen(s); m=strlen(w);
    53     for(int i=0;i<16;i++) scanf("%d",D[0]+i);
    54     for(int i=0;i<n;i++) s[i]=get(s[i]);
    55     for(int i=0;i<m;i++) w[i]=get(w[i]);
    56     reverse(w,w+m);
    57     for(int i=0;i<n;i++){
    58         if(s[i]==0) a[i].i=1;
    59         if(s[i]==1) b[i].i=1;
    60         if(s[i]==2) c[i].i=1;
    61         if(s[i]==3) d[i].i=1;
    62     }
    63     for(int i=0;i<m;i++){
    64         a[i].r=D[0][w[i]];
    65         b[i].r=D[1][w[i]];
    66         c[i].r=D[2][w[i]];
    67         d[i].r=D[3][w[i]];
    68     }
    69     int len=1; while(len<n+m) len<<=1;
    70     FFT(a,len,1); FFT(b,len,1); FFT(c,len,1); FFT(d,len,1);
    71     for(int i=0;i<len;i++){
    72         ans[i]=a[i]*a[i]+b[i]*b[i]+c[i]*c[i]+d[i]*d[i];
    73     }
    74     FFT(ans,len,-1);
    75     int minn=INF;
    76     for(int i=m-1;i<n;i++) 
    77     minn=min(minn,ans[i].num());
    78     cout<<minn<<endl;
    79 }
  • 相关阅读:
    链表和递归
    Async and Await
    Linux下如何对目录中的文件进行统计
    Linux系统添加自定义网卡并更改网卡接口
    运维实战:Linux系统扩展oracle数据库所在的分区
    在Linux命令行发送电子邮件附件的两种方法
    Centos7部署轻量级自动化运维工具pssh
    CentOS7中使用systemctl列出启动失败的服务
    AbstractQueuedSynchronizer详解(转)
    Oracle 11g数据导入到10g
  • 原文地址:https://www.cnblogs.com/xiefengze1/p/9806997.html
Copyright © 2020-2023  润新知