LCS问题的另一种解法

找LCS长度的算法O(nLog n)(AT&T, Unix Group)

       定义:Place(a)={ i |xi=a}
       例:X=<A,B,C,B,D,A,B>
        place(A)={1,6}, place(B)={2,4,7}, place(C)={3}, place(D)={5}，place桶结构如下:

       扫描数组X（序列）一次即可建成桶place，故只要O(n)时间。但请注意桶中每个槽的链表，数值小的排在后面。再定义S(j,k)={ i |若x1…xi与y1…yj的LCS长度为k}

       例：X=<A,B,C,B,D,A,B>，Y=<B,D,C,A,B,A>，当j=4时Y4=BDCA，

              S(4,0)={0}（f与BDCA），

              S(4,1)={1,2}（A、AB与BDCA），
              S(4,2)={3,4,5}（ABC、ABCB、ABCBD与BDCA），
              S(4,3)={6,7} （ABCBDA、ABCBDAB与BDCA），
              S(4,4)=NULL。（这里可以直接跳到最后看该问题的一步一步的运算过程）
       因此对于给定的j，每个S(j,k)恰好产生一个{0,1,2,…n}一个划分(0≤k≤j)。如何由S(j-1,k)生成S(j,k)？(例：根据S(4,k)生成S(5,k))
       设r∈plase(yj)，即有xr=yj。
       例如place(y5)={7,4,2}，有x2=x4=x7=y5=B。如果把yj加到y1…yj-1上(注意xr=yj)，y1…yj-1yj与x1…xr-1xr的LCS，是否会比y1…yj-1与x1…xr-1xr的LCS长度增加?
      CASE1:
      在第j-1轮（上一轮），r-1和r同属于某个集合 。（例：r=7, r-1=6, 有6,7∈S(4,3) ，j-1=4, k=3）

      按定义, x1…xr-1与y1…yj-1的LCS长度为k，x1…xr与y1…yj-1的LCS长度也为k，例：ABCBDA与BDCA、ABCBDAB与BDCA的LCS长度均为3。这意味着：对y1…yj-1，把xr (例子中的B) 加到xr-1之后，LCS并未增长，（即在第j-1轮加上xr对增长LCS没有起到作用）。因此在第j轮时，yj(=xr)加到yj-1之后，y1…yj-1yj与x1…xr-1xr的LCS，就可以比y1…yj-1与x1…xr-1xr的LCS长度增1。e.g. 加上yj后，ABCBDAB与BDCAB的LCS长度变为4。因此，此时应将是S(j-1,k)中所有大于等于r的元素都移到S(j,k)中。

      例：在j=5轮初始时，查表可得place(y5)=place(B) ={7,4,2}。先取r=7，有6,7∈S(4,3) 即r-1和r同属于一个集合，按上述讨论，此时应将 中所有大于等于r的元素都移到 中，于是根据 S(4,3)={6,7} -->S(5,3)={6}，S(5,4)={7}∪S(4,4) 由于S(4,4)=NULL（见上页），故S(5,4)={7}。 
      其次取r=4，有3,4 ∈S(4,2)，亦有r-1和r同属于一个集合，于是根据S(4,2)={3,4,5} --> S(5,2)={3}，S(5,3)={4,5}∪S(5,3)={4,5,6}。最后取r=2，有1,2 ∈S(4,1)，仍然有r-1和r同属于一个集合，于是根据S(4,1)={1,2} --> S(5,1)={1}， S(5,2)={2} ∪S(5,2)={2,3}。

      另外，若plase(yj)有多个r（如上例），则执行时必须按r的值从大向小做，否则出错。

      CASE2:
      在第j-1轮，r-1和r不属于同一个集合，即有r∈S(j-1,k) 而r-1∈S(j-1,k-1) ,说明x1…xr-1与y1…yj-1的LCS长度为k-1，而x1…xr与y1…yj-1的LCS长度为k。即在第j-1轮中，xr已经匹配上，故而在第j轮中，xr不能与yj再次匹配，亦即S(j,k)无须变化。另外, 中的j只是为了逻辑上的理解而引入的，实际执行时，S(j-1,k)与S(j,k)使用同一物理空间（均为Sk）。

       初始时S0={0,1,2,3,…,n}，S1, S2,…Sm均为空。为了使上述的集合分裂和集合保序合并易于进行，Sk中的元素以树形式组织，组织形式如下：

       初始时所有元素i均为叶结点,除S0以外，其他均为“NULL”，将所有节点都加入S0中去，完成初始化。
       算法:

       把元素i (i =1,2,….n) 依次插入place桶中第 H(X[i])个链表的最前部；（O(n)）H(X[i])= X[i]的ASCII编码 - 字符中最小的ASCII编码。

       for j=1 to m do /*对yj计算Sk*/
           {   a←place[H(y[j])] /*桶的place(yj)槽中的首个（值最大）元素*/
               while a≠NULL do/*对place(yj)槽中的所有r从大到小做*/
                 {   r←a->dat; /*r∈plase(yj)*/
                     k←Find(r); /*找到集合Sk*/
                     if k=Find(r-1) then /* CASE1:r-1和r同属于一个集合*/
                        {  split_and_union(k,r,s,l)；/*元素 < r者留在Sk，>= r者放入Sk+1*/
                         }

                      /*若r-1和r不属于同一个集合，则Sk不改变*/
                      a←next[a]

                   }
            }
         k←0; while s[k]≠^ do k←k+1; /*找下标最大的非空集合Sk*/
         LCS-length←k-1; return LCS-length /*找到LCS长度，算法结束*/
       算法分析:
       最坏情况：X=<AAA…A>，Y=<A…A>，此时p=n*m，总时间为O(n*m* Log n)，坏于动态规划法；
       但AT&T的Unix Group证明了该方法的平均时间复杂度为O(n*Log n)，要好于动态规划法。
       算法执行过程：

        X=<A,B,C,B,D,A,B>，Y=<B,D,C,A,B,A>。
        在初始时，S0={0,1,2,3,4,5,6,7}，S1, S2,…Sm均为空。
        在j=1轮，查表可得place(y1)=place(B) ={7,4,2}。
                       取r=7，得S0={0,1,2,3,4,5,6}，S1={7}；
                       取r=4，得S0={0,1,2,3}，S1={4,5,6,7}；
                       取r=2，得S0={0,1 }，S1={2,3,4,5,6,7}；
        在j=2轮，查表可得place(y2)=place(D) ={5}。
                       查r=5，得S0={0,1 }，S1={2,3,4}，S2={5,6,7}。
        在j=3轮，查表可得place(y3)=place(C) ={3}。
                       查r=3，得S0={0,1 }，S1={2}，S2={3,4,5,6,7}。
        在j=4轮，S0={0}，S1={1,2}，S2={3,4,5}。
        在j=5轮，按前述 有S0={0}，S1={1}，S2={2,3}，S3={4,5}，S4={6,7}。
        在j=6轮，S0={0}，S1={1}，S2={2,3}，S3={4,5}，S4={6,7}。

