AC自动机小试 POJ2778

http://poj.org/problem?id=2778

题目翻译：给定m个异常基因片段，在所有长度为N的DNA序列中如‘AAAAA....’,'ATTTTT....','ATCG...','CGTA...',.....中不包含异常序列的DNA个数mod 100000(十万)的值

这题题解是AC自动机和矩阵乘法，首先建立AC自动机

我们先看AC自动机建立过程(会的跳过，不会的英语好的看这篇英文论文：http://wenku.baidu.com/view/8b6f6d85ec3a87c24028c467.html，不会的英文不好的看不懂论文的往下看)

首先定义node为记录类型

{   son[1..4]of int;code:char;fail:longint;danger:boolean;

    son,fail初始为-1或nil或null,danger:=false;

     }

定义函数

function anti(a:char):longint;

begin  case a of

            'A':exit(1);
        'C':exit(2);
        'G':exit(3);
        'T':exit(4);
        end;

end;

1.建立Trie树，即字母树

2.初始化第一层节点fail指针为root节点，将第一层节点加入队列中

3.取队首元素x，将所有x的son入队

4.for y ∈ x.son

       begin

            temp:=x.fail;

            while (temp>0)or(temp<>root) do

                 begin

                     for z ∈ temp.son

                        if z.code=y.code

                           then flag:=true,break;

                     if flag then temp:=z,break;

                     temp:=temp.fail;

                 end;

           if  temp=root

             then

                  for z ∈ root.son

                        if z.code=y.code

                           then temp:=z,break;

            y.fail:=temp;

       end;

5.Complete Build AC Automaton

注意到一个经典的矩阵乘法题，给定一个有向图，问从A点恰好走k步（允许重复经过边）到达B点的方案数mod p的值

答案是将有向图的邻接矩阵乘k次输出G[A，B]

我们可以把A点看做root点，走n步不经过异常节点的方案数mod p(p=100000)的值

那么很明显一点是x节点的son中的异常节点(设为y)是不可能由x转移到的

题中每个节点最多有4个son

如果某个x.son为nil呢

设z=x.fail,当前son为i

由于x.fail为x.son匹配失败时的最长后缀在trie中的位置，那么以fail指针往上(fail指针总是指向上面)寻找z.son[i]不为nil的z节点后终止

然后判断z.son[i]是否为异常节点，不是则g[x,z.son[i]]++或者若z=root且root.son[i]=nil 则g[x,root]++;

按照trie的bfs序列

取当前元素x

for i:=1 to 4 do

  if x.son[i]=-1

    then

      begin

            z:=x.fail

            while (z<>root)and(z.son[i]=-1) do

                 z:=z.fail;

            if  (z=root)and(root.son[i]=-1)

              then g[x,root]++

             else

                if   (z.son[i]<>-1)and(not z.son[i].danger)

                  then g[x,z.son[i]]++;

      end

   else if  not x.son[i].danger then g[x,x.son[i]]++;        

然后就是将g矩阵乘n次，将第root行加起来mod 100000 就是ans

代码pascal

program poj2778;type    matrix=array[0..101,0..101]of int64;const   mo=100000;var     a:string;nodes,leng,i,j,point,m,n,top,x,tail,som,k:longint;ans:int64;        b,g:matrix;        son:array[0..101,0..4]of longint;        off:array[0..101]of longint;        code:array[0..101]of char;        stack:array[0..110]of longint;        comp:array[0..110]of boolean;procedure muti(a,b:matrix;var c:matrix);begin   fillchar(c,sizeof(c),0);        for i:=0 to nodes do          for j:=0 to nodes do            begin              for k:=0 to nodes do                c[i,j]:=(a[i,k]*b[k,j]+c[i,j])mod mo;            end;end;function anticode(a:char):longint;begin   case a of        'A':exit(1);        'C':exit(2);        'G':exit(3);        'T':exit(4);        else exit(0);        end;end;begin        readln(m,n);        fillchar(son,sizeof(son),$ff);        for i:=1 to m do          begin            readln(a);            leng:=length(a);            point:=0;            for j:=1 to leng do              begin                if son[point,anticode(a[j])]<>-1                  then                    point:=son[point,anticode(a[j])]                  else                    begin                      inc(nodes);                      son[point,anticode(a[j])]:=nodes;                      code[nodes]:=a[j];                      point:=nodes;                    end;              end;            comp[point]:=true;          end;        top:=0;        tail:=1;        while top<tail do          begin            inc(top);            x:=stack[top];            for i:=1 to 4 do              if son[x,i]<>-1                then                  begin                    inc(tail);                    stack[tail]:=son[x,i];                    if x=0 then continue;                    som:=off[x];                    repeat                      if son[som,i]<>-1                        then                          begin                            som:=son[som,i];                            break;                          end;                      som:=off[som];                    until som<=0;                    if som=0                      then                        if son[som,i]<>-1                          then                            som:=son[som,i];                    if comp[som]                      then comp[son[x,i]]:=true;                    off[son[x,i]]:=som;                  end;          end;        top:=0;nodes:=0;        while top<tail do          begin            inc(top);            x:=stack[top];            if comp[x] then continue;            if x>nodes then nodes:=x;            for i:=1 to 4 do              if son[x,i]=-1                then                  begin                    som:=off[x];                    while (som>0)and(son[som,i]=-1) do                      begin                        som:=off[som];                      end;                    if (son[som,i]<>-1)and(not comp[son[som,i]])                      then                        inc(g[x,son[som,i]])                      else                        if (som=0)and(son[som,i]=-1)                          then                            inc(g[x,0]);                  end                else                  if (son[x,i]<>-1)and(not comp[son[x,i]])                    then                      inc(g[x,son[x,i]]);          end;        for i:=0 to nodes do          b[i,i]:=1;        while (n>0) do          begin            if (n and 1)>0              then                muti(b,g,b);            n:=n>>1;            muti(g,g,g);          end;        for i:=0 to nodes do          ans:=(ans+b[0,i])mod mo;        writeln(ans);end.

另一篇题解：http://hi.baidu.com/%D2%D5%C1%D6010/blog/item/6db06ccf0a3b440993457e7b.html