Category Archives: 算法学习

打个招phd广告,短期有效

Software Engineering PhD or RA position at HKUST

I am looking for a candidate for the position of either the research assistant or the PhD graduate student in the area of software engineering or programming languages. The student must be strong in either algorithms or math and passionate about programming. The candidate must also be self-motivating, curious, and independent.

The salary of RA is around 9000 HKD per month, and graduate student will be awarded a full scholarship at around 13K per month.

Interested students are welcomed to visit my website: www.cse.ust.hk/~charlesz for more details of the research projects and to send his/her resumes and other supporting documents to [email protected]

几个月的小总结

有趣的是,我还是没能遵守诺言,再忙也抽点空来写写blog。一晃四个月过去了,也学习了一些有趣的东西,又到了该总结一番的时候。
遗憾的是,本来想做一个纯粹的软工research,最终没能实现。这其中有值得总结的地方,首当其冲的问题便是我先入为主的认为一个问题它的复杂度高,所以我只要能把降下来那就可以发paper。其次,因为我的懒惰,所以惧怕去做一些不太兴趣的事情。所以,当我好不容易把一个很实际的问题提炼成图论模型时,我自然不想松手。
事实证明,其实naïve的方法,才是最好的方法:又好理解又好写。根据Occam剃刀原理,其它方法应该被排除了。可是,刚开始,我并没意识到这一点,因为我没有料到一个非常普通的数据结构: sparse bitmap,在实际中性能表现那么好。我承认我起初并非不知道这个结构,而是出于写起来比较麻烦,而没有去做。后来发现GCC里面实现了,我才赶忙拿出来尝试。记得第一次实验完的时候,我心是完全凉了,知道它会好没想到它会这么好。我心里一横,想反正没做多久,干脆抛弃了赶紧找新topic,离deadline还有一个多月,或许还来得及。
就这样,我悄无声息的开始尝试弄一些新的东西出来。后来我有了些小进展,跟老板说,他也觉得奇怪,怎么我突然就不做以前那个了?我没敢说我的实验,因为当初就是因为老板对我的信任,才让我尝试的,我不能告诉他失败了,我只能说contribution不够大,希望做一些results加强,这样更有希望。
对了,还没介绍我一直在做什么。pointer analysis,一个已经研究了四十年的东西,就是静态分析指针该指向何处,任何一个现代的编译器都会分析,只是算法不同分析结果的精度不同,因为广义上它是个undecidable问题。我选它并非是想challenge the whole world,我自知没这个能力,原因仅仅是这个问题很容易理解。参加过今年TIC决赛的同学估计还记得那道A题,便是exactly我所一直在研究的一个问题,然后出了题去恶心大牛的。虽然9月的ICSE我最终没能做出一些东西,但我肯定还会继续探索,因为我在PL差得还太多太多,做一个researcher不能操之过急,当然更不能临阵退缩。
不过,虽然高水平的研究是一个漫长的过程,但HKUST并不会因此而同情你,你得发paper,失败不得。ICSE过后,我又把我那个已经丢在垃圾桶中的算法捡起来,重新改了改并增加了一些证明,弄得还算fancy,投了出去。没报希望,但也必须尝试。
下一步,该往何处?自然,这有段空隙,是该补补基础的时候了。打算花一个月时间,把compiler中的一些基本问题再弄明白些,比如data flow analysis的lattice theory以及其应用。不能以后出去别人问起这些基本问题我心里还打鼓,讲得个是是而非的,那说我是做PL的不是个千古笑话么?
下一篇blog打算从compiler前端开始复习,第一个topic是DFA minimization和predict parser。看能不能把Hopcroft的O(nlgn)算法讲清楚。

后缀数组相关算法(三)

本来打算这篇文章就开始讲SA的应用,可我决定推迟一篇文章,主要原因是:
1. 最近读到的一篇文章中有一个想法我觉得很新颖,这对于希望从事研究类工作的人来说有一定的帮助;
2. 感觉前面的的文章很没有条理,想补充点内容把没说的说完整。

先从SA构造算法的分类开始说起。一般来说,现在有四种思路构造SA,分别如下:
1. Prefix Doubling,这个思路的代表就是MM(倍增)算法;
2. Recursive,代表是DC3算法;
3. Induced Copying,代表是MF算法;
4. 从后缀树构造后缀数组。

