lu {Matrix} | R Documentation |
Computes triangular decompositions of square matrices.
lu(x, ...) ## S4 method for signature 'dgeMatrix': lu(x, warnSing = TRUE, ...) ## S4 method for signature 'dgCMatrix': lu(x, errSing = TRUE, ...)
x |
a matrix of square dimension. No missing values or IEEE special values are allowed. |
warnSing |
(when x is a
"denseMatrix" ) logical specifying if a warning
should be signalled when x is singular. |
errSing |
(when x is a
"sparseMatrix" ) logical specifying if an error
(see stop ) should be signalled when x is
singular. When x is singular, lu(x, errSing=FALSE)
returns NA instead of an LU decomposition. No
warning is signalled and the useR should be careful in that case.
|
... |
further arguments passed to or from other methods. |
lu()
is a generic function with special methods for different types
of matrices. Use showMethods("lu")
to list all the methods
for the lu
generic.
The method for class dgeMatrix
(and all dense
matrices) is based on LAPACK's "dgetrf"
subroutine. It returns
a decomposition also for singular matrices.
The method for class dgCMatrix
(and all sparse
matrices) is based on functions from the CSparse library. It signals
an error (or returns NA
, when errSing = FALSE
, see
above) when the decomposition algorithm fails, as when x
is
(too close to) singular.
An object of class "LU"
, i.e., "denseLU"
or "sparseLU"
, see sparseLU
; this is
a representation of a triangular decomposition of x
.
Golub, G., and Van Loan, C. F. (1989). Matrix Computations, 2nd edition, Johns Hopkins, Baltimore.
Tim Davis (2005) http://www.cise.ufl.edu/research/sparse/CSparse/
Timothy A. Davis (2006) Direct Methods for Sparse Linear Systems, SIAM Series “Fundamentals of Algorithms”.
Class definitions LU
and sparseLU
and function expand
;
qr
, chol
.
##--- Dense ------------------------- x <- Matrix(rnorm(9), 3, 3) lu(x) ##--- Sparse ------------------------ pm <- as(readMM(system.file("external/pores_1.mtx", package = "Matrix")), "CsparseMatrix") str(pmLU <- lu(pm)) # p is a 0-based permutation of the rows # q is a 0-based permutation of the columns ## permute rows and columns of original matrix ppm <- pm[pmLU@p + 1L, pmLU@q + 1L] pLU <- pmLU@L %*% pmLU@U ## equal up to "rounding" ppm[1:14, 1:5] pLU[1:14, 1:5] # product can have extra zeros ## "prove" consistency (up to rounding): i0 <- ppm != pLU & ppm == 0 iN <- ppm != pLU & ppm != 0 stopifnot(all(abs((ppm - pLU)[i0]) < 1e-7), # absolute error for true 0 all(abs((ppm - pLU)[iN]/ppm[iN]) < 1e-9)) # relative error