Tridiagonal matrix algorithm  TDMA (Thomas algorithm)
From CFDWiki
(Added more material) 
(→Variants) 

(16 intermediate revisions not shown)  
Line 30:  Line 30:  
\left[  \left[  
\begin{matrix}  \begin{matrix}  
  {  +  {d_1 } \\ 
  {  +  {d_2 } \\ 
\cdot \\  \cdot \\  
\cdot \\  \cdot \\  
  {  +  {d_n } \\ 
\end{matrix}  \end{matrix}  
\right].  \right].  
</math>  </math>  
  For such systems, the solution can be obtained in <math>O(n)</math> operations instead of <math>O(n^3)</math> required by Gaussian Elimination. A first sweep eliminates the <math>a_i</math>'s, and then an (abbreviated) backward substitution produces the solution. Example of such matrices commonly arise from the discretization of 1D problems (e.g. the 1D  +  For such systems, the solution can be obtained in <math>O(n)</math> operations instead of <math>O(n^3)</math> required by Gaussian Elimination. A first sweep eliminates the <math>a_i</math>'s, and then an (abbreviated) backward substitution produces the solution. Example of such matrices commonly arise from the discretization of 1D problems (e.g. the 1D Poisson problem). 
== Algorithm ==  == Algorithm ==  
+  
+  The following algorithm performs the TDMA, overwriting the original arrays. In some situations this is not desirable, so some prefer to copy the original arrays beforehand.  
+  
Forward elimination phase  Forward elimination phase  
: for k = 2 step until n do <br>  : for k = 2 step until n do <br>  
:: <math>m = {{a_k } \over {b_{k  1} }} </math> <br>  :: <math>m = {{a_k } \over {b_{k  1} }} </math> <br>  
  :: <math> b_k  +  :: <math> b_k = b_k  mc_{k  1} </math> <br> 
  :: <math> d_k  +  :: <math> d_k = d_k  md_{k  1} </math> <br> 
: end loop (k)  : end loop (k)  
Backward substitution phase  Backward substitution phase  
  : <math> x_n = {{d_n  +  : <math> x_n = {{d_n } \over {b_n }} </math> <br> 
: for k = n1 stepdown until 1 do <br>  : for k = n1 stepdown until 1 do <br>  
  :: <math> x_k = {{d_k  +  :: <math> x_k = {{d_k  c_k x_{k + 1} } \over {b_k }} </math> <br> 
: end loop (k)  : end loop (k)  
+  
+  == Assessments ==  
+  
+  This algorithm is only applicable to matrices that are diagonally dominant, which is to say  
+  
+  :<math>\left  b_i \right \vert > \left  a_i \right \vert + \left  c_{i} \right \vert \quad \quad i \in {1,...,n} </math>  
== Variants ==  == Variants ==  
Line 67:  Line 76:  
</math> <br>  </math> <br>  
  In  +  In matrix form, this is 
  In other  +  :<math> 
+  \left[  
+  \begin{matrix}  
+  {b_1} & {c_1} & { } & { } & {a_1} \\  
+  {a_2} & {b_2} & {c_2} & { } & { } \\  
+  { } & {a_3} & {b_3} & \cdot & { } \\  
+  { } & { } & \cdot & \cdot & {c_{n1}}\\  
+  {c_n} & { } & { } & {a_n} & {b_n}\\  
+  \end{matrix}  
+  \right]  
+  \left[  
+  \begin{matrix}  
+  {x_1 } \\  
+  {x_2 } \\  
+  \cdot \\  
+  \cdot \\  
+  {x_n } \\  
+  \end{matrix}  
+  \right]  
+  =  
+  \left[  
+  \begin{matrix}  
+  {d_1 } \\  
+  {d_2 } \\  
+  \cdot \\  
+  \cdot \\  
+  {d_n } \\  
+  \end{matrix}  
+  \right].  
+  </math>  
+  
+  In this case, we can make use of the ShermanMorrison formula to avoid the additional operations of Gaussian elimination and still use the Thomas algorithm. We are now solving a problem of the form  
+  
+  :<math>(A' + uv^T)x = d</math>  
+  
+  where  
+  
+  :<math> u^T = [b_1\ 0\ 0\ ...\ 0\ c_n],\ v^T = [1\ 0\ 0\ ...\ 0\ a_1/b_1].</math>  
+  
+  <math>A'</math> is a slightly different tridiagonal system than above, and the solution to the perturbed system is obtained by solving  
+  
+  :<math>A'y=d,\ A'q=u</math>  
+  
+  and compute <math>x</math> as  
+  
+  :<math>x = y  \{(v^T y)/(1 + (v^T q))\} q</math>  
+  
+  In other situations, the system of equations may be ''block tridiagonal'', with smaller submatrices arranged as the individual elements in the above matrix system. Simplified forms of Gaussian elimination have been developed for these situations.  
== References ==  == References ==  
+  #{{referencebookauthor=Thomas, L.H.year=1949title=Elliptic Problems in Linear Differential Equations over a Networkrest=Watson Sci. Comput. Lab Report, Columbia University, New York.}}  
#{{referencebookauthor=Conte, S.D., and deBoor, C.year=1972title=Elementary Numerical Analysisrest= McGrawHill, New York.}}  #{{referencebookauthor=Conte, S.D., and deBoor, C.year=1972title=Elementary Numerical Analysisrest= McGrawHill, New York.}}  
  +  ''Still TODO: Add more references, more on the variants, make things nicer looking, and maybe more performance type info''  
  ''TODO: Add more references, more on the variants, and maybe more performance type info''  +  [[User:JasondJasond]] 22:59, 7 April 2006 (MDT) 
  [[User:JasondJasond]]  + 
Latest revision as of 12:24, 26 November 2012
Contents 
Introduction
The tridiagonal matrix algorithm (TDMA), also known as the Thomas algorithm, is a simplified form of Gaussian elimination that can be used to solve tridiagonal systems of equations. A tridiagonal system may be written as
where and . In matrix form, this system is written as
For such systems, the solution can be obtained in operations instead of required by Gaussian Elimination. A first sweep eliminates the 's, and then an (abbreviated) backward substitution produces the solution. Example of such matrices commonly arise from the discretization of 1D problems (e.g. the 1D Poisson problem).
Algorithm
The following algorithm performs the TDMA, overwriting the original arrays. In some situations this is not desirable, so some prefer to copy the original arrays beforehand.
Forward elimination phase
 for k = 2 step until n do

 end loop (k)
Backward substitution phase

 for k = n1 stepdown until 1 do

 end loop (k)
Assessments
This algorithm is only applicable to matrices that are diagonally dominant, which is to say
Variants
In some situations, particularly those involving periodic boundary conditions, a slightly perturbed form of the tridiagonal system may need to be solved:
In matrix form, this is
In this case, we can make use of the ShermanMorrison formula to avoid the additional operations of Gaussian elimination and still use the Thomas algorithm. We are now solving a problem of the form
where
is a slightly different tridiagonal system than above, and the solution to the perturbed system is obtained by solving
and compute as
In other situations, the system of equations may be block tridiagonal, with smaller submatrices arranged as the individual elements in the above matrix system. Simplified forms of Gaussian elimination have been developed for these situations.
References
 Thomas, L.H. (1949), Elliptic Problems in Linear Differential Equations over a Network, Watson Sci. Comput. Lab Report, Columbia University, New York..
 Conte, S.D., and deBoor, C. (1972), Elementary Numerical Analysis, McGrawHill, New York..
Still TODO: Add more references, more on the variants, make things nicer looking, and maybe more performance type info Jasond 22:59, 7 April 2006 (MDT)