Solving Linear System of Equations The Undo button for Linear - PowerPoint PPT Presentation

Solving Linear System of Equations

The “Undo” button for Linear Operations Matrix-vector multiplication: given the data 𝒚 and the operator 𝑩 , we can find 𝒛 such that 𝒛 = 𝑩 𝒚 𝑩 𝒚 𝒛 transformation What if we know 𝒛 but not 𝒚 ? How can we “undo” the transformation? 𝑩 !𝟐 𝒚 𝒛 Solve 𝑩 𝒚 = 𝒛 for 𝒚 ?

Image Blurring Example • Image is stored as a 2D array of real numbers between 0 and 1 (0 represents a white pixel, 1 represents a black pixel) • 𝒚𝒏𝒃𝒖 has 40 rows of pixels and 100 columns of pixels • Flatten the 2D array as a 1D array • 𝒚 contains the 1D data with dimension 4000, • Apply blurring operation to data 𝒚, i.e. 𝒄 = 𝑩 𝒚 where 𝑩 is the blur operator and 𝒄 is the blurred image

Blur operator 𝒚 Blur operator 𝑩 blurred "original” Blur operator 𝒄 image image (4000,4000) (4000,) (4000,) 𝒄 = 𝑩 𝒚

”Undo” Blur to recover original image 𝒄 Assumptions: 1. we know the blur operator 𝑩 Solve 2. the data set 𝒄 does not 𝑩 𝒚 = 𝒄 have any noise (“clean for 𝒚 data” 𝒚 What happens if we add some noise to 𝒄 ?

”Undo” Blur to recover original image 𝒛 + 𝑏 ∗ 10 !" (𝑏 ∈ 0,1 ) 𝒛 + 𝑏 ∗ 10 !# (𝑏 ∈ 0,1 ) Solve 𝑩 𝒚 = 𝒛 for 𝒚 How much noise can we add and still be able to recover meaningful information from the original image? At which point this inverse transformation fails? We will talk about sensitivity of the “undo” operation later.

Linear System of Equations How do we actually solve 𝑩 𝒚 = 𝒄 ? We can start with an “easier” system of equations… Let’s consider triangular matrices (lower and upper): 𝑦 ! 𝑀 !! 0 𝑐 … 0 𝑦 " 𝑀 "! 𝑀 "" 𝑐 " … 0 = ⋮ ⋮ ⋮ ⋱ ⋮ ⋮ 𝑦 # 𝑀 #! 𝑀 #" … 𝑀 ## 𝑐 # 𝑦 ! 𝑉 !! 𝑉 !" … 𝑉 !# 𝑐 𝑦 " 0 𝑉 "" … 𝑉 "# 𝑐 " = ⋮ ⋱ ⋮ ⋮ ⋮ ⋮ 𝑦 # … 𝑉 ## 𝑐 # 0 0

Example: Forward-substitution for lower triangular systems 𝑦 ! 2 0 0 0 2 𝑦 " 3 2 0 0 2 = 𝑦 # 1 2 6 0 6 𝑦 $ 1 3 4 2 4 2 𝑦 $ = 2 → 𝑦 $ = 1 3 𝑦 $ + 2 𝑦 % = 2 → 𝑦 % = 2 − 3 = −0.5 2 1 𝑦 $ + 2 𝑦 % + 6 𝑦 & = 6 → 𝑦 & = 6 − 1 + 1 = 1.0 6 1 𝑦 $ + 3 𝑦 % + 4 𝑦 & + 2 𝑦 # = 4 → 𝑦 & = 4 − 1 + 1.5 − 4 = 0.25 𝑦 ! 1 2 𝑦 " −0.5 = 𝑦 # 1.0 𝑦 $ 0.25

Example: Backward-substitution for upper triangular systems 𝑦 ! 2 8 4 2 2 𝑦 " 0 4 4 3 4 = 𝑦 $ 0 0 6 2 4 𝑦 % 1 0 0 0 2 𝑦 # = 1 2 4 − 2 1 = 1 2 𝑦 & = 6 2 4 − 4 1 2 − 3 1 = 1/2 = 1 2 𝑦 % = 4 4 8 2 − 8 1 8 − 4 1 2 − 2 1 = −2 2 𝑦 $ = 2 = −1 2

LU Factorization How do we solve 𝑩 𝒚 = 𝒄 when 𝑩 is a non-triangular matrix? We can perform LU factorization: given a 𝑜×𝑜 matrix 𝑩 , obtain lower triangular matrix 𝑴 and upper triangular matrix 𝑽 such that 𝑩 = 𝑴𝑽 where we set the diagonal entries of 𝑴 to be equal to 1. 1 0 𝑉 !! 𝑉 !" … 𝑉 !# 𝐵 !! 𝐵 !" … 𝐵 !# … 0 𝑀 "! 1 0 𝑉 "" … 𝑉 "# 𝐵 "! 𝐵 "" … 𝐵 "# … 0 = ⋮ ⋮ ⋱ ⋮ ⋮ ⋮ ⋱ ⋮ ⋱ ⋮ ⋮ ⋮ 𝑀 #! 𝑀 #" … 𝑉 ## 𝐵 #! 𝐵 #" … 𝐵 ## … 1 0 0

LU Factorization … 𝑉 !# … 𝐵 !# 1 0 𝑉 !! 𝑉 !" 𝐵 !! 𝐵 !" … 0 𝑀 "! 1 0 𝑉 "" … 𝑉 "# 𝐵 "! 𝐵 "" … 𝐵 "# … 0 = ⋮ ⋮ ⋱ ⋮ ⋮ ⋮ ⋱ ⋮ ⋱ ⋮ ⋮ ⋮ 𝑀 #! 𝑀 #" … 𝑉 ## 𝐵 #! 𝐵 #" … 𝐵 ## … 1 0 0 Assuming the LU factorization is know, we can solve the general system

LU Factorization (with pivoting) 𝑩 = 𝑸𝑴𝑽 Factorize: 𝑸𝑴𝑽 𝒚 = 𝒄 𝒛 𝑴 𝒛 = 𝑸 𝑼 𝒄 (Solve for 𝒛 ) Forward-substitution 𝑽 𝒚 = 𝒛 (Solve for 𝒚 ) Backward-substitution

Example Assume the 𝑩 = 𝑴𝑽 factorization is known, yielding: 1 0 0 0 2 8 4 1 0.5 1 0 0 0 −2 1 2.5 𝑴 = 𝑽 = 0.5 1 1 0 0 0 3 −1 2 0.5 1 0.5 0.5 0 0 0 0.75 2 Determine the solution 𝒚 that satisfies 𝑩𝒚 = 𝒄 , when 𝒄 = 1 4 𝑴 𝑽𝒚 = 𝒄 𝒛 First, solve the lower-triangular Then, solve the upper-triangular system 𝑴 𝒛 = 𝒄 for the variable 𝒛 system 𝑽 𝒚 = 𝒛 for the variable 𝒚 1 0 0 0 2 2 8 4 1 2 0.5 1 0 0 0 −2 1 2.5 2 1 𝒚 = 𝒛 = 0.5 1 1 0 1 0 0 3 −1 −1 0.5 0.5 0.5 1 4 0 0 0 0.75 3

Methods to solve linear system of equations 𝑩 𝒚 = 𝒄 • LU • Cholesky • Sparse

LU Factorization - Algorithm

2x2 LU Factorization (simple example) 𝐵 !! 𝐵 !" 𝑉 !! 𝑉 !" 1 0 = 𝑀 "! 1 𝐵 "! 𝐵 "" 0 𝑉 "" 𝐵 !! 𝐵 !" 𝑉 !! 𝑉 !" = 𝐵 "! 𝐵 "" 𝑀 "! 𝑉 !! 𝑀 "! 𝑉 !" + 𝑉 ""

𝑏 $$ : scalar LU Factorization 𝒃 $% : row vector (1×(𝑜 − 1)) 𝒃 %$ : column vector (𝑜 − 1)×1 𝑏 $$ 𝒃 $% 𝑩 %% : matrix (𝑜 − 1)×(𝑜 − 1) 𝐵 !! 𝐵 !" … 𝐵 !# 𝑏 !! 𝒃 !" 𝑣 !! 𝒗 !" 𝐵 "! 𝐵 "" … 𝐵 "# 1 𝟏 = = 𝒃 "! 𝑩 "" 𝟏 𝑽 "" 𝒎 "! 𝑴 "" ⋮ ⋮ ⋱ ⋮ 𝐵 #! 𝐵 #" … 𝐵 ## 𝒃 %$ 𝑩 %% 𝑏 !! 𝒃 !" 𝑣 !! 𝒗 !" = 𝒃 "! 𝑩 "" 𝑣 !! 𝒎 "! 𝒎 "! 𝒗 !" + 𝑴 "" 𝑽 ""

𝑏 $$ : scalar LU Factorization 𝒃 $% : row vector (1×(𝑜 − 1)) 𝒃 %$ : column vector (𝑜 − 1)×1 𝑏 $$ 𝒃 $% 𝑩 %% : matrix (𝑜 − 1)×(𝑜 − 1) 𝐵 !! 𝐵 !" … 𝐵 !# 𝑏 !! 𝒃 !" 𝑣 !! 𝒗 !" 𝐵 "! 𝐵 "" … 𝐵 "# 1 𝟏 = = 𝒃 "! 𝑩 "" 𝟏 𝑽 "" 𝒎 "! 𝑴 "" ⋮ ⋮ ⋱ ⋮ 𝐵 #! 𝐵 #" … 𝐵 ## 𝒃 %$ 1) First row of 𝑽 is 𝑩 %% the first row of 𝑩 𝑏 !! 𝒃 !" 𝑣 !! 𝒗 !" = 𝒃 "! 𝑩 "" 𝑣 !! 𝒎 "! 𝒎 "! 𝒗 !" + 𝑴 "" 𝑽 "" 2) 𝒎 <= = = Known! % !! 𝒃 <= 3) 𝑵 = 𝑴 "" 𝑽 "" = 𝑩 "" − 𝒎 "! 𝒗 !" First column of 𝑴 is the first column of 𝑩 / 𝑣 == Need another factorization!

Example 1) First row of 𝑽 is the first row of 𝑩 2 8 4 1 3 3 1 2 2) First column of 𝑴 is the first column of 𝑩 / 𝑣 == 𝑵 = 1 2 6 2 3) 𝑴 "" 𝑽 "" = 𝑩 "" − 𝒎 "! 𝒗 !" 1 3 4 2 1 0 0 0 2 8 4 1 0.5 0 0 0 0 0 0 0 𝑴 = 𝑽 = 0.5 0 0 0 0 0 0 0 0.5 0 0 0 0 0 0 0 2 3 3 4 2 0.5 2 6 2 4 𝑴 %% 𝑽 %% = 𝑩 %% − 𝒎 %$ 𝒗 $% = − 2 0.5 3 4 2 4 2 0.5 2 8 4 1 1 −2 1 2.5 𝑵 = 1 −2 4 1.5 1 −1 2 1.5