其中,Induced Copying目前公认是在实践中最快的方法,但是据我所知目前这类算法实现起来都不太优美,它们大多需要用到Sedgewick和Bentley发明的三路划分字符串快速排序算法,这样代码量较大从而在比赛中不容易推广。
本文要说的算法叫KA(Ko and Aluru),也是Recursive Algorithm中的一种,从论文中找到的数据表明要比DC3算法快且节省内存。我并没有实现该算法,并不是说实现难度很大,而是这个算法开始部分的优美让我陶醉,但后面的部分的笨拙(并非说低效)让我觉得这个算法没能一气呵成,所以我暂时不想去写。
1. 基本算法
首先约定,原串是T。对于T的每个后缀T(i),把它划分为两个集合:T(i) > T(i+1)的为L集,T(i) < T(i+1)的为S集,其中T(n)属于S。对T按照单个字符基数排序形成数组A,这样T就以字符集E为基础被分成了很多个Bucket,设head(a)表示字母a的Bucket在A中的起始序号,tail(a)表示字母a的Bucket在A中结束序号,如下图所示:
sa3_1
然后,对S集的所有后缀排序,形成数组B(如何做这个后面说明)。
得到了B,从右向左一次扫描。对每一个i,将suf=B[i]放到A中tail[ suf[0] ]的位置,然后tail[ suf[0] ]减少1。这个步骤没什么好说的,实际上就是把B的结果归位到A中。
神奇的是后面。下面先来看一个定理:
定理1.对于两个相同字母开头的后缀suf1和suf2,如果分别来自S集和T集,那么suf1一定最终排在suf2的后面。
其实很容易证明,suf1[0] == suf2[0],又suf1[0] < suf1[1],suf2[0] > suf2[1],所以suf1[1] > suf2[1],得suf1 > suf2,得证。
于是,利用这个定理,对B数组从左至右扫描,对遇到的每个后缀suf,把该suf放到A数组的tail[ suf[0] ]位置上,并且tail[ suf[0] ]减少1。
后面的思路和DC3类似,要做到O(n)的复杂度,就必须使用S集的排序导出L集的排序,下面首先来看一个过程:
从A数组开始从左至右扫描,对遇到的每个后缀suf,设R[suf]表示这个后缀在T中的起始位置。如果R[suf] – 1是一个T集的后缀,设为suf2,则将其放入head[ suf2[0] ]的位置,并且head[ suf2[0] ]增加1,这个步骤叫做suf1确定suf2。
这个过程一举两得,不但有了L的排序,同时也将S和L归并完成,是KA算法最为让人拍手叫绝的一步(读到此至少我是眼前一亮)。
定理2.以上算法得到的A是T的后缀数组。
要证明它的正确性其实很简单,如下:
1. A的第一个元素肯定是一个S集的后缀。如果是L集的,设为suf1,那么必然存在一个后缀suf2小于suf1。如果suf2仍然是L集的,那么必然有suf3……这个过程不可能永远进行下去,所以必然最小的后缀是来自S集的(是不是感觉和定理1冲突了?);
2. 对于L中每一个后缀suf,如果要被放到A中,必然是R[suf]+1这个后缀已经被放入A。而R[suf]+1这个后缀除非是S集元素,否则它会要求R[suf]+2的后缀先放入A,直到遇到一个S集的元素,所以,每一个元素都是被A中排在其之前的元素确定的。所以可知算法进行到A[i]时,A[i]必然已赋值;
3. 其次,每一个A[i]的赋值必然是正确的,这一步可以通过数学归纳法来证明。已知第1位A是正确赋值的,设前i位A都能正确赋值,那么第i+1位如果错误,可知必然存在j > i + 1,A[j]应该填充A[i+1]。A[j]必然和A[i+1]的首位字母是一样的(按照算法的流程,不一样是不可能争同一个位置的,对吧?),因此就有R[A[j]]+1的后缀小于R[ A[i+1] ]+1的后缀,那么根据假设和算法的流程,必然会先安排A[j]到位置i+1,因此赋值错误是不可能的。
剩下唯一的难点就是对S进行排序。
方法有二。一是和DC3类似,把所有的S按T中出现的顺序提出来,然后将连续的字符通过排序替换成一个新字符,然后递归形成后缀数组。这样做的难点在于,DC3中只需要将连续的三个字符堪称一个新符号,这里却是不规则的,因此普通的基数排序很难进行,所以算法的作者搞了一个很麻烦的方法,想了解的可以参看算法的原文。
另外一个方法就是I&T算法,是Induced Copying的一种,它直接用前文说的快速字符串排序把S搞定。很直接很暴力,但未必不高效。
最后需要注意的是,其实先将L排好序,再导出S的顺序并形成最终的结果也是可以的,步骤和上面类似,但是放入每个bucket的顺序恰好相反。所以,可以视S和L的数量来定该先排序哪个,这样可以改善效率。
2. 进一步提速
KA算法的导出和归并步骤相比DC3要高效许多,要想进一步提速,就需要减少排序S集的时间。更好的排序思路很难找到,那么就通过给S集瘦身来达到目的吧。
现在,只提取T[i]是S后缀但T[i+1]是L后缀的i形成S*集。将S*按照原有的方法排序,然后依然按原来的方法先吧S*里面的元素分配到A中。
剩下的S集元素可以通过S*的元素导出。方法是,再一次从右向左扫描S*数组,如果R[ B[i] ] – 1是S集的,那么把放到开头字母对应的bucket中最后一个位置上。
这一个改进的证明其实很简单,我不再多说。
其实,这些算法的所有证明其实都和一个性质有关:单调性。参透了这一点就算是对KA算法真正的理解了。

原文下载: KA Suffix Arrary 构造算法

后缀数组相关算法(二)

