lecture 11 matrix algorithms and applications
play

Lecture 11. Matrix Algorithms and Applications Introduction to - PowerPoint PPT Presentation

Lecture 11. Matrix Algorithms and Applications Introduction to Computer Science, IIIS, Tsinghua Outline } Introduction to divide and conquer } faster matrix multiplication algorithm } Algorithms for matrix operations } LUP decomposition, matrix


  1. LUP decomposition a w T 1 a w T 1 = × v 1 v A’ / A’-vw T /a … a 1 If we change the order of the 2nd line and the i-th line B Still a lower triangle matrix!

  2. LUP decomposition-summary } LU-decomposition A=LU } When the element a ii =0 } Change the order of the rows by multiplying a permutation matrix Q } Finally A=E 1 …E k U, for a sequence of elementary matrix = } So finally we can write A=PLU } L is a unit lower-triangular matrix } U is an upper-triangular matrix } P is a permutation matrix

  3. LUP-decomposition } After the LUP-decomposition A=PLU, every operation becomes easy: } find inversion } compute determinant det(A)=det(P)det(L)det(U) } Solve linear system (Ax=b) (PLU)x=b ó {Py=b; Lz=y; Ux=z}, takes O(n 2 ) time.

  4. LUP-decomposition } After the LUP-decomposition, every operation becomes easy: } find inversion } compute determinant } Solve linear system (Ax=b) } How can we find the inverse of a triangle matrix? } better than O(n 3 )

  5. LUP-decomposition } After the LUP-decomposition, every operation becomes easy: } find inversion } compute determinant } Solve linear system (Ax=b) } How can we find the inverse of a triangle matrix? − 1 = A − 1 ⎡ − A − 1 CB − 1 ⎤ ⎡ ⎤ A C ⎢ ⎥ ⎢ ⎥ B − 1 0 B 0 ⎣ ⎦ ⎣ ⎦ } So we can recursively compute A -1 and B -1 , (also triangle) } Time: T(n)=2T(n/2)+ Θ (n ω )= Θ (n ω ) (the same as matrix multiplication)

  6. } What’s the running time for Gaussian elimination?

  7. } What’s the running time for Gaussian elimination? } O(n 3 ) a w T 1 a w T 1 = × v 1 v A’ / A’-vw T /a … a 1 } How can we do the LUP-decomposition faster?

  8. How to accelerate LUP decomposition? } Gaussian elimination: } Eliminate elements by order: } O(n 3 ) time.

  9. How to accelerate LUP decomposition? } The Hopcroft-Bunch (1974) version: } Eliminate elements by order:

  10. Hopcroft-Bunch algorithm U X V Y

  11. Hopcroft-Bunch algorithm U X Y -VU -1 X

  12. Hopcroft-Bunch algorithm • U becomes upper triangle! • Find U -1 in O(n ω ) time U X V Y

  13. Hopcroft-Bunch algorithm • Block elimination! U X Y -VU -1 X

  14. Hopcroft-Bunch algorithm U X V Y

  15. Hopcroft-Bunch algorithm U X Y -VU -1 X

  16. Hopcroft-Bunch algorithm • U is always upper-triangle, so it’s easy to find its inverse • We just need to pick a non-zero element at k-th row U X V Y

  17. Hopcroft-Bunch algorithm U X Recursively run this on Y -VU -1 X Y -VU -1 X

  18. Hopcroft-Bunch algorithm Every time we just need to multiply a lower-triangle matrix: 1 U X 1 1 VU -1 V Y 1 1 1 1 1

  19. Running Time: } M(n,p): time for multiplying an n × n and an n × p matrix } I(n): time for inverse an n × n matrix } Time: } Inversions: log n n n n i 1 ) I ( 1 ) I ( 2 ) ... I ( n / 2 ) i I ( 2 − ∑ + + + = 2 4 2 i 1 =

  20. Running Time: } M(n,p): time for multiplying an n × n and an n × p matrix } I(n): time for inverse an n × n matrix } Time: } Inversions: log n n n n i 1 ) I ( 1 ) I ( 2 ) ... I ( n / 2 ) i I ( 2 − ∑ + + + = 2 4 2 i 1 = } Multiplications: every time we need to compute Y-VU -1 X log n n i 1 M ( 2 − , n ) ∑ ≤ i 2 i 1 =

  21. Running Time: } M(n,p): time for multiplying an n × n and an n × p matrix } I(n): time for inverse an n × n matrix What’s M(n,p) ? } Time: } Inversions: log n n n n i 1 ) I ( 1 ) I ( 2 ) ... I ( n / 2 ) i I ( 2 − ∑ + + + = 2 4 2 i 1 = } Multiplications: every time we need to compute Y-VU -1 X log n n i 1 M ( 2 − , n ) ∑ ≤ i 2 i 1 =

  22. Running Time: } M(n,p): O(pn 𝜕 -1 ) } a trivial method, there is good ways } I(n): O(n 𝜕 ) } Time: } Inversions: log n n n n i 1 ) I ( 1 ) I ( 2 ) ... I ( n / 2 ) i I ( 2 − ∑ + + + = 2 4 2 i 1 = } Multiplications: every time we need to compute Y-VU -1 X log n n i 1 M ( 2 − , n ) ∑ ≤ i 2 i 1 = Thus, the total time is still O(n 𝜕 ).

  23. Summary } A reduces to B: we can use an algorithm for B to solve A (A è B) } We have shown: Matrix Multiplication LUP-decomposition Matrix inversion } so LUP-decomposition and matrix inversion has O(n ω ) time algorithm.

  24. Summary } A reduces to B: we can use an algorithm for B to solve A (A è B) } We have shown: Matrix Multiplication LUP-decomposition Matrix inversion } if we can show this reduction, then the three operations has the same time complexity! } there is an oracle for matrix inversion, compute A × B?

  25. Multiplication reduces to Inversion } Very tricky result: } To compute AB, − 1 ⎛ ⎞ ⎛ ⎞ I A I − A AB ⎜ ⎟ ⎜ ⎟ I B I − B = ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ I I ⎝ ⎠ ⎝ ⎠

  26. Summary } A reduces to B: we can use an algorithm for B to solve A (A è B) Matrix Multiplication LUP-decomposition Matrix inversion Compute det(A) Solve Ax=b } All operation can be solved in O(n ω ) time.

  27. Eigenvalues/Eigenvector? } However, there is no such algorithm by FMM for computing eigenvalues.

  28. Example of applications } Matching M: } A set of vertex-disjoint edges } Matched vertices: the vertices associated with an edge in M } Free vertices: unmatched vertices } Perfect matching in bipartite graph: } No free vertices

  29. Another example } Perfect matching in bipartite graph: } Problem: decide whether there is a perfect matching in G

  30. Adjacency matrix } The adjacency matrix A of graph G=(V, E): } A i,j =1 iff (i,j) ∈ E 1 1 1 0 1 0 1 0 0 1 0 1 0 0 1 1

  31. } Perfect matching in bipartite graph: } Whether there is a perfect matching in G } By determinant of adjacency matrix 1 1 1 0 1 0 1 0 0 1 0 1 0 0 1 1

  32. Another example } Perfect matching in bipartite graph: } Whether there is a perfect matching in G } By determinant of adjacency matrix a 11 a 12 a 13 a 14 a 21 a 22 a 23 a 24 a 31 a 32 a 33 a 34 } σ is the permutation of {1,…,n} a 41 a 42 a 43 a 44 } sgn( σ ) is +1 or -1

  33. Another example } Perfect matching in bipartite graph: } Whether there is a perfect matching in G } By determinant of adjacency matrix a 11 a 12 a 13 a 14 a 21 a 22 a 23 a 24 a 31 a 32 a 33 a 34 } σ is the permutation of {1,…,n} a 41 a 42 a 43 a 44 } sgn( σ ) is +1 or -1 Permutation means every row and column contains exactly one element

  34. } A perfect matching in a bipartite graph with n vertices each side can also be seen as a permutation of {1,…,n}

  35. } A perfect matching in a bipartite graph with n vertices each side can also be seen as a permutation of {1,…,n} } Construct an n × n matrix X: x 11 x 12 x 13 0 x 21 0 x 23 0 0 x 32 0 x 34 0 0 x 43 x 44 } x ij is a variable is there is an edge from A i to B j } x ij =0 if there is no edge from A i to B j

  36. } A perfect matching in a bipartite graph with n vertices each side can also be seen as a permutation of {1,…,n} } Construct an n × n matrix X: x 11 x 12 x 13 0 x 21 0 x 23 0 0 x 32 0 x 34 0 0 x 43 x 44 } From the definition of determinant, det(X) ≡ 0 if there is no perfect matching.

  37. } The determinant is a polynomial of degree n } We can randomly pick each nonzero x ij from {1,…,2n}, and compute det(X) } if there is no perfect matching, det(X) is always 0 } otherwise, det(X) is not 0 unless we are unlucky enough to pick the root of the polynomial det(X). } Computing det(X) takes O(n ω ) time. 5 3 8 0 1 0 4 0 0 7 0 3 0 0 2 6

  38. } The determinant is a polynomial of degree n } We can randomly pick each nonzero x ij from {1,…,2n}, and compute det(X) } As in the polynomial identity testing } Computing det(X) takes O(n ω ) time. } If there is no perfect matching, determinant is always 0. } If there are perfect matchings, Pr[det(X) ≠ 0] ≥ 1/2 } Proof in the “Theory of Computation” course. } So we can repeat k time to make error probability ≤ 1/2 k .

  39. Graph algorithms by FMM } Transitive Closure: O(n ω ) } All-pair shortest path (APSP): } undirected unweighted: O(n ω ) } directed unweighted: O(n 𝜈 ) } (1+ ε )-approximate with weights [0…M]: O(n 𝜕 log M/ ε ) } Maximum Matching: O(n ω ) } Dominance Product: O(n (3+ ω )/2 ) } All-pair bottleneck path (APBP): } vertex-capacitated graphs: O(n 𝜈 ) } edge-capacitated graphs: O(n (3+ ω )/2 )

  40. Graph algorithms by FMM } Transitive Closure: O(n ω ) } All-pair shortest path (APSP): } undirected unweighted: O(n ω ) } directed unweighted: O(n 𝜈 ) } (1+ ε )-approximate with weights [0…M]: O(n 𝜕 log M/ ε ) } Maximum Matching: O(n ω ) } Dominance Product: O(n (3+ ω )/2 ) } All-pair bottleneck path (APBP): } vertex-capacitated graphs: O(n 𝜈 ) } edge-capacitated graphs: O(n (3+ ω )/2 ) 0 𝜈 ≈ 2.531 (3+ ω )/2 ≈ 2.687 3 ω≈ 2.373

  41. Simple problem } In a undirected graph G, check whether G contains a triangle. } trivial algorithm: checking every triple of vertices (u,v,w) } time: O(n 3 )

  42. Simple problem } In a undirected graph G, check whether G contains a triangle. } By FMM?

  43. Adjacency matrix of G=(V,E) } A i,j =1 iff (i,j) ∈ E } Consider A 2 : } (A 2 ) i,j = ∑ k A i,k ·A k,j × = A ( a ) B ( b ) C ( c ) = = = i j i j i j

  44. Adjacency matrix of G=(V,E) } A i,j =1 iff (i,j) ∈ E } Consider A 2 : } (A 2 ) i,j = ∑ k A i,k ·A k,j } (A 2 ) i,j >0 equivalent to ∃ k such that A i,k =A k,j =1, } that is (i,k), (k,j) ∈ E

  45. Adjacency matrix of G=(V,E) } A i,j =1 iff (i,j) ∈ E } Consider A 2 : } (A 2 ) i,j = ∑ k A i,k ·A k,j } (A 2 ) i,j >0 equivalent to ∃ k such that A i,k =A k,j =1, } that is (i,k), (k,j) ∈ E } So G has a triangle: ∃ (i,j) such that A i,j =1 and (A 2 ) i,j ≥ 1 3 1 2 7 6 5 4

  46. What about Quadrilateral? } Checking whether a graph contains a 4-cycle: 3 1 2 7 6 5 4

  47. What about Quadrilateral? } Checking whether a graph contains a 4-cycle: } There must be (i,j) such that (A 2 ) i,j ≥ 2: } (A 2 ) 2,6 =2 sinceA 2,7 =A 7,6 =1, A 2,5 =A 5,6 =1 } Thus time bound is also O(n 𝜕 ) 3 1 2 7 6 5 4

  48. Problem solved? } What about the 4-node subgraph is exactly a 4 cycle? } subgraph on the 4 vertices has 4 edges 3 3 1 2 1 2 7 7 ✅ ❌ 6 6 5 4 5 4

  49. Problem solved? } What about the 4-node subgraph is exactly a 4 cycle? } subgraph on the 4 vertices has 4 edges 3 3 1 2 1 2 7 7 ✅ ❌ 6 6 5 4 5 4

  50. Randomized algorithm } [Vassilevska, Wang, Williams, Yu 2014] all 4-node subgraphs except “clique” and “independent set” can be found in O(n 𝜕 ) time with high probability.

  51. Ideas } Consider v u ✓ ( A 2 ) u , v ◆ ∑ 2 ( u , v ) 2 E ✓ ◆ } What if we have a “diamond” } What if we have a “clique”

  52. Ideas } Consider v u ✓ ( A 2 ) u , v ◆ ∑ 2 ( u , v ) 2 E ✓ ◆ } if we have a “diamond” } can only be: } if we have a “clique” } any edge can be (u,v) } so it is counted 6 times

  53. Ideas } So ✓ ( A 2 ) u , v ◆ ∑ = 6 × #( )+#( ) 2 ( u , v ) 2 E ✓ ◆ } But it ’ s still hard to separate them

  54. Ideas } So ✓ ( A 2 ) u , v ◆ X= = 6 × #( )+#( ) ∑ 2 ( u , v ) 2 E } If X mod 6 ≠ 0, then of course the graph contains a ✓ ◆ diamond } If X mod 6 = 0, we randomly delete some vertices } with some probability, the number of diamonds will not be a multiple of 6, then return true } we can run this procedure many times, if the graph contains no diamond, it always return false

  55. What about 4-cycle?

  56. What about 4-cycle? } Consider 2 ✓ ( A 2 ) u , v ◆ ∑ 2 ( u , v ) 62 E ∑ } If we have a diamond… } If we have a 4-cycle…

  57. What about 4-cycle? } Consider 2 ✓ ( A 2 ) u , v ◆ ∑ 2 ( u , v ) 62 E ∑ } If we have a diamond, it is counted once!! } If we have a 4-cycle, it is counted twice

  58. Ideas } So ✓ ( A 2 ) u , v ◆ ∑ X= = 6 × #( )+#( ) 2 ( u , v ) 2 E ✓ ◆ 2 ✓ ( A 2 ) u , v ◆ ∑ Y= = 2 × #( )+#( ) 2 ( u , v ) 62 E ∑ } So by X-Y=6 × (#clique)-2 × (#4-cycle), using similar idea we can check whether a graph contains a 4-cycle.

  59. Dominance Product ( ♦ ) } Original product } C[i,j]= ∑ k A[i,k]·B[k,j] } Dominance Product } M[i,j]=(A s B)[i,j]=|{k | A[i,k] ≤ B[k,j]}|

  60. Dominance Product ( ♦ ) } Dominance Product } M[i,j]=(A s B)[i,j]=|{k | A[i,k] ≤ B[k,j]}| M[i,j]=3 (Dominance Product)

  61. Dominance Product ( ♦ ) } Dominance Product } M[i,j]=(A s B)[i,j]=|{k | A[i,k] ≤ B[k,j]}| } From Matou š ek 1991, Dominance Product for any two n × n matrices can be computed in O(n (3+ ω )/2 ) time } From V. Vassilevska et al 2007, for two sparse matrices, with m 1 and m 2 valid elements respectively, the running time will be ( 1 ) / 2 O ( m m n ) ω − 1 2

  62. Algorithm for Dominance Product? } Dominance Product } M[i,j]=(A s B)[i,j]=|{k | A[i,k] ≤ B[k,j]}|

  63. Algorithm for Dominance Product? } Dominance Product } M[i,j]=(A s B)[i,j]=|{k | A[i,k] ≤ B[k,j]}| } A[i,1] compare with B[1,j].

  64. Algorithm for Dominance Product? } Dominance Product } M[i,j]=(A s B)[i,j]=|{k | A[i,k] ≤ B[k,j]}| } A[i,1] compare with B[1,j]. } Thus, elements of the i-th column of A are only compared with i-th row of B.

  65. Algorithm for Dominance Product? } Dominance Product } M[i,j]=(A s B)[i,j]=|{k | A[i,k] ≤ B[k,j]}| } A[i,1] compare with B[1,j]. } Thus, elements of the i-th column of A are only compared with i-th row of B.

  66. Algorithm for Dominance Product? } Dominance Product } M[i,j]=(A s B)[i,j]=|{k | A[i,k] ≤ B[k,j]}| } A[i,1] compare with B[1,j]. } Thus, elements of the i-th column of A are only compared with i-th row of B. } Thus, define the set L k ={A[i,k]} ∪ {B[k,j]}, and sort L k } Partition the sorted list of each L k into r parts: } each with 2n/r elements Increasing L k

  67. } Construct: } A b [i,k]= 1 if A[i,k] is in b-th part of L k 0 otherwise if i<j, A i ·B j is always } B b [k,j]= 1 if B[k,j] is in b-th part of L k included in the dominance product 0 otherwise

  68. M[i,j]=(A s B)[i,j]=|{k | A[i,k] ≤ B[k,j]}| } For each b=1…r, compute the Boolean product C b =A b ·(B b+1 +B b+2 +...+B r ), then compute ∑ b C b . } Since the corresponding elements in b-th part of A are ≤ the (b+1)-th part of B } What’s left?

  69. M[i,j]=(A s B)[i,j]=|{k | A[i,k] ≤ B[k,j]}| } For each b=1…r, compute the Boolean product C b =A b ·(B b+1 +B b+2 +...+B r ), then compute ∑ b C b . } Then we only need to compare elements in the same parts of the partition } Running Time?

  70. M[i,j]=(A s B)[i,j]=|{k | A[i,k] ≤ B[k,j]}| } For each b=1…r, compute the Boolean product C b =A b ·(B b+1 +B b+2 +...+B r ), then compute ∑ b C b . } Time: O(r·n 𝜕 ) } Then we only need to compare elements in the same parts of the partition } Time: for each part, O(n·(n/r) 2 ), in total: O(rn·(n/r) 2 )=O(n 3 /r) } To balance, we let r=n (3- 𝜕 )/2 , total time: O(n (3+ 𝜕 )/2 )

Recommend


More recommend