动态规划之矩阵链乘 C++实现
原理
在上一次的文章当中,作者讲解了什么是动态规划,以及动态规划的一个举例应用,这次,我们来看看如何运用动态规划来解决矩阵链乘问题。
关于矩阵的乘法,运用如下公式:
C=A×B
其中
cij=∑k=1naikbkj
并不是任意两个矩阵都能相乘,必须满足
A 的列数
na 等于
B 的行数
mb ,即有
Anama,Bnbmb 且
na=mb 。相乘之后,矩阵
C 的行列分别为
ma、nb 。
如果不只是两个矩阵相乘,而是 n 个矩阵相乘呢?同样可以用上面的公式算出结果,但是存在着一个乘积次数的问题,若有 Apq,Bqr 且 na=mb ,则其相乘次数为 pqr 。假设有三个矩阵 ⟨A1,A2,A3⟩ ,其规模分别为 10×100、100×5、5×50 。如果按照 ((A1A2)A3) 的顺序计算,需要做 10∗100∗5+10∗5∗50=7500 次标量乘法。而如果按照 (A1(A2A3)) 的顺序计算,需要做 100∗5∗50+10∗100∗50=75000 次标量乘法。因此,按照第一种顺序计算矩阵的乘法要比第二种顺序块10倍。
矩阵链问题描述如下:给定 n 个矩阵的链 ⟨A1,A2…An⟩ ,矩阵 Ai 的规模为 pi−1×pi(1≤i≤n) ,求完全括号化方案,使得计算乘积 A1A2…An 所需标量乘法次数最少。
令 m[i,j] 表示计算矩阵 Aij 所需标量乘法次数的最小值,那么 A1,A2…An 相乘,最小代价括号化方案的递归求解公式变为:
m[i,j]=⎧⎩⎨0maxi≤k<j{m[i,k]+m[k+1,j]+pi−1pkpj}(i=j)(i<j)
以下附上自底向上方法 O(n3) ,递归版 Ω(2n) 与带备忘的自顶向下法 O(n3) 的实现。
源代码:
#include <iostream>
#include <utility>
#include <vector>
using namespace std;
vector<int> temp_Vec{ 30,35,15,5,10,20,25 }; //A1(30X35), A2(35X15.).....A6(20X25)
//Bottom-up.
pair<vector<vector<int>>, vector<vector<int>>> Matrix_Chain_Order(const vector<int> &temp_Vec) {
auto temp_n = temp_Vec.size() - 1;
vector<vector<int>> temp_VecM, temp_VecS;
temp_VecM.resize(temp_n + 1);
temp_VecS.resize(temp_n + 1);
for(auto &i : temp_VecM) {
i.resize(temp_n + 1);
}
for (auto &i : temp_VecS) {
i.resize(temp_n + 1);
}
for(auto i = 1; i <= temp_n; ++i) {
temp_VecM[i][i] = 0;
}
for(auto l = 2; l <= temp_n; ++l) {
for (auto i = 1; i <= temp_n - l + 1; ++i) {
auto j = i + l - 1;
temp_VecM[i][j] = INT_MAX;
for(auto k = i; k <= j - 1; ++k) {
auto temp_q = temp_VecM[i][k] + temp_VecM[k + 1][j] + temp_Vec[i - 1] * temp_Vec[k] * temp_Vec[j];
if(temp_q < temp_VecM[i][j]) {
temp_VecM[i][j] = temp_q;
temp_VecS[i][j] = k;
}
}
}
}
return make_pair(temp_VecM, temp_VecS);
}
//Memoized of part1.
int Lookup_Chain(vector<vector<int>> &temp_VecM, const vector<int> &temp_Vec, const int &i, const int &j) {
if(temp_VecM[i][j] < INT_MAX) {
return temp_VecM[i][j];
}
if(i == j) {
temp_VecM[i][j] = 0;
}
else {
for(auto k = i; k <= j - 1; ++k) {
auto temp_q = Lookup_Chain(temp_VecM, temp_Vec, i, k) + Lookup_Chain(temp_VecM, temp_Vec, k + 1, j) + temp_Vec[i - 1] * temp_Vec[k] * temp_Vec[j];
if(temp_q < temp_VecM[i][j]) {
temp_VecM[i][j] = temp_q;
}
}
}
return temp_VecM[i][j];
}
//Memoized of part2.
int Memoized_Matrix_Chain(const vector<int> &temp_Vec) {
auto temp_n = temp_Vec.size() - 1;
vector<vector<int>> temp_VecM;
temp_VecM.resize(temp_n + 1);
for (auto &i : temp_VecM) {
i.resize(temp_n + 1);
}
for(auto &i : temp_VecM) {
for(auto &j : i) {
j = INT_MAX;
}
}
return Lookup_Chain(temp_VecM, temp_Vec, 1, temp_n);
}
//Recusive.
int Recursive_Matrix_Chain(const vector<int> &temp_Vec, const int &i, const int &j) {
if (i == j) {
return 0;
}
auto temp_n = temp_Vec.size() - 1;
vector<vector<int>> temp_VecM;
temp_VecM.resize(temp_n + 1);
for (auto &x : temp_VecM) {
x.resize(temp_n + 1);
}
for (auto &x : temp_VecM) {
for (auto &y : x) {
y = INT_MAX;
}
}
for (auto k = i; k <= j - 1; ++k) {
auto temp_q = Recursive_Matrix_Chain(temp_Vec, i, k) + Recursive_Matrix_Chain(temp_Vec, k + 1, j) + temp_Vec[i - 1] * temp_Vec[k] * temp_Vec[j];
if (temp_q < temp_VecM[i][j]) {
temp_VecM[i][j] = temp_q;
}
}
return temp_VecM[i][j];
}
//Print optimal bracket scheme.
void Print_Optimal_Parens(const vector<vector<int>> &temp_VecS, const int &i, const int &j) {
if(i == j) {
cout << "A" << i << " ";
}
else {
cout << "(" << " ";
Print_Optimal_Parens(temp_VecS, i, temp_VecS[i][j]);
Print_Optimal_Parens(temp_VecS, temp_VecS[i][j] + 1, j);
cout << ")" << " ";
}
}
int main() {
auto temp_Pair = Matrix_Chain_Order(temp_Vec);
auto temp_VecM = temp_Pair.first, temp_VecS = temp_Pair.second;
cout << temp_VecM[2][5] << endl; //A2 to A5 by bottom-up.
Print_Optimal_Parens(temp_VecS, 2, 5); //A2 to A5 by bottom-up.
cout << endl << endl;
for(auto &i : temp_VecM) { //All of the value.
for(auto &j : i) {
cout << j << " ";
}
cout << endl;
}
cout << endl << Recursive_Matrix_Chain(temp_Vec, 1, 6) << endl; // A1 to A6 by recursive.
cout << endl << Memoized_Matrix_Chain(temp_Vec) << endl; //A1 to A6 by memoized.
return 0;
}