相隔半个月之久,我又看了许多关于SA的论文,但似乎给我留下深刻映像的并不多。我希望寻找简单而又巧妙的思路也未能成功(不一定复杂度很低),也或许是我已经找到了很好的思考角度但是没有很强的分析能力于是依然无获而终。
本文打算详细说一下DC3算法,早在好几个月前便答应要写。而下一篇将会说一个SA很有意思的应用,来自数据压缩领域(猜到是什么了吗?)。
DC3最早见于论文的名字叫skew算法,不久他们伙同Google的一个人一同撰写了另外一篇文章,内容上和前一篇没有很大区别,但是从更本质的角度阐述了skew算法,即difference cover原理,并提出更一般的DC(v)算法,而原来的skew也就正式更名为DC3算法,因此若只想了解DC3,看哪篇文章都一样。
这个算法和前面不同的是它利用了分治的思想,而该思路最早是由Farach提出并用于Suffix Tree的构造。设想,如果我们能分别得到奇数位置的所有后缀和偶数位置的所有后缀形成的SA,然后再进行一次归并操作,是不是就可以O(nlgn)内构造SA呢?
假设偶数位置的所有后缀形成的SA为SA2,奇数位置为SA1,相应的Rank数组为Rank1和Rank2,原字符串是S(注意,S从1开始计数)。现在SA2已经求出,那么对于所有奇数位置的后缀i,利用文章1中的方法构造二元组< S[i], Rank2[i+1]>,显然对于这个二元组的每一个元,大小比较都有良序定义,则直接基数排序之就得到SA1,避免了一次递归调用。
再来考虑SA2该如何求。既然是分治,因此只要能构造出和原问题类似的子问题即可。令S’是S从位置2开始的后缀。先走题一下,设有两个字符串A和B(为方便讨论,两个都是偶数长度),如果有A < B,那将A和B中的每两个字符看成一个字符后(即每个新的字符是一个二元组,大小比较方法是先比较第一维再比较第二维),依然有A’ < B’,即这样的映射具有保序性。同理,将S’中1、2位置压缩成一个字符,3、4位置压缩,这样S’的所有后缀恰好对应S的所有奇数位置压缩后的后缀。由于有保序性,因此求得S’的SA就得到了S的SA2。
由于每轮只递归一次,|S’| = |S| / 2,因此方程是f(n) = f(n/2) + O(n),即算法的时间复杂度是O(n),空间复杂度也是O(n)。
以上的算法很完美,但想必大家都能看出算法的描述和数字3没有任何关系,所以以上所述其实并不是DC3算法。细心的读者或许还能看出另外一个破绽,我没有讲归并如何进行。设p1是来自奇数位的后缀,p2是偶数位,则p1[1]和p2[1]比较后,剩下的部分还不能马上确定大小关系(因为又分别来自SA2和SA1)。实际上,的确是有办法可以做到线性时间内合并,但非常的不trivial,不值得推荐。
而真正的DC3和上述算法的区别就是最后一步归并变得非常简单。试想我们已经有了所有起始位置mod 3 != 0的SA12和Rank12,而SA0和Rank0是位置mod 3 = 0的后缀。归并这两者时,设p0,p1,p2分别是mod 3 = 0,1,2的后缀,那么p0和p1归并可以先比较第一个字符,然后就分别剩下起始mod 3 = 1和mod 3 = 2的后缀,而它们的大小关系通过Rank12就直接得到。同理,比较p0和p2时,先比较前两个字符,又剩下两个p2和p1的后缀,大小关系也能直接判别。
而SA12的求法和上述算法中的SA2是同理的,只不过要提取出mod 3 = 1,2的所有后缀,并且是每三个字符压缩成一个新符号。而唯一不同的地方是,SA12中要先连续放置mod 3 = 1的后缀,然后再放置mod 3 = 2的后缀,这样当符号压缩了以后才是原问题的子问题。
那difference cover究竟是什么呢?定义S是一个DC samples集合,对任意i和j属于[0,n),存在k使得i+k和j+k都属于S,比如上面的SA12。根据此理论,如果每次压缩4个字符,那么也只要选mod 4 = 2, 3的后缀作为子问题排序即可。当然,更详细的数论的结论建议大家看看论文,简单而有启发。
最后,程序的代码在原始论文的后面有,写得很巧妙,基本上没什么可以改的(我试图自己写一个新的实现,结果始终不如论文作者思路清晰),因此我也不贴这里了。

后缀数组相关算法(一)

        很久便说起要写后缀数组(简称SA)的相关内容,结果一直以来未曾有动静,原因是我既忙又懒。最近上网搜了些帖子,发现SA这东西早已经深入人心了,自然感觉自己已经开始慢慢的脱离于时代。
其实不用把SA看得有多神秘,它只不过是一个字符串的所有后缀排序后的结果。稍微有点常识的人都可以想到,把所有后缀取出来然后跑一个快排就可以拿到SA。没错,如果你认为O(n^2lgn)的复杂度你完全可以接受的话,显然这是一种最快速的编程方法。
自然,无论是竞赛中还是实际的工程项目中,上诉算法都是难以利用的,它复杂度实在太高了。我们尝试改进这个算法,联想到RK算法利用hash的思想将串匹配问题的复杂度降到了O(n),那么如果找到一个hash函数F,使得对任意的串S1和S2,其中S1 < S2,有F(S1) < F(S2),即F具备保序性,那么可以先把这些后缀都映射成整数,然后再快排或基排就可以解决问题。
嗯,这个想法到是好,但为了比原来的算法更优,映射操作的复杂度不能太高,所以这个F怎么构造是个大问题。联想到后缀排序和普通排序相比,有一个很大的不同是后缀之间是有关系的,比如第i个位置开始的后缀第i+1个位置开始的后缀相比就多了S(i)而已。为了利用这个关系,自然想到从后往前取后缀的方法。下面介绍一个实现非常简单的自创O(n^2)算法。
设sa(i)表示排名第i位的后缀的起始位置,rank(i)表示起始位置为i的后缀的排名,显然有sa和rank是互为反函数。假设现在处理到第i个后缀,那么构造二元组P(i) = <S(i), rank(i+1)>。对于之前已经处理过的开始位置为k(k > i)的后缀,有P(i) < P(k) 当且仅当S(i) < S(k) 或者 S(i) == S(k) && rank(i+1) < rank(k+1),由于i+1在i之前已经处理完毕,自然rank(i+1)此时也是有定义的。
使用类似插入排序的过程,使用上面定义的比较大小的方法从sa数组的后面开始找到suffix(i)的位置。在查找的过程中,对于所有P(k) > P(i),令rank(k) ++,因为suffix(i)跑到了它的前面去(读者可以想想看这样做为什么不会错?)。
基本原理很简单,当然代码也非常简洁,C描述如下:

[code:cpp]
void cons_suffix_array()
{
int i, j;

for ( i = n – 1; i > -1; –i ) {
for ( j = n – 2 – i; j > -1; –j ) {
if ( buf[i] > buf[ SA[j] ] || ( buf[i] == buf[ SA[j] ] && Rank[i+1] > Rank[ SA[j]+1 ] ) ) break;
SA[j+1] = SA[j];
++Rank[ SA[j] ];
}
SA[j+1] = i;
Rank[i] = j+1;
}
}

