动态规划是本书介绍的五种算法设计方法中难度最大的一种，它建立在最优原则的基础上。采用动态规划方法，可以优雅而高效地解决许多用贪婪算法或分而治之算法无法解决的问题。在介绍动态规划的原理之后，本章将分别考察动态规划方法在解决背包问题、图象压缩、矩阵乘法链、最短路径、无交叉子集和元件折叠等方面的应用。
3.1算法思想
和贪婪算法一样，在动态规划中，可将一个问题的解决方案视为一系列决策的结果。不同的是，在贪婪算法中，每采用一次贪婪准则便做出一个不可撤回的决策，而在动态规划中，还要考察每个最优决策序列中是否包含一个最优子序列。
例3-1[最短路经]考察图12-2中的有向图。假设要寻找一条从源节点s=1到目的节点d=5的最短路径，即选择此路径所经过的各个节点。第一步可选择节点2，3或4。假设选择了节点3，则此时所要求解的问题变成：选择一条从3到5的最短路径。如果3到5的路径不是最短的，则从1开始经过3和5的路径也不会是最短的。例如，若选择的子路径（非最短路径）是3，2，5(耗费为9)，则1到5的路径为1，3，2，5(耗费为11)，这比选择最短子路径3，4，5而得到的1到5的路径1，3，4，5(耗费为9)耗费更大。
所以在最短路径问题中，假如在的第一次决策时到达了某个节点v，那么不管v是怎样确定的，此后选择从v到d的路径时，都必须采用最优策略。
例3-2[0/1背包问题]考察13.4节的0/1背包问题。如前所述，在该问题中需要决定x1..xn的值。假设按i=1，2，.，n的次序来确定xi的值。如果置x1=0，则问题转变为相对于其余物品（即物品2，3，.，n），背包容量仍为c的背包问题。若置x1=1，问题就变为关于最大背包容量为c-w1的问题。现设r{c，c-w1}为剩余的背包容量。
在第一次决策之后，剩下的问题便是考虑背包容量为r时的决策。不管x1是0或是1，[x2，.，xn]必须是第一次决策之后的一个最优方案，如果不是，则会有一个更好的方案[y2，.，yn]，因而[x1，y2，.，yn]是一个更好的方案。
假设n=3,w=[100,14,10],p=[20,18,15],c=116。若设x1=1，则在本次决策之后，可用的背包容量为r=116-100=16。[x2，x3]=[0,1]符合容量限制的条件，所得值为15，但因为[x2，x3]=[1，0]同样符合容量条件且所得值为18，因此[x2，x3]=[0，1]并非最优策略。即x=[1，0，1]可改进为x=[1，1，0]。若设x1=0，则对于剩下的两种物品而言，容量限制条件为116。总之，如果子问题的结果[x2，x3]不是剩余情况下的一个最优解，则[x1，x2，x3]也不会是总体的最优解。
例3-3[航费]某航线价格表为：从亚特兰大到纽约或芝加哥，或从洛杉矶到亚特兰大的费用为$100；从芝加哥到纽约票价$20；而对于路经亚特兰大的旅客，从亚特兰大到芝加哥的费用仅为$20。从洛杉矶到纽约的航线涉及到对中转机场的选择。如果问题状态的形式为（起点，终点），那么在选择从洛杉矶到亚特兰大后，问题的状态变为（亚特兰大，纽约）。从亚特兰大到纽约的最便宜航线是从亚特兰大直飞纽约，票价$100。而使用直飞方式时，从洛杉矶到纽约的花费为$200。不过，从洛杉矶到纽约的最便宜航线为洛杉矶-亚特兰大-芝加哥-纽约，其总花费为$140（在处理局部最优路径亚特兰大到纽约过程中选择了最低花费的路径：亚特兰大-芝加哥-纽约）。
如果用三维数组（tag，起点，终点）表示问题状态，其中tag为0表示转飞，tag为1表示其他情形，那么在到达亚特兰大后，状态的三维数组将变为（0，亚特兰大，纽约），它对应的最优路径是经由芝加哥的那条路径。
当最优决策序列中包含最优决策子序列时，可建立动态规划递归方程（dynamic-programmingrecurrenceequation），它可以帮助我们高效地解决问题。
例3-4[0/1背包]在例3-2的0/1背包问题中，最优决策序列由最优决策子序列组成。假设f(i,y)表示例15-2中剩余容量为y，剩余物品为i，i+1，.，n时的最优解的值，即：和利用最优序列由最优子序列构成的结论，可得到f的递归式。f(1,c)是初始时背包问题的最优解。可使用（15-2）式通过递归或迭代来求解f(1,c)。从f(n,*)开始迭式，f(n,*)由（15-1）式得出，然后由（15-2）式递归计算f(i,*)(i=n-1，n-2，.，2)，最后由（15-2）式得出f(1,c)。
对于例15-2，若0≤y＜10，则f(3,y)=0；若y≥10，f(3,y)=15。利用递归式（15-2），可得f(2,y)=0(0≤y＜10)；f（2，y）=15（10≤y＜14）；f（2，y）=18（14≤y＜24）和f（2，y）=33（y≥24）。因此最优解f(1,116)=max{f（2，116），f（2，116-w1）+p1}=max{f（2，116），f（2，16）+20}=max{33，38}=38。
现在计算xi值，步骤如下：若f(1,c)=f(2,c)，则x1=0，否则x1=1。接下来需从剩余容量c-w1中寻求最优解，用f(2,c-w1)表示最优解。依此类推，可得到所有的xi(i=1.n)值。
在该例中，可得出f(2,116)=33≠f(1,116)，所以x1=1。接着利用返回值38-p1=18计算x2及x3，此时r=116-w1=16，又由f(2,16)=18，得f(3,16)=14≠f(2,16)，因此x2=1，此时r=16-w2=2，所以f(3,2)=0，即得x3=0。
动态规划方法采用最优原则（principleofoptimality）来建立用于计算最优解的递归式。所谓最优原则即不管前面的策略如何，此后的决策必须是基于当前状态（由上一次决策产生）的最优决策。由于对于有些问题的某些递归式来说并不一定能保证最优原则，因此在求解问题时有必要对它进行验证。若不能保持最优原则，则不可应用动态规划方法。在得到最优解的递归式之后，需要执行回溯（traceback）以构造最优解。
编写一个简单的递归程序来求解动态规划递归方程是一件很诱人的事。然而，正如我们将在下文看到的，如果不努力地去避免重复计算，递归程序的复杂性将非常可观。如果在递归程序设计中解决了重复计算问题时，复杂性将急剧下降。动态规划递归方程也可用迭代方式来求解，这时很自然地避免了重复计算。尽管迭代程序与避免重复计算的递归程序有相同的复杂性，但迭代程序不需要附加的递归栈空间，因此将比避免重复计算的递归程序更快。
3.2应用
3.2.10/1背包问题
1.递归策略
在例3-4中已建立了背包问题的动态规划递归方程，求解递归式（15-2）的一个很自然的方法便是使用程序15-1中的递归算法。该模块假设p、w和n为输入，且p为整型，F(1,c)返回f(1,c)值。
程序15-1背包问题的递归函数
intF(inti,inty)
{//返回f(i,y).
if(i==n)return(y if(y returnmax(F(i+1,y),F(i+1,y-w[i])+p[i]); } 例3-5设n=5，p=[6,3,5,4,6]，w=[2,2,6,5,4]且c=10,求f(1,10)。为了确定f(1,10)，调用函数F(1,10)。递归调用的关系如图15-1的树型结构所示。每个节点用y值来标记。对于第j层的节点有i=j，因此根节点表示F(1,10)，而它有左孩子和右孩子，分别对应F(2,10)和F(2,8)。总共执行了28次递归调用。但我们注意到，其中可能含有重复前面工作的节点，如f(3,8)计算过两次，相同情况的还有f(4,8)、f(4,6)、f(4,2)、f(5,8)、f(5,6)、f(5,3)、f(5,2)和f(5,1)。如果保留以前的计算结果，则可将节点数减至19，因为可以丢弃图中的阴影节点。 正如在例3-5中所看到的，程序15-1做了一些不必要的工作。为了避免f(i,y)的重复计算，必须定义一个用于保留已被计算出的f(i,y)值的表格L，该表格的元素是三元组(i,y,f(i,y))。在计算每一个f(i,y)之前，应检查表L中是否已包含一个三元组(i,y,*)，其中*表示任意值。如果已包含，则从该表中取出f(i,y)的值，否则，对f(i,y)进行计算并将计算所得的三元组(i,y,f(i,y))加入表L。L既可以用散列（见7.4节）的形式存储，也可用二叉搜索树(见11章)的形式存储。 2.权为整数的迭代方法 当权为整数时，可设计一个相当简单的算法（见程序15-2）来求解f(1,c)。该算法基于例3-4所给出的策略，因此每个f(i,y)只计算一次。程序15-2用二维数组f[][]来保存各f的值。而回溯函数Traceback用于确定由程序15-2所产生的xi值。函数Knapsack的复杂性为(nc），而Traceback的复杂性为(n)。 程序15-2f和x的迭代计算 template voidKnapsack(Tp[],intw[],intc,intn,T**f) {//对于所有i和y计算f[i][y] //初始化f[n][] for(inty=0;y<=yMax;y++) f[n][y]=0; for(inty=w[n];y<=c;y++) f[n][y]=p[n]; //计算剩下的f for(inti=n-1;i>1;i--){ f[i][y]=f[i+1][y]; for(inty=w[i];y<=c;y++) f[i][y]=max(f[i+1][y],f[i+1][y-w[i]]+p[i]); f[1][c]=f[2][c]; if(c>=w[1]) f[1][c]=max(f[1][c],f[2][c-w[1]]+p[1]); voidTraceback(T**f,intw[],intc,intn,intx[]) {//计算x for(inti=1;i if(f[i][c]==f[i+1][c])x[i]=0; else{x[i]=1; c-=w[i];} x[n]=(f[n][c])1:0; 3.元组方法（选读） 程序15-2有两个缺点：1)要求权为整数；2)当背包容量c很大时，程序15-2的速度慢于程序15-1。一般情况下，若c＞2n，程序15-2的复杂性为W(n2n)。可利用元组的方法来克服上述两个缺点。在元组方法中，对于每个i，f(i,y)都以数对(y,f(i,y))的形式按y的递增次序存储于表P(i)中。同时，由于f(i,y)是y的非递减函数，因此P(i)中各数对(y,f(i,y))也是按f(i,y)的递增次序排列的。 例3-6条件同例3-5。对f的计算如图15-2所示。当i=5时，f由数对集合P(5)=[(0,0),(4,6)]表示。而P(4)、P(3)和P(2)分别为[(0,0),(4,6),(9,10)]、[(0,0)(4,6),(9,10),(10,11)]和[(0,0)(2,3)(4,6)(6,9)(9,10)(10,11)]。 为求f(1,10)，利用式（15-2）得f(1,10)=max｛f(2,10)，f(2,8)+p1｝。由P(2)得f(2,10)=11、f(2,8)=9(f(2,8)=9来自数对(6，9))，因此f(1,10)=max｛11,15｝=15。现在来求xi的值，因为f(1,10)=f(2,6)+p1，所以x1=1；由f(2,6)=f(3,6-w2)+p2=f(3,4)+p2，得x2=1；由f(3,4)=f(4,4)=f(5,4)得x3=x4=0；最后，因f(5,4)≠0得x5=1。 检查每个P(i)中的数对，可以发现每对(y,f(i,y))对应于变量xi,.,xn的0/1赋值的不同组合。设（a,b）和（c,d）是对应于两组不同xi,.,xn的0/1赋值，若a≥c且b＜d，则(a,b)受(b,c)支配。被支配者不必加入P(i)中。若在相同的数对中有两个或更多的赋值，则只有一个放入P(i)。假设wn≤C，P(n)=[(0,0),(wn,pn)]，P（n）中对应于xn的两个数对分别等于0和1。对于每个i，P(i)可由P(i+1)得出。首先，要计算数对的有序集合Q，使得当且仅当wi≤s≤c且(s-wi,t-pi)为P(i+1)中的一个数对时，（s,t）为Q中的一个数对。现在Q中包含xi=1时的数对集，而P(i+1)对应于xi=0的数对集。接下来，合并Q和P(i+1)并删除受支配者和重复值即可得到P(i)。 如6.4.3节定义的，数字化图像是m×m的像素阵列。假定每个像素有一个0~255的灰度值。因此存储一个像素至多需8位。若每个像素存储都用最大位8位，则总的存储空间为8m2位。为了减少存储空间，我们将采用变长模式（variablebitscheme），即不同像素用不同位数来存储。像素值为0和1时只需1位存储空间；值2、3各需2位；值4，5，6和7各需3位；以此类推，使用变长模式的步骤如下： 1)图像线性化根据图15-3a中的折线将m×m维图像转换为1×m2维矩阵。 2)分段将像素组分成若干个段，分段原则是：每段中的像素位数相同。每个段是相邻像素的集合且每段最多含256个像素，因此，若相同位数的像素超过256个的话，则用两个以上的段表示。 3)创建文件创建三个文件：SegmentLength,BitsPerPixel和Pixels。第一个文件包含在2)中所建的段的长度(减1)，文件中各项均为8位长。文件BitsPerPixel给出了各段中每个像素的存储位数（减1），文件中各项均为3位。文件Pixels则是以变长格式存储的像素的二进制串。 4)压缩文件压缩在3)中所建立的文件，以减少空间需求。 上述压缩方法的效率（用所得压缩率表示）很大程度上取决于长段的出现频率。 例3-8考察图15-3b的4×4图像。按照蛇形的行主次序，灰度值依次为10，9，12，40，50，35，15，12，8，10，9，15，11，130，160和240。各像素所需的位数分别为4，4，4，6，6，6，4，4，4，4，4，4，4，8，8和8，按等长的条件将像素分段，可以得到4个段[10，9，12]、[40，50，35]、[15,12,8,10,9,15,11]和［130,160,240］。因此，文件SegmentLength为2，2，6，2；文件BitsPerSegment的内容为3，5，3，7；文件Pixels包含了按蛇形行主次序排列的16个灰度值，其中头三个各用4位存储，接下来三个各用6位，再接下来的七个各用4位，最后三个各用8位存储。因此存储单元中前30位存储了前六个像素： 101010011100111000110010100011 这三个文件需要的存储空间分别为：文件SegmentLength需32位；BitsPerSegment需12位；Pixels需82位，共需126位。而如果每个像素都用8位存储，则存储空间需8×16=128位，因而在本例图像中，节省了2位的空间。 假设在2)之后，产生了n个段。段标题（segmentheader）用于存储段的长度以及该段中每个像素所占用的位数。每个段标题需11位。现假设li和bi分别表示第i段的段长和该段每个像素的长度，则存储第i段像素所需要的空间为li*bi。在2)中所得的三个文件的总存储空间为11n+ni=1libi。可通过将某些相邻段合并的方式来减少空间消耗。如当段i和i+1被合并时，合并后的段长应为li+li+1。此时每个像素的存储位数为max{bi，bi+1}位。尽管这种技术增加了文件Pixels的空间消耗，但同时也减少了一个段标题的空间。 例3-9如果将例15-8中的第1段和第2段合并，合并后，文件SegmentLength变为5，6，2，BitsPerSegment变为5，3，7。而文件Pixels的前36位存储的是合并后的第一段：001010001001001100111000110010100011其余的像素（例15-8第3段）没有改变。因为减少了1个段标题，文件SegmentLength和BitsPerPixel的空间消耗共减少了11位，而文件Pixels的空间增加6位，因此总共节约的空间为5位，空间总消耗为121位。 我们希望能设计一种算法，使得在产生n个段之后，能对相邻段进行合并，以便产生一个具有最小空间需求的新的段集合。在合并相邻段之后，可利用诸如LZW法（见7.5节）和霍夫曼编码（见9.5.3节）等其他技术来进一步压缩这三个文件。 令sq为前q个段的最优合并所需要的空间。定义s0=0。考虑第i段(i＞0)，假如在最优合并C中，第i段与第i-1，i-2，.，i-r+1段相合并，而不包括第i-r段。合并C所需要的空间消耗等于：第1段到第i-r段所需空间+lsum(i-r+1,i)*bmax(i-r+1,i)+11 其中lsum(a,b)=bj=a lj，bmax(a,b)=max{ba,...,bb}。假如在C中第1段到第i-r段的合并不是最优合并，那么需要对合并进行修改，以使其具有更小的空间需求。因此还必须对段1到段i-r进行最优合并，也即保证最优原则得以维持。故C的空间消耗为： si=si-r+lsum（i-r+1,i）*bmax（i-r+1,i）+11 r的值介于1到i之间，其中要求lsum不超过256(因为段长限制在256之内)。尽管我们不知道如何选择r，但我们知道，由于C具有最小的空间需求，因此在所有选择中，r必须产生最小的空间需求。 假定kayi表示取得最小值时k的值，sn为n段的最优合并所需要的空间，因而一个最优合并可用kay的值构造出来。 例3-10假定在2)中得到五个段，它们的长度为[6，3，10，2，3]，像素位数为[1，2，3，2，1]，要用公式（15-3）计算sn，必须先求出sn-1，.，s0的值。s0为0，现计算s1：s1=s0+l1*b1+11=17kay1=1s2由下式得出： s2=min{s1+l2b2,s0+(l1+l2)*max{b1,b2}}+11=min{17+6,0+9*2}+11=29 kay2=2 以此类推，可得s1.s5=[17，29，67，73，82]，kay1.kay5=[1，2，2，3，4]。因为s5=82，所以最优空间合并需82位的空间。可由kay5导出本合并的方式，过程如下：因为kay5=4，所以s5是由公式（15-3）在k=4时取得的，因而最优合并包括：段1到段(5-4)=1的最优合并以及段2，3，4和5的合并。最后仅剩下两个段：段1以及段2到段5的合并段。 1.递归方法 用递归式（15-3）可以递归地算出si和kayi。程序15-3为递归式的计算代码。l，b，和kay是一维的全局整型数组，L是段长限制（256），header为段标题所需的空间(11)。调用S(n)返回sn的值且同时得出kay值。调用Traceback(kay,n)可得到最优合并。 现讨论程序15-3的复杂性。t(0)=c（c为一个常数）：（n＞0），因此利用递归的方法可得t(n)=O(2n)。Traceback的复杂性为(n)。 程序15-3递归计算s,kay及最优合并 intS(inti) {//返回S(i)并计算kay[i] if(i==0)return0; //k=1时,根据公式（15-3）计算最小值 intlsum=l[i],bmax=b[i]; ints=S(i-1)+lsum*bmax; kay[i]=1; //对其余的k计算最小值并求取最小值 for(intk=2;k<=i&&lsum+l[i-k+1]<=L;k++){ lsum+=l[i-k+1]; if(bmax intt=S(i-k); if(s>t+lsum*bmax){ s=t+lsum*bmax; kay[i]=k;} returns+header; voidTraceback(intkay[],intn) {//合并段 if(n==0)return; Traceback(kay,n-kay[n]); cout<<"Newsegmentbeginsat"<<(n-kay[n]+1)< 2.无重复计算的递归方法 通过避免重复计算si，可将函数S的复杂性减少到(n)。注意这里只有n个不同的si。 例3-11再考察例15-10中五个段的例子。当计算s5时，先通过递归调用来计算s4，.，s0。计算s4时，通过递归调用计算s3，.，s0，因此s4只计算了一次，而s3计算了两次，每一次计算s3要计算一次s2，因此s2共计算了四次，而s1重复计算了16次！可利用一个数组s来保存先前计算过的si以避免重复计算。改进后的代码见程序15-4，其中s为初值为0的全局整型数组。 程序15-4避免重复计算的递归算法 {//计算S(i)和kay[i] //避免重复计算 if(s[i]>0)returns[i];//已计算完 //计算s[i] //首先根据公式（15-3）计算k=1时最小值 s[i]=S(i-1)+lsum*bmax; //对其余的k计算最小值并更新 if(s[i]>t+lsum*bmax){ s[i]=t+lsum*bmax; s[i]+=header; returns[i]; 3.迭代方法 倘若用式（15-3）依序计算s1,.,sn，便可得到一个复杂性为(n)的迭代方法。在该方法中，在si计算之前，sj必须已计算好。该方法的代码见程序15-5，其中仍利用函数Traceback（见程序15-3）来获得最优合并。 程序15-5迭代计算s和kay voidVbits(intl[],intb[],intn,ints[],intkay[]) {//计算s[i]和kay[i] intL=256,header=11; s[0]=0; //根据式（15-3）计算s[i] for(inti=1;i<=n;i++){ //k=1时,计算最小值 intlsum=l, bmax=b[i]; s[i]=s[i-1]+lsum*bmax; if(s[i]>s[i-k]+lsum*bmax){ s[i]=s[i-k]+lsum*bmax; 3.2.3矩阵乘法链 下面举一个得益于选择合适秩序计算A*B*C矩阵的实例：考虑两个3维图像的匹配。图像匹配问题的要求是，确定一个图像需旋转、平移和缩放多少次才能逼近另一个图像。实现匹配的方法之一便是执行约100次迭代计算，每次迭代需计算12×1个向量T： T=A(x,y,z)*B(x,y,z)*C(x,y,z) 其中A，B和C分别为12×3，3×3和3×1矩阵。(x,y,z)为矩阵中向量的坐标。设t表示计算A(x,y,z)*B(x,y,z)*C(x,y,z)的计算量。假定此图像含256×256×256个向量，在此条件中，这100个迭代所需的总计算量近似为100*2563*t≈1.7*109t。若三个矩阵是按由左向右的顺序相乘的，则t=12*3*3+12*3*1=144；但如果从右向左相乘，t=3*3*1+12*3*1=45。由左至右计算约需2.4*1011个操作，而由右至左计算大概只需7.5*1010个操作。假如使用一个每秒可执行1亿次操作的计算机，由左至右需40分钟，而由右至左只需12.5分钟。 在计算矩阵运算A*B*C时，仅有两种乘法顺序（由左至右或由右至左），所以可以很容易算出每种顺序所需要的操作数，并选择操作数比较少的那种乘法顺序。但对于更多矩阵相乘来说，情况要复杂得多。如计算矩阵乘积M1×M2×.×Mq，其中Mi是一个ri×ri+1矩阵(1≤i≤q)。不妨考虑q=4的情况，此时矩阵运算A*B*C*D可按以下方式（顺序）计算： A*((B*C)*D)A*(B*(C*D))(A*B)*(C*D)(A*(B*C))*D 不难看出计算的方法数会随q以指数级增加。因此，对于很大的q来说，考虑每一种计算顺序并选择最优者已是不切实际的。 与求解0/1背包及图像压缩问题一样，本递归方法也须避免重复计算c(i,j)和kay(i,j)，否则算法的复杂性将会非常高。 例3-13设q=5和r=（10,5,1,10,2,10），式中待求的c中有四个c的s=0或1，因此用动态规划方法可立即求得它们的值：c(1,1)=c(5,5)=0;c(1,2)=50;c(4,5)=200。现计算C(2，5)：c(2,5)=min{c(2,2)+c(3,5)+50,c(2,3)+c(4,5)+500,c(2,4)+c(5,5)+100}（15-5）其中c(2,2)=c(5,5)=0；c(2,3)=50；c(4,5)=200。再用递归式计算c(3,5)及c(2,4):c(3,5)=min{c(3,3)+c(4,5)+100,c(3,4)+c(5,5)+20}=min{0+200+100,20+0+20}=40c(2,4)=min{c(2,2)+c(3,4)+10,c(2,3)+c(4,4)+100}=min{0+20+10,50+10+20}=30由以上计算还可得kay(3,5)=4，kay(2,4)=2。现在,计算c(2,5)所需的所有中间值都已求得，将它们代入式（15-5）得： c(2,5)=min{0+40+50,50+200+500,30+0+100}=90且kay(2,5)=2 再用式（15-4）计算c(1,5)，在此之前必须算出c(3,5)、c(1,3)和c(1,4)。同上述过程，亦可计算出它们的值分别为40、150和90，相应的kay值分别为4、2和2。代入式（15-4）得： c(1,5)=min{0+90+500,50+40+100,150+200+1000,90+0+200}=190且kay(1,5)=2 此最优乘法算法的消耗为190，由kay(1,5)值可推出该算法的最后一步，kay(1,5)等于2，因此最后一步为M12×M35，而M12和M35都是用最优法计算而来。由kay(1,2)=1知M12等于M11×M22，同理由kay(3,5)=4得知M35由M34×M55算出。依此类推，M34由M33×M44得出。因而此最优乘法算法的步骤为： M11×M22=M12 M33×M44=M34 M34×M55=M35 M12×M35=M15 计算c(i,j)和kay(i,j)的递归代码见程序15-6。在函数C中，r为全局一维数组变量，kay是全局二维数组变量，函数C返回c(ij)之值且置kay[a][b]=kay(a,b)(对于任何a,b)，其中c(a,b)在计算c(i,j)时皆已算出。函数Traceback利用函数C中已算出的kay值来推导出最优乘法算法的步骤。 设t(q)为函数C的复杂性，其中q=j-i+1（即Mij是q个矩阵运算的结果）。当q为1或2时，t(q)=d，其中d为一常数；而q>2时，t(q)=2q-1k=1t(k)+eq，其中e是一个常量。因此当q＞2时，t(q)＞2t(q-1)+e，所以t(q)=W(2q)。函数Traceback的复杂性为(q)。 程序15-6递归计算c(i,j)和kay(i,j) intC(inti,intj) {//返回c(i,j)且计算k(i,j)=kay[i][j] if(i==j)return0;//一个矩阵的情形 if(i==j-1){//两个矩阵的情形 kay[i][i+1]=i; returnr[i]*r[i+1]*r[r+2];} //多于两个矩阵的情形 //设u为k=i时的最小值 intu=C(i,i)+C(i+1,j)+r[i]*r[i+1]*r[j+1]; kay[i][j]=i; //计算其余的最小值并更新u for(intk=i+1;k intt=C(i,k)+C(k+1,j)+r[i]*r[k+1]*r[j+1]; if(r u=t; kay[i][j]=k; returnu; voidTraceback(inti,intj,int**kay) {//输出计算Mij的最优方法 if(i==j)return; Traceback(i,kay[i][j],kay); Traceback(kay[i][j]+1,j,kay); cout<<"MultiplyM"< cout<<"andM"<<(kay[i][j]+1)<<","< 若避免再次计算前面已经计算过的c（及相应的kay），可将复杂性降低到（q3）。而为了避免重复计算，需用一个全局数组c[][]存储c(i,j)值，该数组初始值为0。函数C的新代码见程序15-7： 程序15-7无重复计算的c(i,j)计算方法 {//返回c(i,j)并计算kay(i,j)=kay[I][j] //检查是否已计算过 if（c[i][j]>）returnc[i][j]; //若未计算,则进行计算 c[i][j]=r[i]*r[i+1]*r[i+2]; returnc[i][j];} for(intk==i+1;k if(t kay[i][j]=k;} c[i][j]=u; c的动态规划递归式可用迭代的方法来求解。若按s=2，3，.，q-1的顺序计算c(i,i+s)，每个c和kay仅需计算一次。 例3-14考察例3-13中五个矩阵的情况。先初始化c(i,i)(0≤i≤5)为0，然后对于i=1,.,4分别计算c(i,i+1)。c(1,2)=r1r2r3=50，c(2,3)=50，c(3,4)=20和c(4,5)=200。相应的kay值分别为1，2，3和4。 当s=2时，可得： c(1,3)=min{c(1,1)+c(2,3)+r1r2r4,c(1,2)+c(3,3)+r1r3r4}=min =150 且kay(1,3)=2。用相同方法可求得c(2,4)和c(3,5)分别为30和40，相应kay值分别为2和3。 当s=3时，需计算c(1,4)和c(2,5)。计算c(2,5)所需要的所有中间值均已知(见(15-5)式)，代入计算公式后可得c(2,5)=90，kay(2,5)=2。c(1,4)可用同样的公式计算。最后，当s=4时，可直接用（15-4）式来计算c(1,5)，因为该式右边所有项都已知。 计算c和kay的迭代程序见函数MatrixChain（见程序15-8），该函数的复杂性为(q3)。计算出kay后同样可用程序15-6中的Traceback函数推算出相应的最优乘法计算过程。 程序15-8c和kay的迭代计算 voidMatrixChain(intr[],intq,int**c,int**kay) {//为所有的Mij计算耗费和kay //初始化c[i][i],c[i][i+1]和kay[i][i+1] for(inti=1;i c[i][i]=0; c[i][i+1]=r[i]*r[i+1]*r[i+2]; c[q][q]=0; //计算余下的c和kay for(ints=2;s for(inti=1;i<=q-s;i++){ //k=i时的最小项 c[i][i+s]=c[i][i]+c[i+1][i+s]+r[i]*r[i+1]*r[i+s+1]; kay[i][i+s]=i; //余下的最小项 for(intk=i+1;k intt=c[i][k]+c[k+1][i+s]+r[i]*r[k+1]*r[i+s+1]; if(t c[i][i+s]=t; kay[i][i+s]=k;} 3.2.4最短路径 假设G为有向图，其中每条边都有一个长度（或耗费），图中每条有向路径的长度等于该路径中各边的长度之和。对于每对顶点(i,j)，在顶点i与j之间可能有多条路径，各路径的长度可能各不相同。我们定义从i到j的所有路径中，具有最小长度的路径为从i到j的最短路径。 例3-15如图15-4所示。从顶点1到顶点3的路径有 1)1,2,5,3 2)1,4,3 3)1,2,5,8,6,3 4)1,4,6,3 由该图可知,各路径相应的长度为10、28、9、27，因而路径3)是该图中顶点1到顶点3的最短路径。 在所有点对最短路径问题（all-pairsshorest-pathsproblem）中，要寻找有向图G中每对顶点之间的最短路径。也就是说，对于每对顶点(i,j)，需要寻找从i到j的最短路径及从j到i的最短路径。因此对于一个n个顶点的图来说，需寻找p=n(n-1)条最短路径。假定图G中不含有长度为负数的环路，只有在这种假设下才可保证G中每对顶点(i,j)之间总有一条不含环路的最短路径。当有向图中存在长度小于0的环路时，可能得到长度为－∞的更短路径，因为包含该环路的最短路径往往可无限多次地加上此负长度的环路。 设图G中n个顶点的编号为1到n。令c(i,j,k)表示从i到j的最短路径的长度，其中k表示该路径中的最大顶点。因此，如果G中包含边 例3-16考察图15-4。若k=0,1,2,3，则c(1,3,k)=∞；c(1,3,4)=28；若k=5,6,7，则c(1,3,k)=10；若k=8,9,10，则c(1,3,k)=9。因此1到3的最短路径长度为9。对于任意k＞0，如何确定c(i,j,k)呢？中间顶点不超过k的i到j的最短路径有两种可能：该路径含或不含中间顶点k。若不含，则该路径长度应为c(i,j,k-1)，否则长度为c(i,k,k-1)+c(k,j,k-1)。c(i,j,k)可取两者中的最小值。因此可得到如下递归式： c(i,j,k)=min{c(i,j,k-1),c(i,k,k-1)+c(k,j,k-1)}，k＞0 //寻找最短路径的长度 //初始化c（i，j，1） for（inti=1；i<=n;i++） for(intj=1;j<=n;j++) c(i,j,0)=a(i,j);//a是长度邻接矩阵 //计算c(i,j,k)(0 for(intk=1;k<=n;k++) for(inti=1;i<=n;i++) if(c(i,k,k-1)+c(k,j,k-1) c(i,j,k)=c(i,k,k-1)+c(k,j,k-1); elsec(i,j,k)=c(i,j,k-1); 图15-5最短路径算法的伪代码 注意到对于任意i，c(i,k,k)=c(i,k,k-1)且c(k,i,k)=c(k,i,k-1)，因而，若用c(i,j)代替图15-5的c(i,j,k)，最后所得的c(i,j)之值将等于c(i,j,n)值。此时图15-5可改写成程序15-9的C++代码。程序15-9中还利用了程序12-1中定义的AdjacencyWDigraph类。函数AllPairs在c中返回最短路径的长度。若i到j无通路，则c[i][j]被赋值为NoEdge。函数AllPairs同时计算了kay[i][j]，其中kay[i][j]表示从i到j的最短路径中最大的k值。在后面将看到如何根据kay值来推断出从一个顶点到另一顶点的最短路径（见程序15-10中的函数OutputPath）。 程序15-9c和kay的计算 voidAdjacencyWDigraph {//所有点对的最短路径 //对于所有i和j，计算c[i][j]和kay[i][j] //初始化c[i][j]=c（i，j，0） for(intj=1;j<=n;j++){ c[i][j]=a[i][j]; kay[i][j]=0; for(i=1;i<=n;i++) //计算c[i][j]=c(i,j,k) Tt1=c[i][k]; Tt2=c[k][j]; Tt3=c[i][j]; if(t1!=NoEdge&&t2!=NoEdge&&(t3==NoEdge||t1+t2 c[i][j]=t1+t2; 程序15-10输出最短路径 voidoutputPath(int**kay,inti,intj) {//输出i到j的路径的实际代码 if(kay[i][j]==0)cout< else{outputPath(kay,i,kay[i][j]); outputPath(kay,kay[i][j],j);} voidOutputPath(T**c,int**kay,TNoEdge,inti,intj) {//输出从i到j的最短路径 if(c[i][j]==NoEdge){ cout<<"Thereisnopathfrom"< return;} cout<<"Thepathis"< cout< outputPath(kay,i,j); cout< 例3-17图15-6a给出某图的长度矩阵a，15-6b给出由程序15-9所计算出的c矩阵，15-6c为对应的kay值。根据15-6c中的kay值，可知从1到5的最短路径是从1到kay[1][5]=4的最短路径再加上从4到5的最短路径，因为kay[4][5]=0，所以从4到5的最短路径无中间顶点。从1到4的最短路径经过kay[1][4]=3。重复以上过程，最后可得1到5的最短路径为：1，2，3，4，5。 3.2.5网络的无交叉子集 在11.5.3节的交叉分布问题中，给定一个每边带n个针脚的布线通道和一个排列C。顶部的针脚i与底部的针脚Ci相连，其中1≤i≤n，数对(i,Ci)称为网组。总共有n个网组需连接或连通。假设有两个或更多的布线层，其中有一个为优先层，在优先层中可以使用更细的连线，其电阻也可能比其他层要小得多。布线时应尽可能在优先层中布设更多的网组。而剩下的其他网组将布设在其他层。当且仅当两个网组之间不交叉时，它们可布设在同一层。我们的任务是寻找一个最大无交叉子集（MaximumNoncrossingSubset，MNS)。在该集中，任意两个网组都不交叉。因(i,Ci)完全由i决定，因此可用i来指定(i,Ci)。 例3-18考察图15-7（对应于图10-17）。(1,8)和(2,7)（也即1号网组和2号网组）交叉，因而不能布设在同一层中。而(1,8)，(7,9)和(9,10)未交叉，因此可布设在同一层。但这3个网组并不能构成一个MNS，因为还有更大的不交叉子集。图10-17中给出的例子中，集合｛(4,2),(5,5),(7,9),(9,10)｝是一个含4个网组的MNS。 设MNS(i,j)代表一个MNS，其中所有的(u,Cu)满足u≤i，Cu≤j。令size(i,j)表示MNS(i,j)的大小(即网组的数目)。显然MNS(n,n)是对应于给定输入的MNS，而size(n,n)是它的大小。 例3-19对于图10-17中的例子，MNS(10,10)是我们要找的最终结果。如例3-18中所指出的，size(10,10)=4，因为(1,8)，(2,7)，(7,9)，(8,3)，(9,10)和(10,6)中要么顶部针脚编号比7大，要么底部针脚编号比6大，因此它们都不属于MNS(7,6)。因此只需考察剩下的4个网组是否属于MNS(7,6)，如图15-8所示。子集｛(3,4),(5,5)｝是大小为2的无交叉子集。没有大小为3的无交叉子集，因此size(7,6)=2。 当i=1时，(1,C1)是MNS(1,j)的唯一候选。仅当j≥C1时，这个网组才会是MNS(1,j)的一个成员. 下一步，考虑i＞1时的情况。若j＜Ci，则(i,Ci)不可能是MNS(i,j)的成员，所有属于MNS(i,j)的(u,Cu)都需满足u＜i且Cu＜j，因此：size(i,j)=size(i-1,j),j 若j≥Ci，则(i,Ci)可能在也可能不在MNS(i,j)内。若(i,Ci)在MNS(i,j)内，则在MNS(i,j)中不会有这样的(u,Cu)：u＜i且Cu＞Ci，因为这个网组必与(i,Ci)相交。因此MNS(i,j)中的其他所有成员都必须满足条件u＜i且Cu＜Ci。在MNS(i,j)中这样的网组共有Mi-1,Ci-1个。若(i,Ci)不在MNS(i,j)中，则MNS(i,j)中的所有(u,Cu)必须满足u＜i；因此size(i,j)=size(i-1,j)。虽然不能确定(i,Ci)是否在MNS(i,j)中，但我们可以根据获取更大MNS的原则来作出选择。因此：size(i,j)=max{size(i-1,j),size(i-1,Ci-1)+1},j≥Ci（15-8） 虽然从（15-6）式到（15-8）式可用递归法求解，但从前面的例子可以看出，即使避免了重复计算，动态规划递归算法的效率也不够高，因此只考虑迭代方法。在迭代过程中先用式（15-6）计算出size(1,j)，然后再用式（15-7）和（15-8）按i=2,3,.,n的顺序计算size(i,j)，最后再用Traceback来得到MNS(n,n)中的所有网组。 例3-20图15-9给出了图15-7对应的size(i,j)值。因size(10,10)=4，可知MNS含4个网组。为求得这4个网组，先从size(10,10)入手。可用（15-8）式算出size(10,10)。根据式（15-8）时的产生原因可知size(10,10)=size(9,10)，因此现在要求MNS(9,10)。由于MNS(10,10)≠size(8,10)，因此MNS(9,10)中必包含9号网组。MNS(9,10)中剩下的网组组成MNS(8,C9-1)=MNS(8,9)。由MNS(8,9)=MNS(7,9)知，8号网组可以被排除。接下来要求MNS(7,9)，因为size(7,9)≠size(6,9)，所以MNS中必含7号网组。MNS(7,9)中余下的网组组成MNS(6,C7-1)=MNS(6,8)。根据size(6,8)=size(5,8)可排除6号网组。按同样的方法，5号网组，3号网组加入MNS中，而4号网组等其他网组被排除。因此回溯过程所得到的大小为4的MNS为{3,5,7,9}。 注意到在回溯过程中未用到size(10,j)(j≠10)，因此不必计算这些值。 为(n2)。 程序15-11寻找最大无交叉子集 voidMNS(intC[],intn,int**size) {//对于所有的i和j，计算size[i][j] //初始化size[1][*] for(intj=0;j size[1][j]=0; for(j=C[1];j<=n;j++) size[1][j]=1; //计算size[i][*],1 for(inti=2;i for(intj=0;j size[i][j]=size[i-1][j]; for(j=C[i];j<=n;j++) size[i][j]=max(size[i-1][j],size[i-1][C[i]-1]+1); size[n][n]=max(size[n-1][n],size[n-1][C[n]-1]+1); voidTraceback(intC[],int**size,intn,intNet[],int&m) {//在Net[0:m-1]中返回MMS intj=n;//所允许的底部最大针脚编号 m=0;//网组的游标 for(inti=n;i>1;i--) //i号net在MNS中 if(size[i][j]!=size[i-1][j]){//在MNS中 Net[m++]=i; j=C[i]-1;} //1号网组在MNS中 if(j>=C[1]) Net[m++]=1;//在MNS中 3.2.6元件折叠 实际上，在元件折叠问题中，还需考虑连接两个栈的线路所需的附加空间。如，在图15-10b中C5和C6间的线路因C6为折叠点而弯曲。这些线路要求在C5和C6之下留有垂直空间，以便能从栈1连到栈2。令ri为Ci是折叠点时所需的高度。栈1所需的高度为5i=1hi+r6，栈2所需高度为8i=6hi+r6+r9。 在标准单元设计中，电路首先被设计成为具有相同高度的符合线性顺序的元件排列。假设此线性顺序中的元件为C1，.，Cn，下一步元件被折叠成如图15-11所示的相同宽度的行。在此图中，12个标准单元折叠成四个等宽行。折叠点是C4，C6和C11。在相邻标准单元行之间，使用布线通道来连接不同的行。折叠点决定了所需布线通道的高度。设li表示当Ci为折叠点时所需的通道高度。在图15-11的例子中，布线通道1的高度为l4，通道2的高度为l6，通道3的高度为l11。位－片栈折叠和标准单元折叠都会引出一系列的问题，这些问题可用动态规划方法来解决。 1.等宽位－片元件折叠 定义r1=rn+1=0。由元件Ci至Cj构成的栈的高度要求为jk=ilk+ri+rj+1。设一个位－片设计中所有元件有相同宽度W。首先考察在折叠矩形的高度H给定的情况下，如何缩小其宽度。设Wi 为将元件Ci到Cn折叠到高为H的矩形时的最小宽度。若折叠不能实现（如当ri+hi＞H时），取Wi=∞。注意到W1可能是所有n个元件的最佳折叠宽度。 当折叠Ci到Cn时，需要确定折叠点。现假定折叠点是按栈左到栈右的顺序来取定的。若第一点定为Ck+1，则Ci到Ck在第一个栈中。为了得到最小宽度，从Ck+1到Cn的折叠必须用最优化方法，因此又将用到最优原理，可用动态规划方法来解决此问题。当第一个折叠点k+1已知时，可得到以下公式： Wi=w+Wk+1（15-9） 由于不知道第一个折叠点，因此需要尝试所有可行的折叠点，并选择满足（15-9）式的折叠点。令hsum(i,k)=kj=ihj。因k+1是一个可行的折叠点，因此hsum(i,k)+ri+rk+1一定不会超过H。 根据上述分析，可得到以下动态规划递归式： 现在来考察另外一个有关等宽元件的折叠问题：折叠后矩形的宽度W已知，需要尽量减小其高度。因每个折叠矩形宽为w，因此折叠后栈的最大数量为s=W/w。令Hi,j为Ci,.,Cn折叠成一宽度为jw的矩形后的最小高度，H1,s则是所有元件折叠后的最小高度。当j=1时，不允许任何折叠，因此：Hi,1=hsum(i,n)+ri,1≤i≤n 另外，当i=n时，仅有一个元件，也不可能折叠，因此：Hn,j=hn+rn,1≤j≤s 在其他情况下，都可以进行元件折叠。如果第一个折叠点为k+1，则第一个栈的高度为 hsum(i,k)+ri+rk+1。其他元件必须以至多(j-1)*w的宽度折叠。为保证该折叠的最优性，其他元件也需以最小高度进行折叠. 因为第一个折叠点未知，因此必须尝试所有可能的折叠点，然后从中找出一个使式（15-11）的右侧取最小值的点，该点成为第一个折叠点。 2.变宽位－片元件的折叠 首先考察折叠矩形的高度H已定，欲求最小的折叠宽度的情况。令Wi如式（15-10）所示，按照与（15-10）式相同的推导过程，可得： Wi=min{wmin(i,k)+Wk+1|hsum(i,k)+ri+rk+1≤H,i≤k≤n}（15-13） 3.标准单元折叠 用wi定义单元Ci的宽度。每个单元的高度为h。当标准单元行的宽度W固定不变时，通过减少折叠高度，可以相应地减少折叠面积。考察Ci到Cn的最小高度折叠。设第一个折叠点是Cs+1。从元件Cs+1到Cn的折叠必须使用最小高度，否则，可使用更小的高度来折叠Cs+1到Cn，从而得到更小的折叠高度。所以这里仍可使用最优原理和动态规划方法。 令Hi,s为Ci到Cn折叠成宽为W的矩形时的最小高度，其中第一个折叠点为Cs+1。令wsum(i,s)=sj=iwj。可假定没有宽度超过W的元件，否则不可能进行折叠。对于Hn,n因为只有一个元件，不存在连线问题，因此Hn,n=h。对于Hi,s（1≤i＜s≤n）注意到如果wsum(i,s)＞W，不可能实现折叠。若wsum(i,s)≤W，元件Ci和Cj+1在相同的标准单元行中，该行下方布线通道的高度为ls+1（定义ln+1=0）。因而：Hi,s=Hi+1,k（15-14） 当i=s＜n时，第一个标准单元行只包含Ci。该行的高度为h且该行下方布线通道的高度为li+1。因Ci+1到Cn单元的折叠是最优的. 为了寻找最小高度折叠，首先使用式（15-14）和（15-15）来确定Hi,s（1≤i≤s≤n）。最小高度折叠的高度为min｛H1,s｝。可以使用回溯过程来确定最小高度折叠中的折叠点。