Lecture 38-39

1. Lecture 38-39 Matrix Multiplication Chains 20.2.2

2. n A(i,k) B(k,j) C(i,j) = k = 1 Matrix Multiplication Chains • Multiply an m x n matrix A and an n x p matrix B to get an m x p matrix C. • We shall use the number of multiplications as our complexity measure. • n multiplications are needed to compute one C(i,j). • mnp multiplications are needed to compute all mp terms of C.

3. Matrix Multiplication Chains • Suppose that we are to compute the product XxYxZ of three matrices X, Y and Z. • The matrix dimensions are: • X:(100 x 1), Y:(1 x 100), Z:(100 x 1) • Multiply X and Y to get a 100 x 100 matrix T. • 100 * 1 * 100 = 10,000 multiplications. • Multiply T and Z to get the 100 x 1 answer. • 100 * 100 * 1 = 10,000 multiplications. • Total cost is 20,000 multiplications. • 10,000 units of space are needed for T.

4. Matrix Multiplication Chains • The matrix dimensions are: • X:(100 x 1) • Y:(1 x 100) • Z:(100 x 1) • Multiply Y and Z to get a 1 x 1 matrix T. • 1 * 100 * 1 = 100 multiplications. • Multiply X and T to get the 100 x 1 answer. • 100 * 1 * 1 = 100 multiplications. • Total cost is 200 multiplications. • 1 unit of space is needed for T.

5. Product Of 5 Matrices • Some of the ways in which the product of 5 matrices may be computed. • A*(B*(C*(D*E))) right to left • (((A*B)*C)*D)*E left to right • (A*B)*((C*D)*E) • (A*B)*(C*(D*E)) • (A*(B*C))*(D*E) • ((A*B)*C)*(D*E)

6. Evaluating all ways to compute the product takes O(4q/q0.5) time. Find Best Multiplication Order • Number of ways to compute the product of q matrices is O(4q/q1.5).

7. An Application • Registration of pre- and post-operative 3D brain MRI images to determine volume of removed tumor.

8. 3D Registration

9. 3D Registration • Each image has 256 x 256 x 256 voxels. • In each iteration of the registration algorithm, the product of three matrices is computed at each voxel … (12 x 3) * (3 x 3) * (3 x 1) • Left to right computation => 12 * 3 * 3 + 12 * 3*1 = 144 multiplications per voxel per iteration. • 100 iterations to converge.

10. 3D Registration • Total number of multiplications is about 2.4 * 1011. • Right to left computation => 3 * 3*1 + 12 * 3 * 1 = 45 multiplications per voxel per iteration. • Total number of multiplications is about 7.5 * 1010. • With 108 multiplications per second, time is 40 minvs. 12.5 min.

11. Matrix Multiplication Chains Determine the best way to compute the matrix product M1xM2 x M3x … xMq. Let the dimensions of Mi be rixri+1. q-1 matrix multiplications are to be done. Decide the matrices involved in each of these multiplications.

12. Decision Sequence • M1xM2 x M3x … xMq • Determine the q-1 matrix products in reverse order. • What is the last multiplication? • What is the next to last multiplication? • And so on.

13. Problem State • M1xM2 x M3x … xMq • The matrices involved in each multiplication are a contiguous subset of the given q matrices. • The problem state is given by a set of pairs of the form (i, j), i <= j. • The pair (i,j) denotes a problem in which the matrix product MixMi+1x … xMj is to be computed. • The initial state is (1,q). • If the last matrix product is (M1xM2x … xMk) x (Mk+1xMk+2x … xMq), the state becomes {(1,k), (k+1,q)}.

14. So the principle of optimality holds and dynamic programming may be applied. Verify Principle Of Optimality Let Mij=Mi xMi+1x … xMj, i <= j. Suppose that the last multiplication in the best way to compute Mijis MikxMk+1,j for some k, i <= k < j. Irrespective of what k is, a best computation of Mij in which the last product is MikxMk+1,j has the property that Mik and Mk+1,j are computed in the best possible way.

15. Recurrence Equations Let c(i,j) be the cost of an optimal (best)way to compute Mi j, i <= j. c(1,q) is the cost of the best way to multiply the given q matrices. Let kay(i,j) = k be such that the last product in the optimal computation of Mij is Mi k xMk+1 j. c(i,i) = 0, 1 <= i <= q. (Mi i=Mi) c(i,i+1) = riri+1ri+2, 1 <= i < q.(Mi i+1=Mi xMi+1) kay(i,i+1) = i.

16. c(i, i+s), 1 < s < q The last multiplication in the best way to compute Mi,i+s is MikxMk+1,i+s for some k, i <= k < i+s. If we knew k, we could claim: c(i,i+s) = c(i,k) + c(k+1,i+s) + rirk+1ri+s+1 Since i <= k < i+s, we can claim c(i,i+s) = min{c(i,k) + c(k+1,i+s) + rirk+1ri+s+1}, where the min is taken over i <= k < i+s. kay(i,i+s) is the k that yields above min.