[/code]

这个算法最大的有点是编码容易,空间复杂度低,而且常数也低,同时也避免了之前算法的最坏情况。但是,O(n^2)的复杂度也的确很难将其应用到大数据量的场合,充其量可以用于快速建模或者Topcoder这种地方。不过,至少我们打开了思路,觉得SA并非intractable的东西,而且,文中对于已有结果的利用的思路也是可圈可点。
再进一步改进。现在的算法是借用的插入排序的模式,虽然后缀比较的复杂度从O(n)降到了O(1),但是排序本身的复杂度上升了,所以总的改进不算大。那么,能不能将以上后缀比较的方法和某个高效的排序算法(如快排、基排)结合产生更快的算法呢?
为了回答这个问题,先来探究上面比较方法的适用条件。显然,P(i)的两维都必须是可比较的,即对于rank的比较必须有意义。这一条件自然是对快排的沉重打击,因为快排在进行到最后一个元素之前其顺序都是混乱的。那基排呢,看得出在排序完成之前也不满足有序这一条件。所以,上面的算法看来是不能再改进了。
但是,基排有一个绝佳的性质就是部分有序,即在进行第k+1位排序时,若是MSB基排已经按前k位有序了,若是LSB基排则已经按后n-k-1位有序了。假设适用MSB,当后缀按照前k位时,构造二元组P(i) = < rank(i), rank(i+k) >,显然此时P(i)的比较是有意义的(因为rank(i)和rank(i+k)有定义),所以可以直接得到在2k下的排序结果。
这样,每次排序的长度增大了一倍而不是1位,使得基排只需要O(lgn)趟就可以完成,总的复杂度为O(nlgn)。这就是发表在第一篇后缀数组论文中著名的Manber算法,04年国家集训队中讲述的倍增算法就是说它。
懒了,不想写了,还是把这篇文章分两部分吧,下次再说其它更好的算法。

最大流的Shortest Augmenting Path Algorithm

我弱了,一直以为SAP算法是说的Edmonds-Karp算法(EK),现在才知道这个专指AhujaOrlin发明的一个非常高效的增广路系的最大流算法,根据实测,比之前常用的Dinic算法还要好一些,而且两者在代码上几乎无区别,所以非常容易实现。

网络流的基本知识我就不介绍了,直接对比Dinic说下SAP算法。Dinic利用了层次图的概念,使得每次找增广路不用从头开始,从而巧妙的使程序的趟数降到了n,复杂度下降为O(n^2m)SAP则是利用了preflow的距离标号思想,进一步的免去了建层次图的过程,同时preflowGap启发也可以利用上,所以即使复杂度上限仍然是O(n^2m),但实际跑起来要快许多。

这个算法怎么实现呢?回忆Dinic算法找阻塞流的过程,每次增广后,当前点退回到那条关键边的发出点w上,从w开始找新的可行弧。在这里,如果w没有可行弧,则退回到w的前一个点,而SAP在这一步上做了改进,先利用preflow中提高标号的方法(即将w的标号变为与wresidual network上相连的所有点中最小的一个加1),再继续寻找可行弧增广。

这样,有了提升标号过程,就不用重建层次图了。当源点的标号>=n时,程序就结束,或者当有Gap启发时,一旦发现某个高度出现断层,也可立即退出。

就这样,终于把POJ那道用我“优化”过无数次的Dinic却仍然过不到的Core Dual Cpu给过了。给我的教训:以后不能太自以为是,没了解过的算法就别太想当然。

最后,推荐大家尝试此算法。

二分图最大匹配的高效算法

本打算介绍Hopcroft-Karp算法(简称HK)和Unit Capacity Network Flow两种算法的,但无赖后者牵涉到的东西太多,就像和网络流有关的问题一样,主要讨论不同实现和启发对性能的影响。鉴于尔后有集中复习网络流的打算,这里便只说说第一种算法。

关于二分图匹配的匈牙利算法,现在一般作为图算法的入门学习内容,其原因就在于其编程简单,适合不求甚解式的强迫记忆,致使大多数人搞不懂如何从匈牙利树和增广轨入手设计算法(ms扯远了….)。匈牙利算法的理论复杂度是O(VE),但实践中我发现一般是跑得很快的,特别是再贪心一个初始匹配后。不过,我们似乎更关心如何降低这个算法的理论上界,于是有了HK算法。

回忆匈牙利算法,算法一共进行V轮,每一轮都从左部的一个未检查点开始查找增广轨。改进的思想是一次性找到多条增广轨,这样就可以减少算法进行的轮数。那么,这样是不是真的能改进效率,先给出算法描述如下:

 

设初始匹配M为空;

每一轮用当前图中的所有未匹配点求出相对于M的不相交最短增广轨集合P,如果P为空,则M此时已经是最大匹配,否则用P种的每条增广轨扩展M

那么,这个算法为什么高效?我的分析如下:

S1^S2表示求S1S2两个集合的对称差。

  1. M是一个匹配,p是关于M的一条最短增广轨,p’是关于M^p的一条最短增广轨,那么 p’ >= p + 2

    证明:几乎是废话,一条增广轨的长度是奇数,两个奇数当然至少相差2

  2. 一共只用进行O( sqrt(V) )轮。

    证明:设当前的最短增广轨长度L=sqrt(n),显然到达这个长度只需要进行L/2phase。设此时的匹配为M,而最大匹配为M*,则|M*| – |M| <= n / (sqrt(n) + 2) <= sqrt(n)。如果后续每阶段只共现一条增广轨,那么也只要sqrt(n)次就好,所以一共需要3/2sqrt(n)轮增广。

上面的证明用到了一个基本定理,即对于任意两个匹配M1M2,如果|M1| > |M2|,则M1^M2有至少|M1| – |M2|条不相交的增广轨。

