动态规划之矩阵链乘 C++实现

动态规划之矩阵链乘 C++实现

原理

在上一次的文章当中，作者讲解了什么是动态规划，以及动态规划的一个举例应用，这次，我们来看看如何运用动态规划来解决矩阵链乘问题。

关于矩阵的乘法，运用如下公式：
$$C=A\times B$$
其中
$$c_{ij}=\sum _{k=1}^n{a}_{ik}{b}_{kj}$$
并不是任意两个矩阵都能相乘，必须满足$A$的列数$n_a$等于$B$的行数$m_b$，即有$A_{n_a{m}_a},B_{n_b{m}_b}$且$n_a=m_b$。相乘之后，矩阵$C$的行列分别为$m_a、n_b$。

如果不只是两个矩阵相乘，而是$n$个矩阵相乘呢？同样可以用上面的公式算出结果，但是存在着一个乘积次数的问题，若有$A_{pq},B_{qr}$且$n_a=m_b$，则其相乘次数为$pqr$。假设有三个矩阵$⟨A_1,A_2,A_3⟩$，其规模分别为$10\times 100、100\times 5、5\times 50$。如果按照$((A_1{A}_2)A_3)$的顺序计算，需要做$10\ast 100\ast 5+10\ast 5\ast 50=7500$次标量乘法。而如果按照$(A_1(A_2{A}_3))$的顺序计算，需要做$100\ast 5\ast 50+10\ast 100\ast 50=75000$次标量乘法。因此，按照第一种顺序计算矩阵的乘法要比第二种顺序块10倍。

矩阵链问题描述如下：给定$n$个矩阵的链$⟨A_1,A_2\dots A_n⟩$，矩阵$A_i$的规模为$p_i-1\times p_i(1\le i\le n)$，求完全括号化方案，使得计算乘积$A_1{A}_2\dots A_n$所需标量乘法次数最少。

令$m\left[i,j\right]$表示计算矩阵$A_{i}j$所需标量乘法次数的最小值，那么$ A_1,A_2\dots A_n$相乘，最小代价括号化方案的递归求解公式变为：
KaTeX parse error: Limit controls must follow a math operator at position 71: …ratorname*{max}\̲l̲i̲m̲i̲t̲s̲_{i\leq k< j}\l…

以下附上自底向上方法$O(n^3)$，递归版$\Omega (2^n)$与带备忘的自顶向下法$O(n^3)$的实现。

源代码：

#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;
}