2 8 4 1 1 −2 1 2.5 𝑵 = 1 −2 4 1.5 1 −1 2 1.5 1 0 0 0 2 8 4 1 0.5 1 0 0 0 −2 1 2.5 𝑴 = 𝑽 = 0 0 0 0 0 0 0.5 1 0.5 0.5 0 0 0 0 0 0 𝑴 %% 𝑽 %% = 𝑩 %% − 𝒎 %$ 𝒗 $% = 4 1.5 1 2.5 1.5 − 2 0.5 1.25 2 8 4 1 1 −2 1 2.5 𝑵 = 3 −1 1 −2 1 −1 1.5 0.25

2 8 4 1 1 2.5 1 −2 𝑵 = 1 −2 3 −1 1 −1 1.5 0.25 1 0 0 0 2 8 4 1 0.5 1 0 0 0 −2 1 2.5 𝑴 = 𝑽 = 0.5 1 1 0 0 0 3 −1 0.5 0 0.5 0.5 0 0 0 0 𝑴 %% 𝑽 %% = 𝑩 %% − 𝒎 %$ 𝒗 $% = 0.25 − −0.5 = 0.75 1 0 0 0 2 8 4 1 0.5 1 0 0 0 −2 1 2.5 𝑴 = 𝑽 = 0.5 1 1 0 0 0 3 −1 0.5 1 0.5 0.5 0 0 0 0.75