17. Recurrence Equations c(i,i+s) = min i <= k < i+s{c(i,k) + c(k+1,i+s) + rirk+1ri+s+1} c(.,.) terms on right side involve fewer matrices than does the c(.,.) term on the left side. So compute in the order s = 2, 3, …, q-1.

18. Recursive Implementation • See text for recursive codes. • Code that does not avoid recomputation of already computed c(i,j)s runs in Omega(2q) time. • Code that does not recompute already computed c(i,j)s runs in O(q3) time. • Implement nonrecursively for best worst-case efficiency.

19. 1 2 3 4 1 2 3 4 kay(i,j), i <= j c(i,j), i <= j Example • q = 4, (10 x 1) * (1 x 10) * (10 x 1) * (1 x 10) • r = [r1,r2 ,r3,r4,r5] = [10, 1, 10, 1, 10]

20. 0 0 0 0 0 0 0 0 kay(i,j), i <= j c(i,j), i <= j s = 0 c(i,i) and kay(i,i) , 1 <= i <= 4 are to be computed.

21. 0 0 0 0 0 0 0 0 kay(i,j), i <= j c(i,j), i <= j s = 1 c(i,i+1) and kay(i,i+1) , 1 <= i <= 3 are to be computed.

22. 0 1 0 100 0 2 0 10 0 3 0 100 0 0 kay(i,j), i <= j c(i,j), i <= j s = 1 • c(i,i+1) = riri+1ri+2, 1 <= i < q.(Mii+1=MixMi+1) • kay(i,i+1) = i. • r = [r1,r2 ,r3,r4,r5] = [10, 1, 10, 1, 10]

23. 0 1 0 100 0 2 0 10 0 3 0 100 0 0 kay(i,j), i <= j c(i,j), i <= j s = 2 • c(i,i+2) = min{c(i,i) + c(i+1,i+2) + riri+1ri+3, • c(i,i+1) + c(i+2,i+2) + riri+2ri+3} • r = [r1,r2 ,r3,r4,r5] = [10, 1, 10, 1, 10]

24. 0 1 1 0 100 20 0 2 0 10 0 3 0 100 0 0 s = 2 • c(1,3) = min{c(1,1) + c(2,3) + r1r2r4, • c(1,2) + c(3,3) + r1r3r4} • r = [r1,r2 ,r3,r4,r5] = [10, 1, 10, 1, 10] • c(1,3) = min{0 + 10 + 10, 100 + 0 + 100} kay(i,j), i <= j c(i,j), i <= j

25. 0 1 1 0 100 20 0 2 3 0 10 20 0 3 0 100 0 0 s = 2 • c(2,4) = min{c(2,2) + c(3,4) + r2r3r5, • c(2,3) + c(4,4) + r2r4r5} • r = [r1,r2 ,r3,r4,r5] = [10, 1, 10, 1, 10] • c(2,4) = min{0 + 100 + 100, 10 + 0 + 10} kay(i,j), i <= j c(i,j), i <= j

26. 0 1 1 1 0 100 20 120 0 2 3 0 10 20 0 3 0 100 0 0 s = 3 • c(1,4) = min{c(1,1) + c(2,4) + r1r2r5, • c(1,2) + c(3,4) + r1r3r5, c(1,3) + c(4,4) + r1r4r5} • r = [r1,r2 ,r3,r4,r5] = [10, 1, 10, 1, 10] • c(1,4) = min{0+20+100, 100+100+1000, 20+0+100} kay(i,j), i <= j c(i,j), i <= j

27. 0 1 1 1 0 100 20 120 0 2 3 0 10 20 0 3 0 100 0 0 Determine The Best Way To ComputeM14 • kay(1,4) = 1. • So the last multiplication is M14=M11 xM24. • M11 involves a single matrix and no multiply. • Find best way to compute M24. kay(i,j), i <= j c(i,j), i <= j

28. 0 1 1 1 0 100 20 120 0 2 3 0 10 20 0 3 0 100 0 0 Determine The Best Way To ComputeM24 • kay(2,4) = 3 . • So the last multiplication is M24=M23 xM44. • M23=M22 xM33. • M44 involves a single matrix and no multiply. kay(i,j), i <= j c(i,j), i <= j

29. The Best Way To ComputeM14 • The multiplications (in reverse order) are: • M14=M11 xM24 • M24=M23 xM44 • M23=M22 xM33

30. 1 2 3 4 1 c(i,j), i <= j 2 3 4 0 100 20 120 0 10 20 0 100 0 Time Complexity • O(q2) values of c(i,j) are to be computed, where q is the number of matrices. • c(i,i+s) = min i <= k < i+s{c(i,k) + c(k+1,i+s) + rirk+1ri+s+1}. • Each c(i,j) is the min of O(q) terms. • Each of these terms is computed in O(1) time. • So all c(i,j) are computed in O(q3) time.

31. 1 2 3 4 1 kay(i,j), i <= j 2 3 4 0 1 1 1 0 2 3 0 3 0 Time Complexity • The traceback takes O(1) time to determine each matrix product that is to be done. • q-1 products are to be done. • Traceback time is O(q).

32. Time Complexity • Total time is O(q3) + O(q) = O(q3).

33. Next: Lecture 40 • All Pairs Shortest Paths • 20.2.3