動態規划之多矩陣乘積

算法原理請參考《算法導論》，因為算法這東西千篇一律，關鍵還是實現和理解，這里只提幾個關鍵點，幫助大家理解。

1. 為什么需要動態規划？

比如矩陣A是p x q大小，矩陣B是q x r大小，很明顯，得到的矩陣C是p x r大小，其中花費的時間必定是p*q*r。這只是兩個矩陣，如果存在N個矩陣需要算其乘積呢？那么就需要用到動態規划了，比如A(p x q), B(q x r), C(r x l)這三個矩陣相乘。如果不規划，那么花費的時間是A*B=p*q*r，然后再乘以C，還需要額外花費p*r*l時間。但有可能B*C先乘，然后再乘以A，這樣花費的時間最少。以此推廣到N個矩陣相乘。這就是動態規划的原因！！！

2. 為了幫助大家理解《算法導論》中的算法，這里我給大家講解下數組p。比如A(a x b), B(b x c), C(c x d)這三個矩陣，那么p數組的內容是p[]={a,b,c,d}，書上因為沒有采用C語言數組，所以下標從1開始，但是這里我給出下標從0開始的數組p，這樣矩陣A可以用p[0] x p[1] 表示，那么推廣到N個矩陣，矩陣A:i的維度是p[i] x p[i+1] （所以，書上給出的偽代碼需要做一些下標調整，包括算法）

下面是代碼部分

matric_free函數，內存回收

template <int p>
void matric_free(int **C) {
    for (int i = 0; i < p; i++)
        delete[] C[i];
    delete[] C;
}

matric_show函數，打印矩陣，默認值為-1時，表示未使用，則不打印。

template <int p, int q>
void matric_show(int **A) {
    for (int i = 0; i < p; i++) {
        for (int j = 0; j < q; j++)
            if (A[i][j] == -1)
                printf("%-04s", " ");
            else
                printf("%-04d ", A[i][j]);
        printf("\n");
    }
}

matric_multiply函數，這個是任意兩個矩陣相乘，從代碼也可以推出來，任意兩個矩陣相乘花費的時間是p*q*r

template <int p, int q, int r>
int **matric_multiply(int A[p][q], int B[q][r]) {
    int **C = new int *[p];
    for (int i = 0; i < r; i++)
        C[i] = new int[r];
    for (int i = 0; i<p; i++)
        for (int j = 0; j < r; j++) {
            C[i][j] = 0;
            for (int k = 0; k < q; k++)
                C[i][j] += A[i][k] * B[k][j];
        }
    return C;
}

matric_optimal函數，通過傳入數組p（構建原理，已經在上面的關鍵點中提到了）,length參數代表數組p長度，實際上，數組p會多一個元素出來，所以真實構建的m和s是length-1，由於下標通通改用了C語言形式，所以代碼我做了調整

template <int length>
void matric_optimal(int *p, int ***m, int ***s) {
    int n = length - 1, q;
    if (n < 2) return;//至少兩個矩陣
    *m = new int *[n];
    *s = new int *[n];
    for (int i = 0; i < n; i++) {//生成m, s對應的二維數組
        (*m)[i] = new int[n];
        (*s)[i] = new int[n];
    }
    //初始化
    for (int i=0;i<n;i++)
        for (int j = 0; j < n; j++) 
            (*m)[i][j] = (*s)[i][j] = -1;//-1表示未使用
    for (int i = 0; i < n; i++)
        (*m)[i][i] = 0;
    for (int l = 1; l <= n; l++) //這里針對C語言數組，做了改變
        for (int i = 0; i < n - l + 1; i++)
        {
            int j = i + l - 1;
            for (int k = i; k <j; k++) {
                q = (*m)[i][k] + (*m)[k + 1][j] + p[i] * p[k+1] * p[j+1];//這里針對C語言數組做了改變，原因很簡單，因為數組下標不能取-1，所以需要修改
                if ((*m)[i][j]==-1 || q < (*m)[i][j]) {
                    (*m)[i][j] = q;
                    (*s)[i][j] = k;//保留k值，方便輸出
                }
            }
        }
}

另一個版本是遞歸版本，是recursive_matric_optimal函數。

int recursive_matric_optimal(int *p, int ***m, int ***s, int i, int j) {
    int q, n;
    if (i == j) return 0;
    if (*m == NULL && *s == NULL) {
        n = j - i + 1;
        if (n < 2) return -1;//至少兩個矩陣
        *m = new int *[n];
        *s = new int *[n];
        for (int i = 0; i < n; i++) {//生成m, s對應的二維數組
            (*m)[i] = new int[n];
            (*s)[i] = new int[n];
        }
        //初始化
        for (int i = 0; i < n; i++)
            for (int j = 0; j < n; j++)
                (*m)[i][j] = (*s)[i][j] = -1;
    }
    for (int k = i; k < j; k++) {
        q = recursive_matric_optimal(p, m, s, i, k) + recursive_matric_optimal(p, m, s, k + 1, j) + p[i] * p[k + 1] * p[j + 1];
        if ((*m)[i][j] == -1 || q < (*m)[i][j]) {
            (*m)[i][j] = q;
            (*s)[i][j] = k;
        }
    }
    return (*m)[i][j];
}

matric_show_optimal函數，生成動態規划后的結果，直觀看得出哪個矩陣先乘，時間最短。

void matric_show_optimal(int **s, int i, int j) {
    if (i == j)
        printf("A:%d", i);
    else {
        printf("(");
        matric_show_optimal(s, i, s[i][j]);
        matric_show_optimal(s, s[i][j] + 1, j);
        printf(")");
    }
}

數據錄入

A0 30x35 
A1 35x15
A2 15x5 
A3 5x10 
A4 10x20
A5 20x25

轉換為數組p

int p[] = { 30,35,15,5,10,20,25 };//A(i)的維數表示p(i) x p(i+1)

結果圖

這個結果對應書上的結果圖（書上旋轉了），是一樣的

 

所有代碼均經過測試，結果正確！！！