poj2778 DNA Sequence （AC自动机+矩阵乘法）

poj2778 DNA Sequence

题目大意：一段DNA序列里面如果带有某些特征序列，那么这段DNA序列就有问题。问给出若干段DNA序列，再给定一个长度，问这个长度的所有可能序列中，有多少是没问题的。
嗯，总之先建机行事。
先建好Trie树，升级成AC自动机。由于m很小，每个有病序列的长度也很小，但是n很大，所以想到用快速幂优化，邻接矩阵表示各个状态之间的关系。
假设矩阵An×n中的第i行，第j列元素Aij表示从i状态走1步到达j状态有几种方法，Bn×n矩阵中的第j行，第k列元素 Bkj表示从k状态走1步到达j状态有几种方法，那么矩阵 Cn×n=A×B中的 Cij=∑nk=1(Aik×Bkj) 即表示从i走到j走两步有几种方法。

通过建好的自动机，可以找到每个状态走一步可以到达的下一个状态，由状态之间的关系可以建立一个邻接矩阵，其中第i行第j列表示从状态i走一步到状态j有几种方法，然后这个矩阵的n次方就是每个状态在经过n步之后的结果，其中第一行的总和就是从树根到达其它状态的n步之后一共有几种方式。在跑自动机的时候，需要筛掉那些到达有病节点的状态。

随手写个样例：
3 4
ACC
AGT
GT
邻接矩阵应该长这样（可能行列位置会有变化）：
2 1 0 0 0 0 1 0
1 1 1 0 1 0 0 0
1 1 0 0 0 0 1 0
0 0 0 0 0 0 0 0
1 1 0 0 0 0 1 0
0 0 0 0 0 0 0 0
1 1 0 0 0 0 1 0
0 0 0 0 0 0 0 0
结果应该是201

另外（为了偷懒？？），用了一点trivial的trick，就是ACTG的ASCII码的二进制表示的倒数第第二个和第三个二进制数刚好是从0~3的（A00 C01 T10 G11），所以可以直接将字符右移一位然后&3就可以完成从A到G的字符分选了。
代码比较冗余，也没有怎么优化过，各种手动开空间，一位数组当二维数组做，乘法找地址，跑起来速度一般（172ms）Orz

#include<iostream>#include<cstdio>#include<cstring>#include<string>#include<queue>#include<cstdlib>using namespace std;struct node{    node *fail;    node *son[4];    int index;    bool vis1, vis2, isdest;};char s[20];node *root;int index = 0;long long *matrix, Size = 1, bytesize = 0, qwordsize = 0;void multiply(long long *a, long long *b, long long *save){    memset(save, 0, bytesize);    for(int i = 0; i < Size; i++)        for (int j = 0; j < Size; j++)        {            long long *aa = a + i*Size, *bb = b + j, *sv = save + i*Size + j;            for (int k = 0; k < Size; k++)                *sv += (*(aa + k) )*( *(bb + k*Size));            *sv %= 100000;        }    memcpy(a, save, bytesize);}void f_pow(long long a[], int b){    long long ans = 0;    long long *base = new long long[qwordsize];    long long *save = new long long[qwordsize];    memset(base, 0, bytesize);    for (int i = 0; i < Size; i++)        base[i*Size + i] = 1;    while (b > 0)    {        if (b & 1)            multiply(base, a, save);        multiply(a, a, save);        b >>= 1;    }    for (int i = 0; i < Size; i++)        ans += base[i];    printf("%d\n", ans % 100000);}void buildtrie(char *s, node *root){    int len = strlen(s), i;    node *temp = root, **son;    for (i = 0; i < len - 1; i++)    {        int ch = (s[i] >> 1) & 0x3;        son = &temp->son[ch];        if (*son == NULL)        {            *son = new node;            memset(*son, 0, sizeof(node));            Size++;            index++;            (*son)->index = index;        }        temp = *son;    }    int ch = (s[i] >> 1) & 0x3;    son = &temp->son[ch];    if (*son == NULL)    {        *son = new node;        memset(*son, 0, sizeof(node));        Size++;        index++;        (*son)->index = index;    }    (*son)->isdest = true;}void buildfail(node *root){    queue<node*> que;    for (int i = 0; i < 4; i++)    {        node *temp = root->son[i];        if (temp != NULL)        {            que.push(temp);            temp->fail = root;        }    }    while (!que.empty())    {        node *temp = que.front();        node **son = temp->son;        que.pop();        for (int i = 0; i < 4; i++)        {            if (son[i] != NULL)            {                que.push(son[i]);                node *t = temp;                do                {                    t = t->fail;                    if (t->son[i] != NULL)                        break;                } while (t != root);                if (t->son[i] != NULL)                    son[i]->fail = t->son[i];                else                    son[i]->fail = root;            }        }    }}bool finddest(node *temp){    if (temp->vis1)        return temp->isdest;    if (temp->isdest)        return true;    if (temp != root)        temp->isdest = finddest(temp->fail);    temp->vis1 = true;    return temp->isdest;}void ac_automation(node *root){    node *temp = root;    queue<node*> que;    que.push(root);    root->vis2 = true;    while (!que.empty())    {        node *temp = que.front();        for (int i = 0; i < 4; i++)        {            node *t = temp;            while (t != root)            {                if (t->son[i] != NULL)                    break;                t = t->fail;            }            if (t->son[i] != NULL)            {                if (!finddest(t->son[i]))                {                    if (!t->son[i]->vis2)                    {                        que.push(t->son[i]);                        t->son[i]->vis2 = true;                    }                    matrix[temp->index*Size+t->son[i]->index]++;                }            }            else            {                if(root->son[i]!=NULL)                    matrix[temp->index*Size+root->son[i]->index]++;                else                    matrix[temp->index*Size]++;            }        }        que.pop();    }}int main(){    int m, n;    scanf("%d%d", &m, &n);    root = new node;    memset(root, 0, sizeof(node));    while (m--)    {        scanf("%s", &s);        buildtrie(s, root);    }    matrix = new long long[Size*Size];    qwordsize = Size*Size;    bytesize = qwordsize << 3;    memset(matrix, 0, bytesize);    buildfail(root);    ac_automation(root);    f_pow(matrix, n);}