现在,HK算法的复杂度也很显然了,就是O( sqrt(V) * E )。从实践来看,的确效率也较匈牙利有所提高,虽然又引入了更大的常数。另外,这个算法的实现也相当简单,每次只要把所有未配点一起bfs一次就好,大概只用写50行左右。不过相对dfs实现的匈牙利,显然编程复杂度提高了很多,而且从目前经验来看比赛中一般没有什么题能卡贪心初始的匈牙利的,所以推荐想开阔视野的朋友看看这个算法。

最后我还想请教下懂这个算法的人。很多书上书该算法和Dinic算法的思想一致,我一直没想明白,有谁能提示一下么?

ps:字符串算法的东西过段时间会继续,这两天被无聊的事情折腾得要死,太久没更新blog也很不爽,所以先放一篇上来。

KMP和Extend KMP算法

构思这篇文章是一段非常痛苦的过程,主要在于苦于无法用简洁的语言来表述我想说的话。KMP算法的介绍网络上可谓一搜一大把,目前引用最广泛的还谓Matrix67的那篇,不过就我当年初学时的记忆而言,还是读得晕乎乎的,但似乎还算有较强的可读性。

KMP本质上是构造了DFA并进行了模拟(DFANFA是什么我就不多说了,参考自动机理论的书籍),因此很显然一旦从模版T构造了自动机D,用D去匹配主串S的过程就是线性的。KMP最引人入胜的地方就在于构造D的自匹配过程,它充分利用了D是一个DAG的性质,使得构造过程也是线性的。

我在第一次学习KMP时就被告知,KMP的速度在平均情况下还不如brute force来得快,原因就在于随机字符串的自重复性并不高。尔后,又听说一些算法实际速度是KMP的好几倍,在平均情况下它们的复杂度是亚线性的,比如HorspoolSunday算法,于是就倍感神奇,原来终日崇拜的KMP并非最强大的。是的,线性并非是串匹配算法的复杂度下界,我们完全可以做得更好,而后面我将会说到很多的最优算法。

那么,能不能将KMP变得更快呢?回忆下,KMPDFA,如果能找到其对应的NFA,然后用惯用的位并行的手段来模拟NFA的执行,那么复杂度大约就是O(n/w),其中w位机器字长。这一思路是可行的,形式化如下:

DNFA,则D(i)1表示T(1),T(2)……T(i)S(1),S(2)……S(j)的后缀。状态转移表的构造也很简单,设该表为P,则如下伪码:

For i = 1 to m

P[ T(i) ] |= 1 << i

End

假设从j进行状态转移到j+1只要如下一个位运算:

D = ( ( D << 1 ) & P[ S(j+1) ] ) | 1

于是,这样就将KMP的大大提高了速度,不仅是理论上的,也是实践中的。

KMP的另外一个研究方向是Extend KMP(以下简称EK),它是说求得T与所有的S(i)的最长公共前缀(LCP),当然,要控制复杂度在线性以内。

EK我第一次听说是07baidu校园招聘的笔试题中,它当时的题目是求最长回文子串,当然这是一个耳熟能详,路人皆知可以用Suffix Array很好解决的问题。事后听一个同学说他写了三个算法:Suffix ArraySuffix TreeEK,当时就不明白EK是什么东西,但又没当面问他,于是这个东西就搁置了很久。知道后来北大的月赛一道题说可以用EK来做,我才终于从03年林希德的文章中开始认识到它,就像KMP一样,这个算法也一下就吸引了我。

Q(i)表示TS(i )的后缀的LCPP(i)表示TT(i)的后缀的LCP,那么和KMP一样,我们试图用P来求得Q,而P可以用自匹配求得,并且和求Q的过程相似。

我以求P为例简要说明一下。P(2)就直接匹配即可,从i = 3开始,如下:

k < iE(k) = k + P(k) – 1,且对所有j < i,有E(k) >= E(j)

那么,当E(k) >= i时,便可以推知T(i) = T( i – k + 1 ),于是如果P( i – k + 1 ) < E(k) – i + 1,那么P(i) = P( i – k + 1),否则P(i) >= P( i – k + 1 ),从E(k)开始向后匹配到E(i),有P(i) = E(i) – i + 1,并且更新 k = i

还有就是E(k) < i,肯定有E(k) = i – 1,不过这个不重要,重要的是直接从i开始做暴力匹配即可得到E(i),则P(i) = E(i) – i + 1,更新k = i

希望我把EK说清楚了,不过这种东西还是自己推导一下有意思,而且记忆周期更长。

最后来罗列下题目。KMP的经典题目是POJ 2185,是要找最小覆盖矩形,如果你认为懂KMP了就去尝试它。EK的经典题是POJ 3376,有一定挑战;当然还有就是上面说的最长回文子串,提醒下用分治+EK来做是其中一种方法。

嗯,下次打算说下Suffix Array,主要是它那个传说中的线性DC3算法,不过我现在还没把握能不能简单的把它说清楚,姑且认为可以吧。

字符串系列算法(一)

看到标题后面的序数,可知又是一个系列文章的开篇。不过,这次我是没有那么多时间去准备材料和学习新算法了,就尽力而为,能写什么就写什么吧。

字符串算法一直就是一个难点(个人观点),一般我在比赛中要遇到了字符串的难题,都会紧张好一阵,然后看到有人做了才会去细想。字符串类型的问题一般涉及两类知识点,一个是算法,一个是很高级的数据结构,总结一下主要会遇到以下这些:

算法方面:RKKMPExtended KMP,位并行算法,最小表示法;

数据结构方面:TrieAho-Corasick自动机,Trie图,Suffix Array

 

