pracma (version 1.8.8)

lu: LU Matrix Factorization

Description

LU decomposition of a positive definite matrix as Gaussian factorization.

Usage

lu(A, scheme = c("kji", "jki", "ijk"))
lufact(A)
lusys(A, b)

Arguments

A

square positive definite numeric matrix (will not be checked).

scheme

order of row and column operations.

b

right hand side of a linear system of equations.

Value

• lu returns a list with components L and U, the two lower and upper triangular matrices such that A=L%*%U.

  lufact returns a list with L and U combined into one matrix LU, the rows used in partial pivoting, and det representing the determinant of A. See the examples how to extract matrices L and U from LU.

  lusys returns the solution of the system as a column vector.

Details

For a given matrix A, the LU decomposition exists and is unique iff its principal submatrices of order i=1,...,n-1 are nonsingular. The procedure here is a simple Gauss elimination with or without pivoting.

The scheme abbreviations refer to the order in which the cycles of row- and column-oriented operations are processed. The ``ijk'' scheme is one of the two compact forms, here the Doolite factorization (the Crout factorization would be similar).

lufact applies partial pivoting (along the rows). lusys uses LU factorization to solve the linear system A*x=b.

References

Quarteroni, A., R. Sacco, and F. Saleri (2007). Numerical Mathematics. Second edition, Springer-Verlag, Berlin Heidelberg.

J.H. Mathews and K.D. Fink (2003). Numerical Methods Using MATLAB. Fourth Edition, Pearson (Prentice-Hall), updated 2006.

See Also

qr

Examples

A <- magic(5)
D <- lu(A, scheme = "ijk")     # Doolittle scheme
D$L %*% D$U
##      [,1] [,2] [,3] [,4] [,5]
## [1,]   17   24    1    8   15
## [2,]   23    5    7   14   16
## [3,]    4    6   13   20   22
## [4,]   10   12   19   21    3
## [5,]   11   18   25    2    9

H4 <- hilb(4)
lufact(H4)$det
## [1] 0.0000001653439

x0 <- c(1.0, 4/3, 5/3, 2.0)
b  <- H4 %*% x0
lusys(H4, b)
##          [,1]
## [1,] 1.000000
## [2,] 1.333333
## [3,] 1.666667
## [4,] 2.000000