前言:今天接着学习动态规划算法,学习如何用动态规划来分析解决矩阵链乘问题。首先回顾一下矩阵乘法运算法,并给出C++语言实现过程。然后采用动态规划算法分析矩阵链乘问题并给出C语言实现过程。
#include <iostream>
using namespace std;
#define A_ROWS 3 /*矩阵A的行为3*/
#define A_COLUMNS 2 /*矩阵A的列为2*/
#define B_ROWS 2 /*矩阵B的行为2*/
#define B_COLUMNS 3 /*矩阵B的列为3*/
void matrix_multiply(int A[A_ROWS][A_COLUMNS],int B[B_ROWS][B_COLUMNS],int C[A_ROWS][B_COLUMNS]);
int main()
{ //初始化矩阵A、B、C的值
int A[A_ROWS][A_COLUMNS]={1,0,
1,2,
1,1};
int B[B_ROWS][B_COLUMNS]={1,1,2,
2,1,2};
int C[A_ROWS][B_COLUMNS]={0};
matrix_multiply(A,B,C); /*C=A*B*/
//输出C矩阵的元素
for(int i=0;i<A_ROWS;i++)
{
for(int j=0;j<B_COLUMNS;j++)
cout<<C[i][j]<<" ";
cout<<endl;
}
return 0;
}
void matrix_multiply(int A[A_ROWS][A_COLUMNS],int B[B_ROWS][B_COLUMNS],int C[A_ROWS][B_COLUMNS])
{
if(A_COLUMNS!=B_ROWS) /*不满足矩阵乘法的原则*/
cout<<"error:incompatible dimensions."<<endl;
else
{
int i,j,k;
for(i=0;i<A_ROWS;i++) //
for(j=0;j<B_COLUMNS;j++)
{
C[i][j]=0;
for(k=0;k<A_COLUMNS;k++)
C[i][j]+=A[i][k]*B[k][j]; /*将A的每一行的每一列与B的每一列的每一行的乘积求和*/
}
}
}
2、矩阵链乘问题描述
给定n个矩阵构成的一个链<A1,A2,A3,.......An>,其中i=1,2,...n,矩阵A的维数为pi-1pi,对乘积 A1A2...An 以一种最小化标量乘法次数的方式进行加全部括号。
注意:在矩阵链乘问题中,实际上并没有把矩阵相乘,目的是确定一个具有最小代价的矩阵相乘顺序。找出这样一个结合顺序使得相乘的代价最低。
3、动态规划分析过程
1)最优加全部括号的结构
动态规划第一步是寻找一个最优的子结构。假设现在要计算AiAi+1....Aj的值,计算Ai...j过程当中肯定会存在某个k值(i<=k<j)将Ai...j分成两部分,使得Ai...j的计算量最小。分成两个子问题Ai...k和Ak+1...j,需要继续递归寻找这两个子问题的最优解。
有分析可以到最优子结构为:假设AiAi+1....Aj的一个最优加全括号把乘积在Ak和Ak+1之间分开,则Ai..k和Ak+1..j也都是最优加全括号的。
2)一个递归解
设m[i,j]为计算机矩阵Ai...j所需的标量乘法运算次数的最小值,对此计算A1..n的最小代价就是m[1,n]。现在需要来递归定义m[i,j],分两种情况进行讨论如下:
当i==j时:m[i,j] = 0,(此时只包含一个矩阵)
当i<j 时:从步骤1中需要寻找一个k(i≤k<j)值,使得m[i,j] =min{m[i,k]+m[k+1,j]+pi-1pkpj} (i≤k<j)。
3)计算最优代价
虽然给出了递归解的过程,但是在实现的时候不采用递归实现,而是借助辅助空间,使用自底向上的表格进行实现。设矩阵Ai的维数为pi-1pi,i=1,2.....n。输入序列为:p=<p0,p1,...pn>,length[p] = n+1。使用m[n][n]保存m[i,j]的代价,s[n][n]保存计算m[i,j]时取得最优代价处k的值,最后可以用s中的记录构造一个最优解。书中给出了计算过程的伪代码,摘录如下:
MAXTRIX_CHAIN_ORDER(p)
n=length[p]-1; /*矩阵的最大下标*/
for i=1 to n
do m[i][i]=0; /*对角线元素*/
for t=2 to n /*t 是矩阵链的长度*/
do for i=1 to n-t+1
j=i+t-1 /*i和j的距离为t-1*/
m[i][j]=MAXLIMIT; /*初始化为最大值*/
for k=i to j-1 /*划分的位置逐渐的后移*/
q=m[i][k]+m[k+1][i]+qi-1qkqj;
if q<m[i][j]; /*比较,寻找最小值*/
then m[i][j]=q; /*保存当前矩阵链的乘法次数的最小值*/
s[i][j]=k; /*保存当前矩阵链的最优切割位置*/
return m and s; /*返回最优次数和最优切割位置*/
MATRIX_CHAIN_ORDER具有循环嵌套,深度为3层,运行时间为O(n3)。如果采用递归进行实现,则需要指数级时间Ω(2n),因为中间有些重复计算。递归是完全按照第二步得到的递归公式进行计算,递归实现如下所示:
int recursive_matrix_chain(int *p,int i,int j,int m[N+1][N+1],int s[N+1][N+1]) //矩阵链i-j,m乘法次数最小值,s最优的切割位置 { if(i==j) m[i][j]=0; else { int k; m[i][j]=MAXVALUE; for(k=i;k<j;k++) { int temp=recursive_matrix_chain(p,i,k,m,s)+recursive_matrix_chain(p,k+1,j,m,s)+p[i-1]*p[k]*p[j]; if(temp<m[i][j]) { m[i][j]=temp; s[i][j]=k; } } } return m[i][j]; }
对递归算计的改进,可以引入备忘录,采用自顶向下的策略,维护一个记录了子问题的表,控制结构像递归算法。完整程序如下所示:
4)构造一个最优解
第三步中已经计算出来最小代价,并保存了相关的记录信息。因此只需对s表格进行递归调用展开既可以得到一个最优解。书中给出了伪代码,摘录如下:
PRINT_OPTIMAL_PARENS(s,i,j) if i==j then print "Ai" else print "("; PRINT_OPTIMAL_PARENS(s,i,s[i][j]); PRINT_OPTIMAL_PARENS(s,s[i][j]+1,j); print ")";
4、编程实现
采用C++语言实现这个过程,现有矩阵A1(30×35)、A2(35×15)、A3(15×5)、A4(5×10)、A5(10×20)、A6(20×25),得到p=<30,35,15,5,10,20,25>。实现过程定义两个二维数组m和s,为了方便计算其第一行和第一列都忽略,行标和列标都是1开始。完整的程序如下所示:
#include <iostream> using namespace std; #define N 6 #define MAXVALUE 1000000 void matrix_chain_order(int *p,int len,int m[N+1][N+1],int s[N+1][N+1]); void print_optimal_parents(int s[N+1][N+1],int i,int j); int main() { int p[N+1]={30,35,15,5,10,20,25}; int m[N+1][N+1]={0}; int s[N+1][N+1]={0}; int i,j; matrix_chain_order(p,N+1,m,s); cout<<"m value is:"<<endl; for(i=1;i<=N;++i) { for(j=1;j<=N;++j) cout<<m[i][j]<<'\t'; cout<<endl; } cout<<"s value is:"<<endl; for(i=1;i<=N;++i) { for(j=1;j<=N;++j) cout<<s[i][j]<<'\t'; cout<<endl; } cout<<"The result is:"<<endl; print_optimal_parents(s,1,N); return 0; } void matrix_chain_order(int *p,int len,int m[N+1][N+1],int s[N+1][N+1]) { int i,j,k,t; for(i=0;i<=N;++i) m[i][i]=0; for(t=2;t<=N;t++) //当前矩阵链的长度 { for(i=1;i<=N-t+1;i++) //从第一矩阵开始计算,计算长度为t的最小代价 { j=i+t-1; //长度为t时候的最后一个元素 m[i][j]=MAXVALUE; //初始化为最大代价 for(k=i;k<=j-1;k++) { int temp=m[i][k]+m[k+1][j]+p[i-1]*p[k]*p[j]; if(temp<m[i][j]) { m[i][j]=temp; s[i][j]=k; } } } } } //s中存放着括号当前的位置 void print_optimal_parents(int s[N+1][N+1],int i,int j) { if(i==j) cout<<"A"<<i; else { cout<<"("; print_optimal_parents(s,i,s[i][j]); print_optimal_parents(s,s[i][j]+1,j); cout<<")"; } }
程序测试结果如下所示:
5、总结
动态规划解决问题关键是分析过程,难度在于如何发现其子问题的结构及子问题的递归解。这个需要多多思考,不是短时间内能明白。在实现过程中遇到问题就是数组,数组的下标问题是个比较麻烦的事情,如何能够过合理的去处理,需要一定的技巧。
冷静思考,勇敢面对,把握未来!
转载:http://www.cnblogs.com/Anker/archive/2013/03/10/2952475.html