当然,上面列出的只是个人知道的并且在比赛中比较容易使用的东西,一些如后缀树或者BM算法,还远没有在Contestants中普及,究其原因主要是编程复杂度和思维难度上的问题,所以我也暂时回避它们(嗯,让大家失望了),这个系列也就只针对上述这些算法做讨论。

由于前几天才做了POJ2778,所以就先从这道题说起。

题目大意是给了很多个基因(总长度<100),求有多少种不同的长度为nDNA单链不包含其中任何一个基因,答案mod 100000n < 2000000000

拿到题第一反应自然是DP,然后就发现n很大,没有办法直接做。联想到DP是可以用矩阵来优化的(比如Fibnacci数列),这样n就看作是lgn了,于是现在的任务便是构造这么一个矩阵。

想到矩阵,再进一步,就想到图,貌似现在的问题就是要建图,然后求起点到终点有多少种不同的受限走法。下面引出上述的AC自动机和Trie图。

AC自动机可以看作是KMP算法在多模板时候next函数的扩展,因此它的主要用途就是多模匹配。要构造它,首先用这些模板构造一个Trie树,然后就是对Trie树上的每一个节点,求出它的next函数值,并且对于匹配成功的节点,标记为终结点。可以看出,AC自动机是DFA,这样一旦构造成功,后面的模式匹配过程就是线性时间的事情了。

如何求上述的next函数?很简单,从这棵Trie树的根节点开始,作一个BFS,对于每一个节点x(该节点的字符为c),用Father[x]表示它的父节点,Trie[x][k]表示x的第k个儿子,那么执行如下过程:

if ( Father[ x ] == ROOT ) next[ x ] = ROOT;

else {

y = next[ Father[x] ];

while ( y != ROOT && Trie[y][x] == NULL ) y = next[ y ];

next[ x ] = Trie[ y ][ c ];

}

 

这样,对于总长为L的模板,只要O( L )的时间就能够完成AC自动机的构造。

还有一点很重要,就是在next回溯过程中,对于终止态的传递要附带完成。举个例子,abcababc,后者是前者的前缀,因此前者的c字符也是一个终止状态。

好了,AC自动机就这么简单的构造出来了,可是这道题我们不用它,因为我们显然不是要做多模匹配,另外一个数据结构Trie图才是本文的主角。

Trie图其实和AC自动机很像,也是从Trie演化而来,只是对于图中的每个点,都有|E|个后代(E是字母表),其中有些是Trie原本有的,有些是添加的。添加的那些指向哪些节点的求法和AC自动机基本类似,所以我就不再赘述,大家想想。

附带提一句,如果还没弄明白这两个东西,可以参考zhuzeyuan 2004Maigo 2006的论文,其中前一篇论文把AC自动机叫做前缀树,大家注意一下就好。

好了,既然Trie图有了,后面的算法就很自然了。起点是根ROOT,终点是所有的非终结点,求从ROOT到这些点长为n的路径且任何一条路上不能包含终结点的条数,就是答案,我就不多说了,复杂度O( k^3 * lgn ),其中kTrie图的状态数。

开篇先说到这里,下次我想主要说一下KMP以及Extended KMP的算法本身以及应用。这是一个十分神奇的算法,要想学明白相当不容易,其应用更是让人惊叹,我先在此打个广告给大家留一点悬念,^_^

ps:应要求,写了一个简短的Trie图的说明,见
这里

在南京

南京赛区的现场比赛时间是2007年10月28号,至今已经快过去半年。

我本来不打算再去回顾这个历史,只希望它在我生命中作为一个符号而存在,永远不要去思考细节。可是,半年的挣扎,我始终无法做到肃清脑海中的阴影,对于那段时间周遭的嘲笑,我早就已经学会像木头一样去接受,可是,我却不能接受自己成为一块木头,我必须要完完全全的走出来。

可能会有很多人不理解,一次小小的比赛就能挫折到如此,可谓茶余饭后的笑料。不过,那种钻心的难受是我之前走过的20年路里面少有的,而且随着时间的增长伤口却没有愈合。我已经下了很久的决心一定要把这些事努力的回忆起来,记录下来,可是总能找到十足的借口加以推托,我已经被折磨了许久,我不想再耗下去。

第一天

这一天的主要任务是参加开幕式和下午的热身比赛。我们队本来是不打算参加开幕式,因为恰好和TopCoder四川的Tour比赛冲突。由于这次Tour的奖品非常丰厚,而且加之四川的有实力的人基本上掐指可算,因此参加比赛并拿奖的心情一直驱使着我们一到南航便一路狂奔到附近的网吧。非常不凑巧,那间网吧连Topcoder的速度非常慢,我们提前了十分钟到达,却依然没能在规定时间内注册上,于是,大家都很泄气的在那里观看了一下比赛的情况便走人。

回去只能去参加开幕式,就是些众所周知的领导讲话等形式,我们自然也没什么心情听。不过那时的我还比较乐观,还没有把奖品当一回事,于是一直和队友分析上报队伍的情况。我跟队友说了些什么我不想再说了,总之,一切都太乐观了。

开幕式很快结束,匆匆吃过午饭便赶到了比赛现场。现场的气氛的确非常的压抑,至少在比赛开始之前我是一直都没有很好的心情。不过还好,当打开试题时我紧张的心情一下就好了很多,精神也能迅速集中,这倒不是因为题目很简单,而恰恰相反,这次的热身题目有相当的难度,只有一个较简单的题,而且第一个过掉题目的队出现已经是20分钟左右的事情。我们三人拿到题目也迅速分析了题目的类型,A一开始没有思绪,B是一个推公式的概率,于是我的两个队友马上就扑了上去准备拿它开刀。我默默的看着C,一时由于数据规模太大而没有思路,于是和A题换着想了好几遍,终于还是找到了突破口就是dilworth定理和线段树,只是需要推导一个O(1)的判定条件。此时队友还在很卖力的推导B的公式,而排名上已经有20几个队过了B。于是我打算上C的代码,还好比较顺利,只用了15分钟左右写代码和调试并通过neotium给的一些数据。此时刚在几分钟前有人2Y了C,于是我抱着必死的心情忐忑的提交了C。由于第一次使用文件输入输出,忘记用fopen输入数据,所以第一次提交就TLE掉,发现问题马上修改,再提交,等待了一分钟居然返回Yes,我自己都惊了一下,这是我没有预料到的,于是很兴奋的告诉队友AC了,还搞的旁边的队伍好奇的探头来看Yes是什么样子,这样,我紧张的情绪一下完全就消失了。

