POJ2778 - DNA Sequence - AC自动机+矩阵快速幂

DNA Sequence

题目链接

分类：data structures strings

1.题意概述

• 有M(1<=M<=10)种DNA序列是表示有“病症”的，现在要你构造长度为N(1<=n<=2000000000)的DNA序列，并且它不包含任何一种有疾病的DNA序列。

2.解题思路

• 对于样例病毒序列$\{ACG,C\}$我们构造AC自动机（Tire树）的转移状态如图：

  

  也就是从状态0出发，一共有四步走法：

  • 走A到状态1——安全
  • 走C到状态4——危险
  • 走T到状态0——安全
  • 走G到状态0——危险

  所以当n=1时候，答案就是3，当n=2时候，就是从状态0出发，走2步，形成一个长度为2的字符串，只要路径上没有经过危险节点，答案就这几种走法，我们可以再考虑建立迭代矩阵：$M_{i,j}$表示从节点i到节点j只走一步的走法数，那么迭代n次($M^n$)就表示从节点i到节点j走n步的走法，这个可以用矩阵快速幂解决，答案就是$\sum\limits_{i=1}^nM^i$，但是这里n很大，我们可以考虑矩阵加一维度表示前i次的和，具体思路详见矩阵快速幂专题！

3.AC代码

class Matrix {
public: 
    ll a[140][140];
    int n;

    void Init(int key) {
        memset(a, 0, sizeof a);
        if (key) 
            rep(i, 0, n) a[i][i] = 1;
    }

    Matrix operator+ (const Matrix &b) {
        Matrix c;
        c.n = n;
        rep(i, 0, n) rep(j, 0, n) c.a[i][j] = (a[i][j] + b.a[i][j]) % mod;
        return c;
    }

    Matrix operator+ (int x) {
        Matrix p = *this;
        rep(i, 0, n) p.a[i][i] = (p.a[i][i] + x % mod) % mod;
        return p;
    }

    Matrix operator- (const Matrix &b) {
        Matrix c;
        c.n = n;
        rep(i, 0, n) rep(j, 0, n) c.a[i][j] = (a[i][j] - b.a[i][j] + mod) % mod;
        return c;
    }

    Matrix operator* (const Matrix &b) {
        Matrix c;
        c.n = n;
        c.Init(0);
        rep(i, 0, n) rep(j, 0, n) rep(k, 0, n)
            c.a[i][j] = (c.a[i][j] + a[i][k] * b.a[k][j]) % mod;
        return c;
    }

    Matrix Power(int t) {
        Matrix ans, p = *this;
        ans.n = p.n;
        ans.Init(1);
        while (t) {
            if (t & 1) ans = ans * p;
            p = p * p;
            t >>= 1;
        }
        return ans;
    }

    void Print() {
        rep(i, 0, n) {
            rep(j, 0, n) {
                if (j == 0) printf("%I64d", a[i][j]);
                else printf(" %I64d", a[i][j]);
            }
            puts("");
        }
    }   
};
struct Trie {
    int next[140][4], fail[140];
    bool end[140];
    int root, L;

    int Newnode() {
        rep(i, 0, 4) next[L][i] = -1;
        end[L++] = 0;
        return L - 1;
    }

    void Init() {
        L = 0;
        root = Newnode();
        memset(end, 0, sizeof end);
    }

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

    void Insert(char *s) {
        int len = strlen(s);
        int now = root;
        rep(i, 0, len) {
            int index = getIndex(s[i]);
            if(next[now][index] == -1)
                next[now][index] = Newnode();
            now = next[now][index];
        }
        end[now] = 1;
    }

    void Build() {
        queue<int> Q;
        fail[root] = root;
        rep(i, 0, 4) {
            if(next[root][i] == -1)
                next[root][i] = root;
            else {
                fail[next[root][i]] = root;
                Q.push(next[root][i]);
            }
        }
        while(!Q.empty()) {
            int now = Q.front();
            Q.pop();
            if (end[fail[now]])
                end[now] = 1;
            rep(i, 0, 4) {
                if(next[now][i] == -1)
                    next[now][i] = next[fail[now]][i];
                else {
                    fail[next[now][i]] = next[fail[now]][i];
                    Q.push(next[now][i]);
                }
            }
        }
    }

    Matrix GetMatrix() {
        Matrix ret;
        ret.n = L; ret.Init(0);
        rep(i, 0, L) {
            rep(j, 0, 4) {
                if (!end[i] && !end[next[i][j]])
                    ret.a[i][next[i][j]]++;
            }
        }
        return ret;
    }
} AC;
char ch[10];
inline void solve() {
    int n, l;
    while (~scanf("%d%d\n", &n, &l)) {
        AC.Init();
        rep(i, 0, n) {
            gets(ch);
            AC.Insert(ch);
        }
        AC.Build();
        Matrix f = AC.GetMatrix();
    //  f.Print();
        f = f.Power(l);
        ll ans = 0;
        rep(i, 0, f.n) ans = (ans + f.a[0][i]) % mod;
        printf("%I64d\n", ans);
    }
}