LU decomposition
From CFDWiki
(towards a uniform notation for linear systems : A*Phi = B) 
(Added material, used \vec for vectors, still needs more factorization info) 

Line 1:  Line 1:  
  ==  +  == Description == 
  +  Consider the system of equations <math>A\vec{x}=\vec{b}</math>, where <math>A</math> is an <math>n\times n</math> nonsingular matrix. <math>A</math> may be decomposed into an lower triangular part <math>L</math> and an upper triangular part <math>U</math> that will lead us to a direct procedure for the solution of the original system. This decomposition procedure is quite useful when more than one righthand side (more than one <math>\vec{b}</math>) is to be used.  
  +  
  +  
  +  The algorithm is relatively straightforward  first, we determine the upper and lower triangular parts:  
  +  :<math>A = LU.</math>  
  +  
  +  
  +  
  +  
  +  
  {{  +  Then, 
    +  
+  :<math> A \vec{x} = (LU) \vec{x} = L(U\vec{x})=L\vec{y},</math>  
+  
+  where <math>y=Ux</math>. Once we solve the system  
+  
+  :<math>L\vec{y} = \vec{b},</math>  
+  
+  we will be able to find the solution to the original system by solving  
+  
+  :<math>U\vec{x} = \vec{b}</math>  
+  
+  The first solution is a foward substitution, while the second solution is a backward substitution. Both can be done efficiently once the factorization is available. The forward substitution may be expressed as  
+  
+  :<math>  
+  y_i = {1 \over {L_{ii} }}\left( {b_i  \sum\limits_{j = 1}^{i1} {L_{ij} y_j } } \right),  
+  i = 1,2,\ldots,n</math>  
+  
+  and the backward substitution process may be expressed as  
+  
+  :<math>  
+  x_i = {1 \over {U_{ii} }}\left( {y_i  \sum\limits_{j = i + 1}^n {U_{ij} x_j } } \right),  
+  i = n,n1,\ldots,1.</math>  
+  
+  == Algorithm ==  
+  ''Add factorization here''<br>  
+  Forward substitution <br>  
+  : for k:=1 step until n do <br>  
+  :: for i:=1 step until k1 <br>  
+  ::: <math>b_k=b_kL_{ki}b_{i}</math> <br>  
+  :: end loop (i) <br>  
+  :: <math>b_{k}=b_{k}/L_{kk}</math><br>  
+  : end loop (k)  
+  Backward substitution  
+  : for k:=n stepdown until 1 do <br>  
+  :: for i:=k+1 step until n <br>  
+  ::: <math>b_k=b_kU_{ki}b_{i}</math> <br>  
+  :: end loop (i) <br>  
+  :: <math>x_{k}=b_{k}/U_{kk}</math><br>  
+  : end loop (k)  
+  == Important Considerations ==  
+  As with [[Gaussian elimination]], LU decomposition is probably best used for relatively small, relatively nonsparse systems of equations (with small and nonsparse open to some interpretation). For larger and/or sparse problems, it would probably be best to either use an iterative method or use a direct solver package (e.g. [http://www.cse.psu.edu/~raghavan/software.html DSCPACK]) as opposed to writing one of your own.  
    +  If one has a single lefthandside matrix and many righthand side vectors, then LU decomposition would be a good solution procedure to consider. At the very least, it should be faster than solving each system separately with Gaussian elimination. 
  + 
Revision as of 01:10, 17 December 2005
Description
Consider the system of equations , where is an nonsingular matrix. may be decomposed into an lower triangular part and an upper triangular part that will lead us to a direct procedure for the solution of the original system. This decomposition procedure is quite useful when more than one righthand side (more than one ) is to be used.
The algorithm is relatively straightforward  first, we determine the upper and lower triangular parts:
Then,
where . Once we solve the system
we will be able to find the solution to the original system by solving
The first solution is a foward substitution, while the second solution is a backward substitution. Both can be done efficiently once the factorization is available. The forward substitution may be expressed as
and the backward substitution process may be expressed as
Algorithm
Add factorization here
Forward substitution
 for k:=1 step until n do
 for i:=1 step until k1

 end loop (i)

 for i:=1 step until k1
 end loop (k)
Backward substitution
 for k:=n stepdown until 1 do
 for i:=k+1 step until n

 end loop (i)

 for i:=k+1 step until n
 end loop (k)
Important Considerations
As with Gaussian elimination, LU decomposition is probably best used for relatively small, relatively nonsparse systems of equations (with small and nonsparse open to some interpretation). For larger and/or sparse problems, it would probably be best to either use an iterative method or use a direct solver package (e.g. DSCPACK) as opposed to writing one of your own.
If one has a single lefthandside matrix and many righthand side vectors, then LU decomposition would be a good solution procedure to consider. At the very least, it should be faster than solving each system separately with Gaussian elimination.