C语言实现：

#include <stdio.h>
#include <stdlib.h>
#include <string.h>

#define		STR_LEN		512
#define		CHAR_LEN 	26
#define		H(x)        ((x)-'A')

typedef struct place_node{
	int dat;
	struct place_node* next;
}place_node,*place_node_ptr;

typedef struct place_head{
	place_node_ptr first;
}place_head,*place_head_ptr;

typedef struct {
	int father;
	int next;
}leaf_node,*leaf_ptr;

typedef struct set_node{
	int first_leaf;
}set_node,*set_ptr;

void init_place(place_head_ptr place,char *str,int str_len)
{
	int i;place_node_ptr p;
	for(i=0;i<CHAR_LEN;i++)
		(place+i)->first=NULL;
	for(i=0;i<str_len;i++){
		p=(place_node_ptr)malloc(sizeof(place_node));
		p->dat=i+1;
		p->next=(place+H(str[i]))->first;
		(place+H(str[i]))->first=p;
	}
}

void init_set_leaf(set_ptr s,int s_len,leaf_ptr l,int l_len)
{
	int i;
	for(i=0;i<l_len;i++){
		(l+i)->father=0;
		(l+i)->next=i+1;
	}

	for(i=0;i<s_len;i++)
		(s+i)->first_leaf=-1;

	(l+l_len-1)->next=-1;
	s->first_leaf=0;
}

int find(int x,leaf_ptr leaf)
{
	return (leaf+x)->father;
}

void split_and_union(int k,int r,leaf_ptr leaf,set_ptr set)
{
	int i,j;
	j=(set+k+1)->first_leaf;
	(set+k+1)->first_leaf=r;
	(leaf+r-1)->next=-1;
	for(i=r;(leaf+i)->next != -1;i++)
		(leaf+i)->father=k+1;
	(leaf+i)->father=k+1;
	(leaf+i)->next=j;
}

int main()
{
	int i,j,k;
	int sa_len,sb_len,s_len,l_len;
	leaf_ptr leaf;set_ptr set;
	char str_a[STR_LEN],str_b[STR_LEN];
	place_head_ptr place;place_node_ptr p;
	scanf("%s%s",str_a,str_b);

	sa_len=strlen(str_a);
	sb_len=strlen(str_b);

	place=(place_head_ptr)malloc(CHAR_LEN*sizeof(place_head));
	init_place(place,str_a,sa_len);

	l_len=sa_len+1;s_len=(sa_len > sb_len?sb_len:sa_len)+1;
	leaf=(leaf_ptr)malloc(l_len*sizeof(leaf_node));
	set=(set_ptr)malloc(s_len*sizeof(set_node));

	init_set_leaf(set,s_len,leaf,l_len);

	for(i=1;i<=sb_len;i++)
	{
		for(p=(place+H(str_b[i-1]))->first;p;p=p->next)
		{
			j=p->dat;
			k=find(j,leaf);
			if(k == find(j-1,leaf)){//in the same set
				split_and_union(k,j,leaf,set);
			}
		}
	}

	for(k=0;(set+k)->first_leaf != -1;k++);
	printf("lcs'len = %d\n",k-1);

	return 0;
}

        that'all~~~