POJ 2778 DNA Sequence

最新推荐文章于 2021-03-22 22:58:59 发布

原创最新推荐文章于 2021-03-22 22:58:59 发布 · 485 阅读

1 ·

CC 4.0 BY-SA版权

文章标签：

#build #insert #algorithm #struct #c #c++

POJ 同时被 3 个专栏收录

75 篇文章

订阅专栏

动规

33 篇文章

订阅专栏

AC自动机

6 篇文章

订阅专栏

本文介绍如何使用矩阵快速幂方法优化大规模DNA序列匹配问题。通过构建AC自动机并结合动态规划，实现高效的模式串匹配算法。文章详细阐述了算法原理及实现细节。

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

DNA Sequence

自动机的又一题。这个题目和前面的那个题目POJ1625很相似，不过POJ1625的数据范围没有这么大，如果按照那个题目那样递推过来的话，肯定会超时。对于这么大的数据，一般会想到加速计算，多采用二分法，例如：快速幂取模，矩阵快速幂。而这个题目就是用矩阵快速幂来加速计算的。和前面一样我们可以很容易的写出状态转移方程，其实是一样的方程，这里就不再写了。我们把这个方程写成方程组的形式，那么我们就可以将这些方程组写成矩阵相乘的形式，那么我们就得到了通项公式，然后就是矩阵快速幂了。如果用递归的形式，貌似会RE，我写了一次结果RE了，换成非递归的就可以过了。当然如果取模次数太多的话，也是很容易超时的，取模的速度也是很慢的= . =，开始的时候貌似就在这个上面超了......

/*
author: csuchenan
LANG  : c++
PROG  : POJ2778
ALGORITHM: AC trie + DP + martix
*/

#include <cstdio>
#include <cstring>
#include <queue>
using namespace std ;
const int mod = 100000 ;

int N ;
int L ;

struct Martix{
    long long a[110][110] ;
    int kbit ;
};

int next[110][4] ;
int fail[110] ;
bool flag[110] ;
int cnt ;

void debug(Martix E){
    for(int i = 0 ; i <=E.kbit ; i ++){
        printf("%d:" , i) ;
        for(int j = 0 ;j <= E.kbit ; j ++){
            printf("%d " , E.a[i][j]) ;
        }
        printf("\n") ;
    }
}
inline void initMartix(Martix &E){
    memset(E.a , 0 , sizeof(E.a)) ;
    for(int i = 0 ; i < 110 ; i ++)
        E.a[i][i] = 1 ;
}

inline void init(){
    memset(next , 0 , sizeof(next)) ;
    memset(fail , 0 , sizeof(fail)) ;
    memset(flag , 0 , sizeof(flag)) ;
    cnt = 0 ;
}

int hash(char c){
    switch(c){
        case 'A' : return 0 ;
        case 'C' : return 1 ;
        case 'G' : return 2 ;
        case 'T' : return 3 ;
    }
}

void insert(char * str){
    char *p = str ;
    int b;
    int c(0) ;

    while(*p){
        b = hash(*p) ;
        if(!next[c][b])
            next[c][b] = ++ cnt ;
        c = next[c][b] ;
        p ++ ;
    }
    flag[c] = 1 ;
}

void build_ac(){
    queue<int> Q ;
    int cur(0) ;
    int child  ;
    fail[cur] = 0 ;
    Q.push(cur) ;

    while( !Q.empty() ){
        cur = Q.front() ;
        Q.pop() ;

        for(int i = 0 ; i < 4 ; i ++){
            child = next[cur][i] ;
            if(child){
                Q.push(child) ;

                if(!cur)
                    fail[cur] = 0 ;
                else{
                    int tmp = fail[cur] ;
                    for(;tmp && !next[tmp][i] ; tmp = fail[tmp])
                        ;
                    if(next[tmp][i])
                        fail[child] = next[tmp][i] ;
                    else
                        fail[child] = 0 ;
                }
                if(flag[fail[child]])
                    flag[child] = 1 ;
            }
            else{
                next[cur][i] = next[fail[cur]][i] ;
            }
        }
    }
}

Martix build_martix(){
    Martix ans ;
    ans.kbit = cnt ;
    memset(ans.a , 0 , sizeof(ans.a)) ;
    for(int i = 0 ; i <= cnt ; i ++){
        if(flag[i])
            continue ;
        for(int k = 0 ; k < 4 ; k ++){
            int cur = next[i][k] ;
            if(flag[cur])
                continue ;
            ans.a[ cur ][ i ] ++ ;
        }
    }
    return ans ;
}

Martix mul(Martix A , Martix B){
    int l = A.kbit ;
    Martix ans ;
    ans.kbit = l ;
    for(int i = 0 ; i <= l ; i ++){
        for(int j = 0 ; j <= l ; j ++){
            int total ;
            total = 0 ;
            for(int k = 0 ; k <= l ; k ++){
                total =( total + A.a[i][k]*B.a[k][j])%mod ;
            }
            ans.a[i][j] = total ;
        }
    }
    return ans ;
}

Martix pow_martix(Martix A , int k){
    Martix ans  ;
    initMartix(ans) ;
    ans.kbit = A.kbit;
    while(k){
        if(k&1)
            ans = mul(ans , A) ;
        A = mul(A , A) ;
        k=k>>1 ;
    }
    return ans ;
}
/*
Martix pow_martix(Martix A , int k){
    Martix ans ;
    ans.kbit = A.kbit ;
    initMartix(ans) ;

    if(k==0){
        return ans ;
    }

    ans = pow_martix(A , k>>1) ;
    ans = mul(ans , ans) ;
    if(k&1)
        ans = mul(ans , A) ;
    return ans ;
}*/

void read(){
    char str[15] ;
    scanf("%d%d" , &N , &L) ;
    for(int k = 1 ; k <= N ; k ++){
        scanf("%s" , str) ;
        insert(str) ;
    }
}

inline void print(Martix A){
    int ans(0);
    for(int k = 0 ; k <= A.kbit ; k ++){
        if(flag[k])
            continue ;
        ans = ans + A.a[k][0] ;
    }
    printf("%d\n" , ans%mod) ;
}
int main(){
    //init() ;
    read() ;
    build_ac() ;
    Martix E = build_martix() ;
    Martix ans = pow_martix(E , L) ;
    print(ans);
    return 0 ;
}