从poj2778到hdu2243 AC自动机+矩阵快速幂+等比矩阵求和

于 2019-04-08 19:52:14 发布
阅读量174 收藏 1
点赞数
分类专栏： acm
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.youkuaiyun.com/qq_35885746/article/details/89103450
版权
本文介绍如何使用AC自动机结合矩阵快速幂解决两类字符串匹配问题：一是求不含特定子串的DNA序列数量；二是求包含至少一个给定子串的所有可能字符串的数量。文章详细解释了AC自动机构建过程及矩阵快速幂的应用。

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

poj2778：给你m个病毒串，问不包含病毒串的长度n的DNA片段有几个。

显然是一道AC自动机的题目，但难点在于如何求出答案。

从根出发，在AC自动机上跑，经过k次转移到达某个结点，这个结点所代表的病毒串可以看作长度为k的字符串的后缀，顺着字典树往下跑可以转移到新的长k+1字符串。

在建自动机过程中，把每个病毒串的结尾打上标记，由于AC自动机后缀节点的特性，如果某个节点的后缀节点被打上标记，则此节点也要被打上标记。

题目的数据范围是10个长度10的病毒串，所以Trie树中最多有101个结点，那么AC自动机整个转移就可以构建一个101*101的邻接矩阵，矩阵i行j列的权值是结点i转移到结点j的方案数。

学过离散数学的同学都知道，进行k次转移，从结点i转移到结点j的方案数就是这个矩阵k次幂后点（i，j）的值。

这样一来n很大的问题就可以通过矩阵快速幂来解决。

代码如下：

#include<cstdio>
#include<queue>
#include<algorithm>
#include<cstring>
using namespace std;
const int maxn=105;
const int mod=100000;
int id[130];
struct AC
{
    int next[maxn][4],fail[maxn],flag[maxn];
    int root,L;
    int newnode()
    {
        for(int i=0;i<4;i++)
            next[L][i]=-1;
        flag[L++]=0;
        return L-1;
    }
    void init()
    {
        L=0;
        root=newnode();
    }
    void insert(char buf[])
    {
        int len=strlen(buf);
        int now=root;
        for(int i=0;i<len;i++)
        {
            if(next[now][id[buf[i]]]==-1)
                next[now][id[buf[i]]]=newnode();
            now=next[now][id[buf[i]]];
        }
        flag[now]=1;
    }
    void bfs()
    {
        queue<int> q;
        fail[root]=root;
        for(int i=0;i<4;i++)
        {
            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();
            flag[now]|=flag[fail[now]];
            for(int i=0;i<4;i++)
            {
                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]);
                }
            }
        }
    }
}ac;
struct Mat
{
    long long m[105][105];
    void init()
    {
        memset(m,0,sizeof m);
    }
}ans,base;
Mat mul(Mat a,Mat b)
{
    Mat res;
    res.init();
    for(int i=0;i<ac.L;i++)
        for(int j=0;j<ac.L;j++)
            for(int k=0;k<ac.L;k++)
            {
                res.m[i][j]+=a.m[i][k]*b.m[k][j];
                res.m[i][j]%=mod;
            }
    return res;
}
char s[12];
int main()
{
    id['A']=0,id['C']=1,id['T']=2,id['G']=3;
    int m;
    long long n;
    while(scanf("%d%lld",&m,&n)!=EOF)
    {
        ac.init();
        for(int i=1;i<=m;i++)
        {
            scanf("%s",s);
            ac.insert(s);
        }
        ac.bfs();
        ans.init(),base.init();
        for(int i=0;i<ac.L;i++)
            ans.m[i][i]=1;
        for(int i=0;i<ac.L;i++)
        {
            for(int j=0;j<4;j++)
            {
                int now=ac.next[i][j];
                if(!ac.flag[now])
                    base.m[i][now]++;
            }
        }
        while(n)
        {
            if(n&1)
                ans=mul(ans,base);
            base=mul(base,base);
            n>>=1;
        }
        long long res=0;
        for(int i=0;i<ac.L;i++)
        {
            res+=ans.m[0][i];
            res%=mod;
        }
        printf("%lld\n",res);
    }
    return 0;
}

 

 

hdu2243：给你m个字符串，问长度L以内的所有字符串中，有多少个字符串起码包含一个给定的字符串。

在解决了poj2778以后，这道题的答案就是26^1+26^2+26^3+...+26^n-(A^1+A^2+A^3+...+A^n)，A为AC自动机的所有转移所构成的邻接矩阵。

