Русская Википедия:Задача о порядке перемножения матриц
Задача о порядке перемножения матриц — классическая задача динамического программирования, в которой дана последовательность матриц <math> A_1, A_2, \dots, A_n </math> и требуется минимизировать количество скалярных операций для вычисления их произведения. Матрицы предполагаются совместимыми по отношению к матричному умножению (то есть количество столбцов <math> A_{i - 1}</math> совпадает с количеством строк <math> A_i </math> матрицы).
Подробное описание задачи
Произведение матриц — ассоциативная операция, которая принимает на вход две матрицы размером k×m и m×n и возвращает матрицу размером k×n, потратив на это kmn операций умножения[1].
Когда матрицы велики по одному измерению и малы по другому, количество скалярных операций может серьёзно зависеть от порядка перемножений матриц. Допустим, нам даны 3 матрицы <math> A_1, A_2, A_3 </math> размерами соответственно 10×100, 100×5 и 5×50. Существует 2 способа их перемножения (расстановки скобок): <math>((A_1A_2)A_3)</math> и <math>(A_1(A_2A_3))</math>. В первом случае нам потребуется 10·100·5 + 10·5·50 = 7500 скалярных умножений, а во втором случае 100·5·50 + 10·100·50 = 75000 умножений — разница налицо. Поэтому может быть выгоднее потратить некоторое время на предобработку, решив, в каком порядке лучше всего умножать, чем умножать сразу в лоб.
Таким образом, даны n матриц: <math>A_1: \, p_0 \times p_1</math>, <math>A_2: \, p_1 \times p_2</math>, …, <math>A_n: \, p_{n-1} \times p_{n}</math>. Требуется определить, в каком порядке перемножать их, чтобы количество операций умножения было минимальным.
Решение задачи
Разберём 2 способа решения задачи, чтобы показать, насколько выгодно динамическое программирование в данном случае.
Перебор всех вариантов расстановки скобок
Оценим, сколько же нужно перебрать вариантов расстановки. Обозначим через <math>P(n)</math> количество способов расстановки скобок в последовательности, состоящей из n матриц. Когда матрица одна, то расставлять нечего, вариант один. Если <math>n \ge 2</math>, то количество вариантов, которым можно расставить скобки является произведением количества вариантов, которым можно расставить скобки в составляющих результирующую матрицу произведениях (т.е. если <math> A_3 = A_1*A_2</math>, то количество вариантов, которым мы можем получить матрицу <math>A_3</math> равно произведению количества способов получить матрицу <math>A_1</math> на количество способов получить <math>A_2</math>). Разбиение на матрицы, <math>A_1</math> и <math>A_2</math> может производиться на границе k-й и (k + 1)-й матриц для <math>k = 1, 2, \dots, n-1</math>. Получаем рекуррентное соотношение: <math> P(n) = \left \{ \begin{array}{ll}
1,\ n=1 \\
\sum_{k=1}^{n-1} P(k)P(n-k),\ n>=2
\end{array}
\right.
</math>
Решением аналогичного рекуррентного соотношения является последовательность чисел Каталана, возрастающая как <math>\Omega \left ( \frac{4^n}{n^{1.5}} \right )</math>. Зависимость получается экспоненциальная, непригодная для практического применения в программах. Разберём более перспективный способ.
Динамическое программирование
Сведение задачи к подзадачам
Обозначим результат перемножения матриц <math>A_iA_{(i+1)}... A_j</math> через <math>A_{i..j}</math>, где i<=j. Если i<j, то существует такое k, которое разбивает <math>A_{i..j}</math> между матрицами <math>A_k</math> и <math>A_{k+1}</math>, i<=k<j. То есть для вычисления <math>A_{i..j}</math> надо сначала вычислить <math>A_{i..k}</math>, потом <math>A_{k+1..j}</math> и затем их перемножить. Выбор k является аналогом расстановки скобок между матрицами. Выбрав некоторое k мы свели задачу к двум аналогичным подзадачам для матриц <math>A_{i..k}</math> и <math>A_{k+1..j}</math>.
Рекурсивное решение
Обозначим через m[i, j] минимальное количество скалярных умножений для вычисления матрицы <math>A_{i..j}</math>. Получаем следующее рекуррентное соотношение: <math> m[i,j] = \left \{ \begin{array}{ll}
0, & i=j \\
\min\limits_{i \le k < j}\{m[i,k] + m[k+1,j] + p_{i-1}p_kp_j \}, & i < j
\end{array}
\right.
</math>
Объясняется оно просто: для того, чтобы найти произведение матриц <math>A_{i..j}</math> при i=j не нужно ничего делать — это и есть сама матрица <math>A_i</math>. При нетривиальном случае мы перебираем все точки разбиения матрицы <math>A_{i..j}</math> на матрицы <math>A_{i..k}</math> и <math>A_{k+1..j}</math>, ищем кол-во операций, необходимое чтобы их получить и затем перемножаем для получения матрицы <math>A_{i..j}</math>.(Оно будет равно кол-ву операций, потраченное на решение подзадач + стоимость умножения матриц <math>A_{i..k}A_{k+1..j}</math>). Считаем, что размеры матриц заданы в массиве <math>p</math> и размер матрицы <math>A_i</math> равен <math>p_{i-1} \times p_i</math>. Как обычно рекурсивный метод нельзя использовать напрямую — он будет экспоненциальным из-за большого кол-ва перекрывающихся подзадач.
Динамическое программирование
Будем запоминать в двумерном массиве m результаты вычислений для подзадач, чтобы избежать пересчета для уже вычислявшихся подзадач. После вычислений ответ будет в m[1,n](Сколько перемножений требуется для последовательности матриц от 1 до n — то есть ответ на поставленную задачу).Сложность алгоритма будет O<math>\left ( n^3 \right )</math>, так как у нас <math>{n \choose 2}</math> вариантов выбора i, j : <math>1 \le i \le j \le n</math> и <math>O(N)</math> точек разделения для каждого варианта. Детали станут понятны из реализации.
Реализация
В методе main – пример из начала статьи. Если запустить, выведет 7500, как и ожидается.
public class MatrixChain {
/*
* Возвращает ответ на задачу об оптимальном перемножении матриц, используя
* динамическое программирование.
* Асимптотика решения - O(N^3) время и O(N^2) память.
*
* @param p массив размеров матриц, см.статью
* @return минимальное количество скалярных умножений, необходимое для решения задачи
*/
public static int multiplyOrder(int[] p) {
int n = p.length - 1;
int[][] dp = new int[n][n];
for (int i = 0; i < n; ++i) {
dp[i][i] = 0;
}
for (int l = 1; l < n; ++l) {
for (int i = 0; i < n - l; ++i) {
int j = i + l;
dp[i][j] = Integer.MAX_VALUE;
for (int k = i; k < j; ++k) {
dp[i][j] = Math.min(dp[i][j],
dp[i][k] + dp[k + 1][j] + p[i] * p[k + 1] * p[j + 1]);
}
}
}
return dp[0][n - 1];
}
public static void main(String[] args) {
int[] test = { 10, 100, 5, 50 };
System.out.println(MatrixChain.multiplyOrder(test));
}
}
Приведены только методы, которые непосредственно выполняют необходимые расчеты
dataStore — объект класса, который хранит все данные
Его атрибуты:
public List<List<int>> m; //матрица m
public List<List<int>> s; //матрица s
public List<List<int>> result; //результат всех перемножений
public List<List<List<int>>> source; //массив из 2-мерных матриц (A0,A1,...,An) которые нужно перемножить
public List<int> sizes = new List<int>(); //размеры матриц (записаны таким образом - 12,10,7,4 => значит 3 матрицы размерами 12x10,10x7,7x4)
public string order = new string('a', 0); //правильное расположение скобок
Функциональные участки кода:
//© Paulskit 27.03.2010
//метод который находит матрицу m и s (там же под них и выделяется память)
private void matrixChainOrder(){
int n = dataStore.sizes.Count - 1;
//выделяем память под матрицы m и s
dataStore.m = new List<List<int>>();
dataStore.s = new List<List<int>>();
for (int i = 0; i < n; i++){
dataStore.m.Add(new List<int>());
dataStore.s.Add(new List<int>());
//заполняем нулевыми элементами
for (int a = 0; a < n; a++) {
dataStore.m[i].Add(0);
dataStore.s[i].Add(0);
}
}
//выполняем итерационный алгоритм
int j;
for (int l = 1; l < n; l++) {
for (int i = 0; i < n - l; i++) {
j = i + l;
dataStore.m[i][j] = int.MaxValue;
for (int k = i; k < j; k++) {
int q = dataStore.m[i][k] + dataStore.m[k + 1][j] +
dataStore.sizes[i] * dataStore.sizes[k + 1] * dataStore.sizes[j + 1];
if (q < dataStore.m[i][j]) {
dataStore.m[i][j] = q;
dataStore.s[i][j] = k;
}
}
}
}
}
//метод - простое перемножение 2-х матриц
private List<List<int>> matrixMultiply(List<List<int>> A, List<List<int>> B) {
int rowsA = A.Count;
int columnsB = B[0].Count;
//column count of A == rows count of B
int columnsA = B.Count;
//memory alloc for "c"
List<List<int>> c = new List<List<int>>();
for (int i = 0; i < rowsA; i++) {
c.Add(new List<int>());
for (int a = 0; a < columnsB; a++) {
c[i].Add(0);
}
}
//do multiplying
for (int i = 0; i < rowsA; i++) {
for (int j = 0; j < columnsB; j++) {
for (int k = 0; k < columnsA; k++) {
c[i][j] += A[i][k] * B[k][j];
}
}
}
//return value
return c;
}
//метод, который непосредственно выполняет перемножение в правильном порядке
//первоначально вызывается таким образом
//dataStore.result = matrixChainMultiply(0, dataStore.sizes.Count - 2);
private List<List<int>> matrixChainMultiply(int i, int j) {
if (j > i) {
List<List<int>> x = matrixChainMultiply(i, dataStore.s[i][j]);
List<List<int>> y = matrixChainMultiply(dataStore.s[i][j] + 1, j);
return matrixMultiply(x, y);
}
else return dataStore.source[i];
}
//метод печатающий строку с правильной расстановкой скобок
private void printOrder(int i, int j){
if (i == j) {
order += "A" + i.ToString();
} else {
order += "(";
printOrder(i, dataStore.s[i][j]);
order += "*";
printOrder(dataStore.s[i][j] + 1, j);
order += ")";
}
}
Примечания
К данной задаче сводится задача оптимизации свободной энергии молекулы РНК в биоинформатике (здесь пара скобок в строке мономеров РНК определяет их спаривание). Подобное динамическое программирование реализовано в алгоритмах Nussinov или Zucker.
Литература
- ↑ Существуют и более быстрые, чем kmn, алгоритмы умножения заполненных матриц, но они применяются крайне редко — прирост скорости наблюдается только на матрицах 100×100 и крупнее. Разреженные же матрицы умножают особыми алгоритмами, оптимизированными под ту или иную форму матрицы.
