DNA Sequence POJ - 2778 邻接矩阵 trie图 矩阵快速幂

Posted 2021-03-22 lordxx

首先构造trie图。
我们明确一点的是，给出trie图，那么所有点的转移方式都是唯一可以确定的。即使是没有这个字符，他也会指向根节点。
我们根据离散数学的知识可以知道。计算有向图的邻接矩阵，然后k次方，就能够计算出从某一个点到另一个点，有多少条长度为k的路径。
故，我们构造出来trie图，拿出该图的邻接矩阵，就能计算路径数目。--（注意改图是有向图）--
trie图的构造不说了，模板。

邻接矩阵的构造根据trie图来的。我们在trie图上找到每一个节点，查看他的相邻节点，即A,G,C,T四个点指向的节点。如果能走到，矩阵为1，不能为0。
什么情况不能走到，有两种情况，1是它自身，就是结尾的节点。2是，他的父亲节点是不能走到的，例如给出AGCCT,GC AGC的C节点fail指针指向了GC，表明以C结尾的字符串，有一个GC是危险的，不能走到的，那么AGCCT这个字符串字串AGC就是危险的，AGCC，AGCCT都是危险的。故不能走到。

再继续判断该算法的合理性。我们其实害怕会统计重复。我们思考一下会重复统计么。以远点出发，表示当前串，一开始为空，我们继续以AGCCT,GC为例子。
我们走到GC节点和走到AGC节点，是不会重复统计的。虽然都是GC结尾，但是因为trie图状态唯一确定，走到单独的GC节点，那么该GC节点的前一位，一定不会是A，只能是GCT，因为如果是A，那么根据trie图的性质，肯定会走到AGC这条链上去。故统计不会重复。

（危险的点，他的行列都是0，故不可能有任何点能通过该点，故即使是快速幂k次方后，仍然不会有任何节点通过危险节点能走到其他节点，也不会走到他的子节点，即有向边指向的邻接点。）

#include<iostream>
#include<cstring>
#include<algorithm>
#include<cmath>
#include<cstdlib>
#include<climits>
#include<stack>
#include<vector>
#include<queue>
#include<set>
#include<bitset>
#include<map>
//#include<regex>
#include<cstdio>
#pragma GCC optimize(2)
#define up(i,a,b)  for(int i=a;i<b;i++)
#define dw(i,a,b)  for(int i=a;i>b;i--)
#define upd(i,a,b) for(int i=a;i<=b;i++)
#define dwd(i,a,b) for(int i=a;i>=b;i--)
//#define local
typedef long long ll;
typedef unsigned long long ull;
const double esp = 1e-6;
const double pi = acos(-1.0);
const int INF = 0x3f3f3f3f;
const int inf = 1e9;
using namespace std;
ll read()
{
    char ch = getchar(); ll x = 0, f = 1;
    while (ch<'0' || ch>'9') { if (ch == '-')f = -1; ch = getchar(); }
    while (ch >= '0' && ch <= '9') { x = x * 10 + ch - '0'; ch = getchar(); }
    return x * f;
}
typedef pair<int, int> pir;
#define lson l,mid,root<<1
#define rson mid+1,r,root<<1|1
#define lrt root<<1
#define rrt root<<1|1
const int N = 105;
int trie[N][5];
int fail[N];
int flag[N];
int vis[N];
int cnt = 0;
ll n, m;
char buf[N];
int calc(char c)
{
    if (c == 'A')return 0;
    if (c == 'C')return 1;
    if (c == 'T')return 2;
    if (c == 'G')return 3;
}
void insert(char *s)
{
    int len = strlen(s);
    int root = 0;
    up(i, 0, len)
    {
        int c = calc(s[i]);
        int v = trie[root][c];
        if (!v) {
            trie[root][c] = ++cnt;
        }
        root = trie[root][c];
    }
    vis[root] = 1;
    flag[root] = 1;
}
void getfail()
{
    fail[0] = 0;
    int root = 0;
    queue<int>que;
    up(i, 0, 4)
    {
        if (trie[root][i])
        {
            fail[trie[root][i]] = 0;
            que.push(trie[root][i]);
        }
    }
    while (!que.empty())
    {
        int top = que.front();
        que.pop();
        int fail_ = fail[top];
        if (vis[fail_])vis[top] = 1;
        up(i, 0, 4)
        {
            int v = trie[top][i];
            if (!v)
            {
                trie[top][i] = trie[fail_][i];
                continue;
            }
            fail[v] = trie[fail_][i];
            if (vis[top]&&v)vis[v]=1;
            que.push(v);
        }
    }
}
const int mod = 100000;
struct matrix{
    int len;
    ll mat[N][N];
    void init()
    {
        this->len = cnt+1;
        memset(mat, 0, sizeof(mat));
    }
    void ones()
    {
        up(i, 0, len)
        {
            mat[i][i] = 1;
        }
    }
    matrix operator *( matrix &temp) {
        matrix res;
        res.init();
        up(i, 0, len)
        {
            up(k, 0, len)
            {
                if (!mat[i][k])continue;
                up(j, 0, len)
                {
                    res.mat[i][j] = (res.mat[i][j]+mat[i][k] * temp.mat[k][j]%mod) % mod;
                }
            }
        }
        return res;
    }
    void quick_pow(ll k)
    {
        matrix now = (*this);
        matrix res;
        res.init(); res.ones();
        while (k)
        {
            if (k & 1ll)res = res * now;
            now = now * now;
            k >>= 1ll;
        }
        (*this) = res;
    }
    void print()
    {
        up(i, 0, len)
        {
            up(j, 0, len)cout << mat[i][j] << " ";
            cout << endl;
        }
    }
}adj;
void getmatrix()
{
    upd(i, 0, cnt)
    {
        up(j, 0, 4)
        {
            int v = trie[i][j];
            if (!vis[i] && !vis[v])
            {
                adj.mat[i][v] ++;
            }
        }
    }
}
int main()
{
    m = read(), n = read();
    up(i, 0, m)
    {
        scanf("%s", buf);
        insert(buf);
    }
    getfail();
    adj.init();
    getmatrix();
    adj.quick_pow(n);
//  adj.print();
    ll ans = 0;
    upd(i, 0, cnt)
    {
        ans += adj.mat[0][i];
        ans %= mod;
    }
    printf("%lld
", ans);
    return 0;
}