队友继续在搞B,交了两次都是WA,于是他们继续改进公式。我现在已经非常轻松了,所以很快A也有了思路,发现就是一个最小权闭合子图,用最大流求解即可。这个时候,队友神奇的将B给过掉并且由于我已经知道了A的算法,同时发现机器能够连接到外网,于是起了邪念并且一发不可收拾,跑到POJ去下了一个最大流的代码,然后改了改便提交了,WA了一次后很轻松的AC掉它,至此做掉所有题目并且总排名在第8(最后结束时好像落后到10几),我自己认为比较满意(虽然我知道有很多队伍都是没有像我们一样认真在对待,但是我对自己今天的状态还是比较满意)这次热身结果,于是再次做了一个错误性的决定,提前走人。

现在回头来看,自己当时真的很傻。虽然我不相信RP这种东西,可是在比赛现场贴代码便已经反映出了我一个很大的问题,便是对自己的不自信和对模版的依赖程度太高。可是,我自以为是的表现让我没有去发现自己的问题,假象蒙蔽了双眼。

由于我和neotium没有带饭票,于是赶车回宾馆去拿。在回来等公交车的时候,发生了一件非常神奇的事情:我当时觉得无聊,然后正好手中有个硬币,不知为何我突然想到了具体数学上有一章说的一个游戏,意思就是每个人说一个人头和字的序列,然后开始执硬币看谁的说的序列最先出现。当时neotium还问我是不是知道什么结论才和他玩,我就说这个东西是有个好像是Conwey给出的结论,但是更简单的方法可以用Markov过程来做。然后我说Markov现在很火,估计明天可能会考,结果,正如大家所知道的,南京赛区的正式比赛中就恰好有这么一个完全一样的题,只不过用了0-9的数字代替头和字来形成序列,不过这道题我也是事后才知道,比赛的时候早就完全顾不上去做这么个题了。

回南航吃了饭,遇到川大的Windy7926778一行,和他们简单讨论了今天的题目,于是便坐车回了宾馆。当天晚上我真的非常轻松,因为我本来南京此行的目的不是要拿到多么辉煌的成绩,而今天的表现又证明我很在状态,于是,晚上和领队还有队友聊得非常开心,并没有觉得有什么事情正在悄悄的临近。

第二天

如期,我们到达了比赛现场。时隔一天我再次坐在这里,我并没有感到有如昨天那般的紧张,反而有些兴奋,这是个非常反常的状态,比紧张更可怕。

比赛开始了,我依然是看ABC三题。A我觉得题目描述有问题,和neotium商量了很久也没搞懂题意,于是就先放了。B我读完就觉得是个线段树的陈题,连半点思索的空间都没留便立马开始上机敲代码,而这时比赛才刚刚开始了十分钟而已,我的太过于乐观的心态没能让我把问题分析清楚就开始漫长的一个小时的coding。当我终于反应过来题目很难时(事实上最终也没有队过这题),已经是1个小时零几分的事情了,这个时候已经有比较多的队过了最后一题,还有零散的队伍过了其它题目。我的这一极其错误而没有被制止的行为(事后分析,因为队内当时只有我懂线段树这个东西,我只要比较肯定的说出是用这个方法,就不会有人帮我Review,我想这是一支队伍的悲哀所在),给最后的败局奠定了坚实的基础。

我马上放弃了B,转而去做一个最长弱重复子串的题目。这是一个非常经典的问题,并且是我事前在POJ做过的一个原题(但是我在比赛时却完全没有反应过来,以为是个新题目),一个编码很简单的算法就是后缀数组,另外一个按道理可以过的方法就是RK。因为已经落后了,我便没有想那么多,很快条件反射一样想到了二分枚举答案,自然接下来就联想到用RK来验证。离队友给我说题目不到一分钟的时间,我便自信满满的说方法绝对没问题,很好编码,我可以很快敲出来(那个时候,再一次没有人站出来反驳我)。队友自然很高兴,于是他们便又把机器让我,在一边思考已经有好几个队过了的最后一题。RK我自然比较熟悉,不过由于比较慌张,所以敲代码时居然犯下了一个非常低级的C语言编程错误,忘记将宏定义中的变量用括号括起来,这便是我队最后失利的主要原因之一,因为找到这个错误时比赛只剩下两个小时不到。这份代码,我前后不知道打印了多少份,开始是以为只用了一个基来取模判断精确度不高,后来一直改到三个基数取模仍然答案出不来。凭借我对RK算法的了解,我觉得要是真能把数据做到让如此的小概率事件发生,出题人去破个hash加密算法绝对能力有余。可是,这个时候比赛已经进行到了第三个小时,周围的队伍已经陆续升起了气球,前面的合肥工大的队伍更是已经有三个气球了。我的心情完全乱了(事实说明,我冷静分析把握大局的能力的确很差),摆在面前的午餐哪里有心思去吃,可是又实在没想出来那个字符串的题哪里有问题,而且队友负责的最后一个题也被告知很难没什么想法(估计他也完全被吓住了)。我一时不知道该怎么办,打算让队友来帮我找错,其实就是想让他帮我随便乱改一下说不定有效果。他不知道算法当然不敢碰,只是看到我一个地方为了省时间在for循环时写的i < N / 2好像有把握改改不会影响算法主体(我的N定义是 1 << 18,看出我为什么错了吧),于是随手编译一测试便通过了自己做的几个数据。我大喜,感觉是事情一下有了章法,看到了曙光,可是我连续交几次都是WA, 再次收到沉重的打击(至于另外一个问题,我便是比赛结束也没有发现,很可惜,估计是我后来又用的取一个基做模导致的)。