𝑏 $$ : scalar LU Factorization 𝒃 $% : row vector (1×(𝑜 − 1)) 𝒃 %$ : column vector (𝑜 − 1)×1 𝑏 $$ 𝒃 $% 𝑩 %% : matrix (𝑜 − 1)×(𝑜 − 1) 𝐵 !! 𝐵 !" … 𝐵 !# 𝑏 !! 𝒃 !" 𝑣 !! 𝒗 !" 𝐵 "! 𝐵 "" … 𝐵 "# 1 𝟏 = = 𝒃 "! 𝑩 "" 𝟏 𝑽 "" 𝒎 "! 𝑴 "" ⋮ ⋮ ⋱ ⋮ 𝐵 #! 𝐵 #" … 𝐵 ## 𝒃 %$ 1) First row of 𝑽 is 𝑩 %% the first row of 𝑩 𝑏 !! 𝒃 !" 𝑣 !! 𝒗 !" = 𝒃 "! 𝑩 "" 𝑣 !! 𝒎 "! 𝒎 "! 𝒗 !" + 𝑴 "" 𝑽 "" 2) 𝒎 <= = = Known! % !! 𝒃 <= 3) 𝑵 = 𝑴 "" 𝑽 "" = 𝑩 "" − 𝒎 "! 𝒗 !" First column of 𝑴 is the first column of 𝑩 / 𝑣 == Need another factorization!

Cost of solving linear system of equations

Cost of solving triangular systems * 𝑐 & − ∑ '(&)! 𝑉 &' 𝑦 ' 𝑦 * = 𝑐 * /𝑉 ** 𝑦 & = , 𝑗 = 𝑜 − 1, 𝑜 − 2, … , 1 𝑉 &&