Matrix (version 1.7-1)

lu-methods: Methods for LU Factorization

Description

Computes the pivoted LU factorization of an \(m \times n\) real matrix \(A\), which has the general form $$P_{1} A P_{2} = L U$$ or (equivalently) $$A = P_{1}' L U P_{2}'$$ where \(P_{1}\) is an \(m \times m\) permutation matrix, \(P_{2}\) is an \(n \times n\) permutation matrix, \(L\) is an \(m \times \min(m,n)\) unit lower trapezoidal matrix, and \(U\) is a \(\min(m,n) \times n\) upper trapezoidal matrix.

Methods for denseMatrix are built on LAPACK routine dgetrf, which does not permute columns, so that \(P_{2}\) is an identity matrix.

Methods for sparseMatrix are built on CXSparse routine cs_lu, which requires \(m = n\), so that \(L\) and \(U\) are triangular matrices.

Usage

lu(x, ...)
# S4 method for dgeMatrix
lu(x, warnSing = TRUE, ...)
# S4 method for dgCMatrix
lu(x, errSing = TRUE, order = NA_integer_,
  tol = 1, ...)
# S4 method for dsyMatrix
lu(x, cache = TRUE, ...)
# S4 method for dsCMatrix
lu(x, cache = TRUE, ...)
# S4 method for matrix
lu(x, ...)

Value

An object representing the factorization, inheriting from virtual class LU. The specific class is denseLU unless x inherits from virtual class sparseMatrix, in which case it is sparseLU.

Arguments

x

a finite matrix or Matrix to be factorized, which must be square if sparse.

warnSing

a logical indicating if a warning should be signaled for singular x. Used only by methods for dense matrices.

errSing

a logical indicating if an error should be signaled for singular x. Used only by methods for sparse matrices.

order

an integer in 0:3 passed to CXSparse routine cs_sqr, indicating a strategy for choosing the column permutation \(P_{2}\). 0 means no column permutation. 1, 2, and 3 indicate a fill-reducing ordering of \(A + A'\), \(\tilde{A}' \tilde{A}\), and \(A' A\), where \(\tilde{A}\) is \(A\) with “dense” rows removed. NA (the default) is equivalent to 2 if tol == 1 and 1 otherwise. Do not set to 0 unless you know that the column order of \(A\) is already sensible.

tol

a number. The original pivot element is used if its absolute value exceeds tol * a, where a is the maximum in absolute value of the other possible pivot elements. Set tol < 1 only if you know what you are doing.

cache

a logical indicating if the result should be cached in x@factors[["LU"]]. Note that caching is experimental and that only methods for classes extending generalMatrix or symmetricMatrix will have this argument.

...

further arguments passed to or from methods.

Details

What happens when x is determined to be near-singular differs by method. The method for class dgeMatrix completes the factorization, warning if warnSing = TRUE and in any case returning a valid denseLU object. Users of this method can detect singular x with a suitable warning handler; see tryCatch. In contrast, the method for class dgCMatrix abandons further computation, throwing an error if errSing = TRUE and otherwise returning NA. Users of this method can detect singular x with an error handler or by setting errSing = FALSE and testing for a formal result with is(., "sparseLU").

References

The LAPACK source code, including documentation; see https://netlib.org/lapack/double/dgetrf.f.

Davis, T. A. (2006). Direct methods for sparse linear systems. Society for Industrial and Applied Mathematics. tools:::Rd_expr_doi("10.1137/1.9780898718881")

Golub, G. H., & Van Loan, C. F. (2013). Matrix computations (4th ed.). Johns Hopkins University Press. tools:::Rd_expr_doi("10.56021/9781421407944")

See Also

Classes denseLU and sparseLU and their methods.

Classes dgeMatrix and dgCMatrix.

Generic functions expand1 and expand2, for constructing matrix factors from the result.

Generic functions Cholesky, BunchKaufman, Schur, and qr, for computing other factorizations.

Examples

 
library(stats, pos = "package:base", verbose = FALSE)
library(utils, pos = "package:base", verbose = FALSE)

showMethods("lu", inherited = FALSE)
set.seed(0)

## ---- Dense ----------------------------------------------------------

(A1 <- Matrix(rnorm(9L), 3L, 3L))
(lu.A1 <- lu(A1))

(A2 <- round(10 * A1[, -3L]))
(lu.A2 <- lu(A2))

## A ~ P1' L U in floating point
str(e.lu.A2 <- expand2(lu.A2), max.level = 2L)
stopifnot(all.equal(A2, Reduce(`%*%`, e.lu.A2)))

## ---- Sparse ---------------------------------------------------------

A3 <- as(readMM(system.file("external/pores_1.mtx", package = "Matrix")),
         "CsparseMatrix")
(lu.A3 <- lu(A3))

## A ~ P1' L U P2' in floating point
str(e.lu.A3 <- expand2(lu.A3), max.level = 2L)
stopifnot(all.equal(A3, Reduce(`%*%`, e.lu.A3)))