剩下的问题就是如何计算等比数列求和以及等比矩阵求和。

对于等比数列求和，有:

S[n]=a+a^2+a^3+...a^n

(n==0)S[n]=1

当n为偶数时 S[n]=( a^(n/2)+1 )*S[n/2]

当n为奇数时 S[n]=( a^(n/2+1)+1 )*S[n/2]+( a^(n/2+1) )

类似的，等比矩阵求和也可以套用此方法，只需把a换成A矩阵，1换成E单位矩阵即可。

代码如下：

#include<cstdio>
#include<queue>
#include<algorithm>
#include<cstring>
using namespace std;
const int maxn=35;
struct AC
{
    int next[maxn][26],fail[maxn],flag[maxn];
    int root,L;
    int newnode()
    {
        for(int i=0;i<26;i++)
            next[L][i]=-1;
        flag[L++]=0;
        return L-1;
    }
    void init()
    {
        L=0;
        root=newnode();
    }
    void insert(char buf[])
    {
        int len=strlen(buf);
        int now=root;
        for(int i=0;i<len;i++)
        {
            if(next[now][buf[i]-'a']==-1)
                next[now][buf[i]-'a']=newnode();
            now=next[now][buf[i]-'a'];
        }
        flag[now]=1;
    }
    void bfs()
    {
        queue<int> q;
        fail[root]=root;
        for(int i=0;i<26;i++)
        {
            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();
            flag[now]|=flag[fail[now]];
            for(int i=0;i<26;i++)
            {
                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]);
                }
            }
        }
    }
}ac;
struct Mat
{
    unsigned long long m[35][35];
    void init()
    {
        memset(m,0,sizeof m);
    }
}E,x;
Mat mul(Mat a,Mat b)
{
    Mat res;
    res.init();
    for(int i=0;i<ac.L;i++)
        for(int j=0;j<ac.L;j++)
            for(int k=0;k<ac.L;k++)
            {
                res.m[i][j]+=a.m[i][k]*b.m[k][j];
            }
    return res;
}
Mat add(Mat a,Mat b)
{
    for(int i=0;i<ac.L;i++)
    {
        for(int j=0;j<ac.L;j++)
        {
            a.m[i][j]+=b.m[i][j];
        }
    }
    return a;
}
Mat mpow(Mat a,unsigned long long b)
{
    Mat res,base=a;
    res.init();
    for(int i=0;i<ac.L;i++)
        res.m[i][i]=1;
    while(b)
    {
        if(b&1)
            res=mul(res,base);
        base=mul(base,base);
        b>>=1;
    }
    return res;
}
unsigned long long mpow(unsigned long long a,unsigned long long b)
{
    unsigned long long res=1,base=a;
    while(b)
    {
        if(b&1)
            res=res*base;
        base=base*base;
        b>>=1;
    }
    return res;
}
Mat sum(unsigned long long k)
{
    if(k==1) return x;
    Mat t=sum(k/2);
    Mat res;
    if(k&1)
    {
        Mat cur=mpow(x,k/2+1);
        res=add(mul(add(E,cur),t),cur);
    }
    else
    {
        Mat cur=mpow(x,k/2);
        res=mul(add(E,cur),t);
    }
    return res;
}
unsigned long long sum(unsigned long long a,unsigned long long k)
{
    if(k==1)
        return a;
    unsigned long long t=sum(a,k/2);
    if(k&1)
    {
        unsigned long long cur=mpow(a,k/2+1);
        t+=cur*t;
        t+=cur;
    }
    else
    {
        unsigned long long cur=mpow(a,k/2);
        t+=cur*t;
    }
    return t;
}
char s[12];
int main()
{
    int m;
    long long n;
    while(scanf("%d%lld",&m,&n)!=EOF)
    {
        ac.init();
        for(int i=1;i<=m;i++)
        {
            scanf("%s",s);
            ac.insert(s);
        }
        ac.bfs();
        E.init(),x.init();
        for(int i=0;i<ac.L;i++)
            E.m[i][i]=1;
        for(int i=0;i<ac.L;i++)
        {
            for(int j=0;j<26;j++)
            {
                int now=ac.next[i][j];
                if(!ac.flag[now])
                    x.m[i][now]++;
            }
        }
        x=sum(n);
        unsigned long long ans=0;
        for(int i=0;i<ac.L;i++)
        {
            ans+=x.m[0][i];
        }
        ans=sum(26ull,n)-ans;
        printf("%llu\n",ans);
    }
    return 0;
}

 

评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

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

抵扣说明：

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

余额充值