Poj 2778 [AC自动机,矩阵乘法]

It’s well known that DNA Sequence is a sequence only contains A, C, T and G, and it’s very useful to analyze a segment of DNA Sequence,For example, if a animal’s DNA sequence contains segment ATC then it may mean that the animal may have a genetic disease. Until now scientists have found several those segments, the problem is how many kinds of DNA sequences of a species don’t contain those segments.

Suppose that DNA sequences of a species is a sequence that consist of A, C, T and G,and the length of sequences is a given integer n.

题意:有m种DNA序列是有疾病的,问有多少种长度为n的DNA序列不包含任何一种有疾病的DNA序列。(仅含A,T,C,G四个字符)
dp[i][j]=sum(dp[i-1][k])
枚举第i位填了’A’,’C’,’G’,’T’中的哪一个,自动机从k走到了哪里。
把i看作行,j看做列,那么dp[i][]和dp[i-1][]之间的关系显然和dp[i+1][]和dp[i-1][]的关系相同。
因此,可以用构造转移矩阵方法优化。。。矩阵就是dp【i][j]表示从I到j结点走一步(加一个字母)有几种方法,注意不能走到任何模式串的尾巴也不能它来转移,所以相当于尾结点的行列都不要,然后n次方,就是加了n个字母(题目要求了长度),初始化时候,只把0 0设置为1(没有字母),最后所有(0,i)的值加起来,就是0号结点转移到I的。
注意m=0特判。。。

#include <cstdio>
#include <iostream>
#include <cstring>
#include <queue>
using namespace std;
int n,m,node;
const int maxn = 15;
const long long MOD = 100000;
char ss[maxn][maxn];
int val[maxn*10],last[maxn*10],ch[maxn*10][maxn],fail[maxn*10],kmax;

int getpos(char a){
    if(a == 'A')    return 0;
    if(a == 'T')    return 1;
    if(a == 'C')    return 2;
    if(a == 'G')    return 3;
}

void Insert(int v){
    int len = strlen(ss[v]);
    int tmp = 0;
    for(int i = 0; i < len ; i++){
        int pos = getpos(ss[v][i]);
        if(ch[tmp][pos] == 0){
            memset(ch[node],0,sizeof(ch[node]));
            val[node] = 0;
            ch[tmp][pos] = node++;
        }
        tmp = ch[tmp][pos];
    }
    val[tmp] = v;
}

void getfail(){
    queue<int> q;
    int u = 0;
    for(int i = 0; i < 4; i++){
        if(ch[u][i]){
            int tmp = ch[u][i];
            fail[tmp] = 0;
            q.push(tmp);
            last[tmp] = 0;
        }
    }

    while(!q.empty()){
        int tmp = q.front();
        q.pop();
        for(int i = 0; i < 4; i++){
            if(ch[tmp][i]){
                int v = fail[tmp];
                while(v && ch[v][i]==0) v = fail[v];
                fail[ch[tmp][i]] = ch[v][i];
                q.push(ch[tmp][i]);
               // if(ch[tmp][i]==2)   printf("asd   ch[v][i] = %d \n",ch[v][i]);
               last[ch[tmp][i]] = val[ch[v][i]]?ch[v][i] : last[ch[v][i]];
            }
        }
    }
}

void init(){
    node = 1;
    kmax = 0;
    memset(ch[0],0,sizeof(ch[0]));
    for(int i = 1; i <= m ; i++){
        scanf("%s",ss[i]);
        Insert(i);
    }

}
struct Matrix{
    long long m[maxn*10][maxn*10];
}dp,ini;

void print(Matrix a){
    for(int i = 0 ; i < 15 ; i++){
        for(int j = 0 ; j < 15 ; j++)
            printf("%lld ",a.m[i][j]);
        printf("\n");
    }
    cout <<"\n";
}

void make_matrix(){
   // cout <<"node = " << node <<endl;
    memset(dp.m,0,sizeof(dp.m));
    for(int i = 0 ; i < node; i++){
        if(val[i]||val[last[i]])  continue;
        for(int j = 0 ; j < 4; j++){
            //printf("i = %d  j = %d \n",i,j);
            //printf("last[%d]  = %d    val = %d\n",ch[i][j],last[ch[i][j]],val[ch[i][j]]);
            if(last[ch[i][j]]||val[ch[i][j]])  continue;
            int v = i;

            if(ch[v][j]!=0){
                dp.m[i][ch[v][j]] += 1;
                kmax = max(kmax,max(i,ch[v][j]));
            }
            else{
                while(v && ch[v][j] == 0 )    v = fail[v];
                if(last[ch[v][j]]||val[ch[v][j]])  continue;
                dp.m[i][ch[v][j]] += 1;
                kmax = max(kmax,max(i,ch[v][j]));
            }
        }
    }
   // print(dp);
}


Matrix Mul(Matrix a,Matrix b){
    Matrix p;
    for(int i = 0 ;i <= kmax ; i++)
        for(int j = 0 ; j <= kmax; j++)
            p.m[i][j] = 0;
    for(int i =  0 ;i <= kmax ;i++){
        for(int j = 0 ; j <= kmax ; j++){
                for(int k = 0 ; k <= kmax ; k++)
                    p.m[i][j] = (p.m[i][j]+a.m[i][k]*b.m[k][j])%MOD;
        }
    }
    return p;
}

Matrix pow(Matrix a, int n){
    Matrix p;
     for(int i = 0 ;i <= kmax ; i++)
        for(int j = 0 ; j <= kmax; j++)
            p.m[i][j] = 0;
    for(int i = 0 ; i <= kmax ; i++)  p.m[i][i] = 1;
    while(n){
        if(n & 1)   p = Mul(p,a);
        n >>= 1;
        a = Mul(a,a);
    }
   // print(p);
    return p;
}

void sov(){
     for(int i = 0 ;i <= kmax ; i++)
        for(int j = 0 ; j <= kmax; j++)
            ini.m[i][j] = 0;
    ini.m[0][0] = 1;
    //cout <<"n = "<<n<<endl;
    ini =  Mul(ini,pow(dp,n));
    int ans = 0;
    for(int i = 0 ; i < node; i++)
        ans = (ans + ini.m[0][i])%MOD;
    printf("%d\n",ans);
}
/*
3 10
AT
GATC
TGAC
*/
int qpow(int n){
    int a = 1,p = 4;
    while(n){
        if(n & 1)   a = (a*p)%MOD;
        n >>= 1;
        p = (p*p)%MOD;
    }
    return a;
}

int main(){
    while(~scanf("%d%d",&m,&n)){
        if(m == 0){
           // cout <<"Asd"<<endl;
            printf("%d\n",qpow(n));
            continue;
        }
        init();
        getfail();
        make_matrix();
        sov();
    }
    return 0;
}
  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值