Triangular Solvers
Linear systems A * x == b are easy to solve if A is an upper/lower triangular matrix. We provide a generic function to perform this operation:
x= upper_trisolve(A, b);
The matrix A must be triangular matrix otherwise the function can throw an exception. (For dense matrices the lower part is currently ignored but this might change in the future.) If A has a unit diagonal, the diagonal entries can and must be omitted if the system is solved by:
x= unit_upper_trisolve(A, b);
The implicit diagonal decreases the stress on the memory bandwidth and avoids expensive divisions.
On matrices with non-unit diagonals, the divisions can be circumvented by inverting the diagonal once with invert_diagonal(A) and then using:
x= inverse_upper_trisolve(A, b);
Especially if A is used as preconditioner of an iterative method, the substitution of divisions by multiplications can lead to a significant speed-up.
Likewise, the functions for lower triangular matrices are defined:
x= lower_trisolve(A, b); x= unit_lower_trisolve(A, b); x= inverse_lower_trisolve(A, b);
Return to Other Matrix Functions Table of Content Proceed to Introduction Krylov-Subspace Methods
Triangular Solvers -- MTL 4 -- Peter Gottschling and Andrew Lumsdaine
-- Gen. with
rev. 7542
on Sat Aug 11 2012 by doxygen 1.7.6.1 -- © 2010 by SimuNova UG.