后来,我放下这道题目,和队友一起分析最后一题。由于3个半小时还完全没气球的状况已经让我感到非常的压抑,所以思路便根本打不开。后来不知谁说了个不知是不是动态规划的话,让我一下如梦初醒,立即判断这就是一个非常普通的DP(事实证明我这次没有判断错误),可是我觉得自己今天手背,就让neotium去写,我在一边给他看着。由于题目的确太简单,他很快就写好了,测试了样例没问题便直接交上去。

等待返回结果的这段时间是非常漫长的,我多么希望来的是个Yes,说不定士气就可以被重整。可是,命运就是这样,No仿佛如约而至,一下我的心理防线就开始有点崩溃了。我转过头去,看着窗外,一瞬间就像我已经忘记自己在比赛,看着远处的景色发神,此刻我真的已经绝望。好在我的队友没有和我一样,他立即开始着手改代码。代码真的没什么可以改的,我是看着他一个字一个字的敲进去的,如此简单的逻辑,不可能有错。事实也证明我的猜想,我们的代码的确没有错,在他也很无奈的提交了数次后终于也支持不住了。此刻我真的想哭,昨天的我绝对没想到今天在封board的时候我们居然还没有出题。大家郁闷了一阵,于是提出重新读题看有没有读错的地方。很快,neotium发现有个地方他之前理解错了,于是很快改掉便提交,系统终于在漫长的等待了近一分钟之后返回了Yes。

我们相视一笑,neotium说了句我们终于没有报零,我苦笑,我重来没有认为自己是以零的突破为目标的选手,我不至于这么弱。可是没办法,比赛始终是无情的,我的第一次ICPC区域赛生命只剩下50分钟可以走了。

我让neotium继续帮我看之前还没过的字符串题目,我便去开了个新题(事实证明这又是一个非常冒险的策略,是一次傻瓜式的冒险)。新题又是一个字符串,讲的是一个长字符串可以有多少种方法被许多给定的短字符串给组合出来。规模很大,长串有10,000的长度,短串总长为4000,每个串最多长100。最浅显的DP自然就是O( 4000 * 100 * 10,000 )的复杂度,不可能对。可是,在那种完全失利的情况下,我还能想出什么好方法?我能想到Trie?我自然又往RK上去想了,其结果大家肯定也很容易猜到。

比赛最后2分钟,我疯狂的对这两个字符串题目做修改然后提交,一共提了5次,可是上帝不与我同在,无一幸免。

于是,比赛结束,打开board一看,倒数第8名,恰好和昨天的热身赛结果做了颠倒。我想,这对于一个真正在过去的日子里付出了很大的努力来做ICPC的所有人都是难以接受的事实,我无奈的坐在那里,眼望着别的队伍手握气球在照相,我还能干什么,我已经没有事情可以做了,哨声一响我便成了孤儿,我唯一的事情便只有收拾书包。

那一刻的失落,可能只有经历过赛场的人才能有所体会。我们从成都带来的一切希望顷刻变为幻想,只有我和队友之前打趣说的拿铁牌我们就走路回家的玩笑话最有可能成为现实。其实,当时我的心情并没有现在说的那么难受,反而我是三人中最先开朗起来的,还不断的安慰同伴,就连我自己都不知道这原来是假象,在听Sun的讲座时我突然就崩溃了,一下子变得根本没了理智,毫无救药。

以至后来的颁奖典礼,我现在已经记不得是什么情景,只是还记得老师让我们三人上主席台合影留恋,我死都不去,我害怕面对台下那一双双冰冷的目视着你的眼睛,一时已经无法再有任何勇气可以将自己暴露在空气中。

坐了车,回了宾馆,老师自然心里对我们是完全的失望。不但如此,我仿佛开始觉得所有人的目光都开始变得犀利起来。

为了逃避现实,我们去了附近的网吧玩CS,夜深回去睡觉,没有人多说什么。

第三天

本来打算在南京旅游的行程也完全取消了,那个时候已经没有理智作出任何的判断。我们买的火车票是晚上6点,于是花了些时间去紫金山天文台看了看。整个过程我们三人已经开不出什么玩笑,而我更是连照片都害怕照,于是草草结束了行程,到市内随便找了个公园坐着,可是我无时不在回忆着比赛的场景,一幕一幕的重复,无休止。

终于,上了火车,离开了这个有着无限回忆的地方。

总结

关于我自己,我只想说是态度害了自己,好的心态是成功的一半绝对没错。我从来没有认为我们的队伍烂到如此地步,拿到这样的成绩的确和自己的努力太不成比例,这也是造成我无法接受结果的一个原因。

最后,我想说,ICPC的赛场是一个局面极端容易失控的地方,每个选手都应该在赛前做一系列的充分准备,不要只局限于嘴上讨论遇到什么问题该怎么解决,在比赛时乱起来是没有章法的。还有就是,错一次,即使和我一样在南京遇到与此的挫折,也不一定就可以完全学会如何控制自己的情绪,比如我,在之后的成都赛区又一次遇到了类似的问题,我依然处理的很无奈。第三点非常重要的事情,如果你不是一个非常厉害可以一手遮天的人,那么请不要在平时训练时和队友分头负责做哪些类型的题目,一定得有个交集,不然将导致最后互相听不懂对方说的什么东西,这有违团队理念,必将受到制裁。