Editing Dynamic programming (section)

=== Matrix chain multiplication ===
{{unreferenced section|date=May 2013}}
{{Main|Matrix chain multiplication}}
<!--Show how the placement of parentheses affects the number of scalar multiplications required when multiplying a bunch of matrices.
Show how to write a dynamic program to calculate the optimal parentheses placement.
This is such a long example that it might be better to make it its own article.-->

Matrix chain multiplication is a well-known example that demonstrates utility of dynamic programming. For example,  engineering applications often have to multiply a chain of matrices. It is not surprising to find matrices of large dimensions, for example 100×100. Therefore, our task is to multiply matrices {{tmath|A_1, A_2, .... A_n}}. Matrix multiplication is not commutative, but is associative; and we can multiply only two matrices at a time. So, we can multiply this chain of matrices in many different ways, for example:

: {{math|((A<sub>1</sub> × A<sub>2</sub>) × A<sub>3</sub>) × ... A<sub>n</sub>}}

: {{math|A<sub>1</sub>×(((A<sub>2</sub>×A<sub>3</sub>)× ... ) × A<sub>n</sub>)}}

: {{math|(A<sub>1</sub> × A<sub>2</sub>) × (A<sub>3</sub> × ... A<sub>n</sub>)}}

and so on. There are numerous ways to multiply this chain of matrices. They will all produce the same final result, however they will take more or less time to compute, based on which particular matrices are multiplied. If matrix A has dimensions m×n and matrix B has dimensions n×q, then matrix C=A×B will have dimensions m×q, and will require m*n*q scalar multiplications (using a simplistic [[matrix multiplication algorithm]] for purposes of illustration).

For example, let us multiply matrices A, B and C. Let us assume that their dimensions are m×n, n×p, and p×s, respectively. Matrix A×B×C will be of size m×s and can be calculated in two ways shown below:

# Ax(B×C)    This order of matrix multiplication will require nps + mns scalar multiplications.  
# (A×B)×C    This order of matrix multiplication will require mnp + mps scalar calculations.

Let us assume that m = 10, n = 100, p = 10 and s = 1000. So, the first way to multiply the chain will require 1,000,000 + 1,000,000 calculations. The second way will require only 10,000+100,000 calculations. Obviously, the second way is faster, and we should multiply the matrices using that arrangement of parenthesis.

Therefore, our conclusion is that the order of parenthesis matters, and that our task is to find the optimal order of parenthesis.

At this point, we have several choices, one of which is to design a dynamic programming algorithm that will split the problem into overlapping problems and calculate the optimal arrangement of parenthesis. The dynamic programming solution is presented below.

Let's call m[i,j] the minimum number of scalar multiplications needed to multiply a chain of matrices from matrix i to matrix j (i.e. A<sub>i</sub> × .... × A<sub>j</sub>, i.e. i<=j). We split the chain at some matrix k, such that i <= k < j, and try to find out which combination produces minimum m[i,j].

The formula is:

        '''if''' i = j, m[i,j]= 0
        '''if''' i < j, m[i,j]= min over all possible values of k {{nowrap|(m[i,k]+m[k+1,j] + <math>p_{i-1}*p_k*p_j</math>)}} 
where ''k'' ranges from ''i'' to ''j''&nbsp;&minus;&nbsp;1.

*{{tmath|p_{{(}}i-1{{)}}}} is the row dimension of matrix i,
*{{tmath|p_k}} is the column dimension of matrix k,
*{{tmath|p_j}} is the column dimension of matrix j.

This formula can be coded as shown below, where input parameter "chain" is the chain of matrices, i.e. {{tmath|A_1, A_2, ... A_n}}:

 '''function''' OptimalMatrixChainParenthesis(chain)
     n = length(chain)
     '''for''' i = 1, n
         m[i,i] = 0    ''// Since it takes no calculations to multiply one matrix''
     '''for''' len = 2, n
         '''for''' i = 1, n - len + 1
             j = i + len -1
             m[i,j] = infinity      ''// So that the first calculation updates''
             '''for''' k = i, j-1
                 {{nowrap|1=q = m[i, k] + m[k+1, j] + <math>p_{i-1}*p_k*p_j</math>}}
                 '''if''' q < m[i, j]     ''// The new order of parentheses is better than what we had''
                     m[i, j] = q    ''// Update''
                     s[i, j] = k    ''// Record which k to split on, i.e. where to place the parenthesis''

So far, we have calculated values for all possible {{math|''m''[''i'', ''j'']}}, the minimum number of calculations to multiply a chain from matrix ''i'' to matrix ''j'', and we have recorded the corresponding "split point"{{math|''s''[''i'', ''j'']}}. For example, if we are multiplying chain {{math|A<sub>1</sub>×A<sub>2</sub>×A<sub>3</sub>×A<sub>4</sub>}}, and it turns out that {{math|1=''m''[1, 3] = 100}} and {{math|1=''s''[1, 3] = 2}}, that means that the optimal placement of parenthesis for matrices 1 to 3 is {{tmath|(A_1\times A_2)\times A_3}} and to multiply those matrices will require 100 scalar calculations.

This algorithm will produce "tables" ''m''[, ] and ''s''[, ] that will have entries for all possible values of i and j. The final solution for the entire chain is m[1, n], with corresponding split at s[1, n]. Unraveling the solution will be recursive, starting from the top and continuing until we reach the base case, i.e. multiplication of single matrices.

Therefore, the next step is to actually split the chain, i.e. to place the parenthesis where they (optimally) belong. For this purpose we could use the following algorithm:

 '''function''' PrintOptimalParenthesis(s, i, j)
     '''if''' i = j
         print "A"i
     '''else'''
         print "(" 
         PrintOptimalParenthesis(s, i, s[i, j]) 
         PrintOptimalParenthesis(s, s[i, j] + 1, j) 
         print ")"

Of course, this algorithm is not useful for actual multiplication. This algorithm is just a user-friendly way to see what the result looks like.

To actually multiply the matrices using the proper splits, we need the following algorithm:
<syntaxhighlight lang="javascript">
   function MatrixChainMultiply(chain from 1 to n)       // returns the final matrix, i.e. A1×A2×... ×An
      OptimalMatrixChainParenthesis(chain from 1 to n)   // this will produce s[ . ] and m[ . ] "tables"
      OptimalMatrixMultiplication(s, chain from 1 to n)  // actually multiply

   function OptimalMatrixMultiplication(s, i, j)   // returns the result of multiplying a chain of matrices from Ai to Aj in optimal way
      if i < j
         // keep on splitting the chain and multiplying the matrices in left and right sides
         LeftSide = OptimalMatrixMultiplication(s, i, s[i, j])
         RightSide = OptimalMatrixMultiplication(s, s[i, j] + 1, j)
         return MatrixMultiply(LeftSide, RightSide) 
      else if i = j
         return Ai   // matrix at position i
      else 
         print "error, i <= j must hold"

    function MatrixMultiply(A, B)    // function that multiplies two matrices
      if columns(A) = rows(B) 
         for i = 1, rows(A)
            for j = 1, columns(B)
               C[i, j] = 0
               for k = 1, columns(A)
                   C[i, j] = C[i, j] + A[i, k]*B[k, j] 
               return C 
      else 
          print "error, incompatible dimensions."
</syntaxhighlight>