首页 > 代码库 > 动态规划之矩阵链乘
动态规划之矩阵链乘
问题提出：
对于如下矩阵：
其中各矩阵A[i]下标为
![clip_image001[18]](/d/file/p/2024/07/31/58ea59e7e7df15ddc2f6c88b06f65d92.png "clip_image001[18]"){kind=small}
计算其乘积的结果，以及我们需要计算其最小标量乘法次数。
问题分析：
首先我们需要明确的是何为标量：标量即为没有方向的量，而有方向的量即为矢量。（严谨的定义自己百度去）
那么标量乘法就变成了最基本的数字相乘。
其次对于两个矩阵相乘，需满足下示公式所示的形式：（左边矩阵的列数与右边矩阵的行数必须一致）
上述条件可从矩阵相乘的定义中看出：
在计算机，我们可以用一个二维数组来表示矩阵。
一个m行n列的矩阵与一个n行p列的矩阵相乘，会得到一个m行p列的矩阵，具体步骤如下：
其中第i行第j列位置上的数为第一个矩阵的第i行上的n个数与第二个矩阵第j列上的n个数对应相乘后所得的n个乘积之和。
下图展示了一个简单的例子：
从上图我们可以得出两个矩阵相乘城的次数，即为取得结果矩阵每个点所要做的标量乘法次数之和
而每个点的标量乘法次数为左边矩阵的列数n，结果矩阵有m*p个点，那么总数即为m*n*p；
另外我们也应该知道矩阵乘法符合乘法结合律，但不符合交换律；也就是说，我们可以为矩阵链乘设置括号，括号放置的位置不同，产生标量乘法的总量也是不一样的哦，这一点我们必须明确，对这点不清楚的，可以自己举个例子试试看。如果没有这点，那这题就没意义了。
前面科普了一大堆知识，现在进正题了：
对于矩阵链乘，我们可知道不同的括号化方案，会产生不同的计算代价；
算法导论书中指出了求完全括号化方案，所谓矩阵乘积链为完全括号化需满足如下性质：
要么是单一矩阵，要么是两个完全括号化的矩阵乘积链的积，且已加外括号。说实话这个定义让人有点摸不着头脑，待我仔细解释下：
首先：该乘积链中只有一个矩阵，即单一矩阵，该乘积链必是完全括号化，且该矩阵外部不用加括号。
再者：必须是两个括号化的矩阵链的乘积，我们可以把两个括号化的矩阵链看作两个结果矩阵（即两个单一矩阵），该两个矩阵相乘，并加上外括号把他们括住，即这种形式也叫着完全括号化的。
啰嗦了这么多，感觉简单的说法就是：单一矩阵（不加括号）是完全括号化，两个矩阵（单一矩阵）加上外括号也是完全括号化的，完全括号化的乘积链可以当作单一矩阵处理。
下面上一道课后的证明题：
对n个元素的表达式进行完全括号化时，恰好需要n-1对括号。
证明：
注：采用数学归纳法
已知仅有一个元素时，不需要括号，即0对括号。
假设k个元素时，需要k-1对括号；
当有k+1个元素时，我们可以看作k个元素的完全括号化的乘积链外加1个新加的元素。已知我们可以把完全括号化的乘积链看作单一矩阵，那么现有两个单一矩阵，那么需再加1对括号，才能再次完全括号化，而由假设知，k个元素的完全括号化需k-1对括号，那么k+1个元素的完全括号化方案就需要k-1+1对括号；
命题得证！
额，又科普了部分知识，只不过这些科普对后续理解问题都是有用的。
应用动态规划的方法的步骤：（出自书本）
1、刻画一个最优解的结构特征
2、递归地定义最优解的值
3、计算最优解的值，通常采用自底向上的方法
4、利用计算出的信息构造一个最优解的值
第一步：刻画一个最优解的结构特征
对于最开始提出的问题，我们看作为矩阵链乘设计最优的括号化方案。而设置括号的最终目的就是为了打乱计算顺序，通过结合律来找到最小计算代价，就比如下面一个问题33*4*25，是个正常人就会先计算后面两项，因为这样好计算，同样计算矩阵，我们也希望计算的少，计算的快啊，设置括号就变的尤为重要了。
对于完全括号化，每一个括号的作用是将一个链乘转化成两个完全括号化的链乘的积，前面有提到。
因此我最优解的结构特征就如下了：
对于第m个矩阵到第n个矩阵，他的最优解存在于对它的一次划分中，它的划分有n-m种情况，我们对这些情况中代价取最小，不就获得了最优解
而对于每种情况，我们需要划分后的两部分矩阵链都是最优解，乘积后才可能最小啊，（这部分大家可以用“剪切--粘贴”去反证），这样问题就转换成了两个子部分矩阵链的最优括号化方案的问题，以此不断的递归下去。
第二步：递归地定义最优解的值
我们的划分点k可从m到n-1，代价如下所示，每次的代价如下，而最优解只有一个，即它们的最小值：
其一般化的公式如下：
从这里可以看出，这是一个递归的问题，看着很类似上篇所述的钢条切割问题啊，仔细比对比对吧。
若采用递归策略，可能会有很多的重复子问题，无疑提高了计算的复杂度。也可以采取带备忘的递归来降低复杂度。
第三步：计算最优解的值，通常采用自底向上的方法
另外从公式中可以看出，上述计算结果要依赖于m,n（m<n）的所有可能组合。区间[m，n]组合要依赖于它们之间的子区间组合。所以我们采取自下而上的策略，逐渐扩大m,n之间的区间长度，最终算到m=0,n=n，从而获得最优的括号化方案，
具体方法为：matrix_chain_mutiply.cpp中的void matrix_divide_option(int *num,int** cost,int** point,int length)方法。
第四步、利用计算出的信息构造一个最优解的值
具体解析见matrix_chain_mutiply.cpp最后两个方法
代码如下：
matrix.h文件：
#include<iostream>#ifndef MATRIX_H#define MATRIX_H/*矩阵类row为其行数col为该矩阵的列数，value表示一个二维数组，主要存储该矩阵各点的值*/class matrix{ public: matrix(); //默认的构造函数 matrix(int row,int col); //给定行数和列数，构造一个矩阵 virtual ~matrix(); //析构函数 matrix(const matrix& other); //拷贝构造函数 matrix operator*(const matrix& a); //矩阵的乘法 matrix& operator=(const matrix&); //重写赋值函数 matrix init(); //初始化矩阵 friend std::ostream& operator << (std::ostream&,const matrix &); //重写输出 protected: private: double **value; //主要用于存储矩阵各点的值 int col; //列数 int row; //行数};#endif // MATRIX_H
matrix.cpp
#include<iostream>#include "matrix.h"#include<cstdlib>#include<ctime>matrix::matrix() //默认构造函数{ col = row = 0; value = nullptr;}matrix::matrix(int row, int col) //根据给定行，列数来初始化一个矩阵，并分配相应的内存{ this->row = row; this->col = col; value = new double*[row]; for(int i=0;i<row;i++) { value[i] = new double[col]; }}matrix::matrix(const matrix& other) //复制构造函数，注意此处要用matrix&，重写复制构造函数为了防止值复制造成{ row = other.row; col = other.col; value = new double*[row]; //此处使用new分配内存 for(int i=0;i<row;i++) { value[i] = new double[col]; } for(int i=0;i<row;i++) { for(int j=0;j<col;j++) { value[i][j]=(other.value)[i][j]; } }}matrix::~matrix(){ for(int i=0;i<row;i++) //因为构造函数中new分配有内存，所以在析构函数中要对内存进行释放。 { delete[] value[i]; } delete[] value;}matrix& matrix::operator=(const matrix& other) //重写了赋值函数{ if(this == &other) //防止a=a的情况出现 { return *this; } for(int i=0;i<row;i++) //先释放原有区域的内存，从此处我们可以上面的判断的重要性 { //若没有上述判断，a=a在赋值的同时，也将内存释放了，造成数据丢失 delete[] value[i]; } delete[] value; row = other.row; //将成员变量进行复制 col = other.col; value = new double*[row]; //重新分配内存 for(int i=0;i<row;i++) { value[i] = new double[col]; } for(int i=0;i<row;i++) //对二维数组进行复制 { for(int j=0;j<col;j++) { value[i][j]=(other.value)[i][j]; } } return *this;}matrix matrix::operator*(const matrix& n) //矩阵的乘法{ matrix result(row,n.col); for(int i=0;i<row;i++) { for(int j=0;j<n.col;j++) { (result.value)[i][j] = 0; for(int k=0;k<n.row;k++) { (result.value)[i][j] += value[i][k]*(n.value)[k][j]; } } } return result;}matrix matrix::init() //初始化矩阵{ matrix result(row,col); srand((unsigned)time(NULL)); for(int i=0;i<row;i++) { for(int j=0;j<col;j++) { (result.value)[i][j] = rand()%100; } } return result;}std::ostream& operator <<(std::ostream& os,const matrix& m)//将矩阵输出{ for(int i=0;i<m.row;i++) { for(int j=0;j<m.col;j++) { os << (m.value)[i][j] << ‘\t‘; } os << std::endl; } return os;}
matrix_chain_mutiply.cpp(核心)
#include<iostream>#include "matrix.h"const int SIZE = 5; //SIZE的大小决定了矩阵的数量（SIZE-1），可自行改变const int max_size = 2147483647;void matrix_divide_option(int *,int **,int **,int); //确定最优的括号化方案void show(int **); //输出二维数组void print(int**,int,int); //输出矩阵链乘的括号化方案matrix matrix_mutiply(int**,int,int,matrix*); //获得矩阵链乘的最终结果矩阵int main(){ using namespace std; cout << "input " << SIZE << " number(press enter to end your input)：" << endl; int *num = new int[SIZE]; //用于存储矩阵A[i]下标为num[i]*num[i+1]; int i=0; cin >> num[i]; //获得矩阵的行列数 while(cin.get() != ‘\n‘) { i++; cin >> num[i]; //读取数组维数； } cout << endl << "what you input is:" << endl; for(int j=0;j<SIZE;j++) { cout << num[j] << "-----"; //展示输入的是否正确 } cout << endl << endl; matrix* test = new matrix[SIZE-1]; //为矩阵数组分配内存 for(int j=1;j<SIZE;j++) { matrix m(num[j-1],num[j]); //初始化矩阵大小 test[j-1] = m.init(); //初始化该矩阵各点的值 } cout << "矩阵信息如下：" << endl; for(int j=0;j<SIZE-1;j++) { cout << test[j] << endl; //输出矩阵信息 } int **cost = new int*[SIZE-1]; //存储矩阵相乘时，各标量乘法次数，即代价表 int **point = new int*[SIZE-1]; //存储矩阵括号化方案时的分割点 for(int j=0;j<SIZE-1;j++) { cost[j] = new int[SIZE-1]; //分配内存 point[j] = new int[SIZE-1]; } matrix_divide_option(num,cost,point,SIZE-1); //获取划分方案，获得最小的代价 cout << "各区间的计算代价：" << endl; show(cost); cout << endl; cout << "区间划分点所在的位置：" << endl; show(point); cout << endl; cout << "输出完全括号化方案：" << endl; print(point,0,SIZE-2); cout << endl << endl; cout << "输出结果矩阵：" << endl; matrix result = matrix_mutiply(point,0,SIZE-2,test); cout << result; for(int j=0;j<SIZE-1;j++) { delete[] cost[j]; delete[] point[j]; } delete[] cost; delete[] point; delete[] test; delete[] num; return 0;}/*此方法主要获得计算代价，以及括号化方案的划分点为此问题解决的核心方法。num存储的矩阵的下标，num[i-1]为第i个矩阵的行数，num[i]表示第i个矩阵的列数因为数组存储时，最初始的下标为0，所以下述的定义就变的有点不好理解！仔细想想也容易明白。cost----二维数组，cost[i][j]代表的就是第（i+1）个矩阵到第（j+1）个矩阵的链乘的最小计算代价point----二维数组，point[i][j]代表的就是第（i+1）个矩阵到第（j+1）个矩阵的链乘的最佳划分处，假如为k，则k就是第k+1个矩阵，（i+1）到(k+1)矩阵属于左边的完全括号化矩阵，其余为右侧完全括号化矩阵length代表的是矩阵的个数，即num的长度-1*/void matrix_divide_option(int *num,int** cost,int** point,int length){ //首先将代价表中的值重置 //下标i到下标j之间的矩阵链乘的最小值 for(int i=0;i<length;i++) { for(int j=0;j<length;j++) { if(i>=j) { cost[i][j] = 0; //单个矩阵不存在乘积，代价0，同时也将不合理的下标组合置为0 } else { cost[i][j] = max_size; //先将各合理的组合代价设为最大，以便后续比较使用。 } point[i][j] = -1; //将所有组合的划分点置为-1，以便后续记录。 } }/* l代表矩阵链中矩阵的个数，我们这里从2个矩阵起步.由公式可知，我们要算遍所有的合理组合，同时大区间的组合 还要依赖于其子区间的合理组合，因为完全括号化的性质，最终都归根于2个矩阵的乘积，所以我们从2个矩阵起步， 逐步向上计算，最终就能计算（0，n）的组合的最优括号化的代价。*/ for(int l=2;l<=length;l++) { for(int i=0;i<=length-l;i++) //此处要注意的i+l-1不能超过数组的个数 { int j = i+l-1; //矩阵个数要为l，则j-i需满足j-i+1=l，所以有该等式 int q = 0; for(int k=i;k<j;k++) //注意k要在区间内移动，利用下述if条件，获得所有可能值中的最小值， { q = cost[i][k] + cost[k+1][j] + num[i]*num[k+1]*num[j+1]; if(q < cost[i][j]) { cost[i][j] = q; point[i][j] = k; } } } }}void show(int **p){ using namespace std; for(int i=0;i<SIZE-1;i++) { for(int j=0;j<SIZE-1;j++) { cout << p[i][j] << ‘\t‘; } cout << endl; }}void print(int** point,int i,int j){ if(i == j) { std::cout << "A" << i+1; } else { std::cout << "("; print(point,i,point[i][j]); //输出划分点左侧的完全括号化表达式 print(point,point[i][j]+1,j); //输出划分点右侧的完全括号化表达式 std::cout << ")"; }}/*矩阵链乘的方法与括号化方案表达式输出的方法是同样的道理*/matrix matrix_mutiply(int** point,int i,int j,matrix* m){ if(i==j) return m[i]; else { matrix m1 = matrix_mutiply(point,i,point[i][j],m); matrix m2 = matrix_mutiply(point,point[i][j]+1,j,m); return m1*m2; }}
该问题的解决代码在上述三个文件中。
同时上述代码也存在一定的问题，关于内存是否泄漏以及代码是否符合规范，本人不是非常熟悉。望各位大神，看后觉得有什么问题，直接回复，我也希望提高，毕竟code能力有待提高。
以上问题也属于典型的动态规划问题，也属于递归向迭代的转换。
总结：
应用动态规划的方法的步骤：
1、刻画一个最优解的结构特征
2、递归地定义最优解的值
3、计算最优解的值，通常采用自底向上的方法
4、利用计算出的信息构造一个最优解的值
动态规划之矩阵链乘