Examples for 'Matrix::sparseLU-class'


Sparse LU decomposition of a square sparse matrix

Aliases: sparseLU-class

Keywords: classes

### ** Examples

## Extending the one in   examples(lu), calling the matrix  A,
## and confirming the factorization identities :
A <- as(readMM(system.file("external/pores_1.mtx",
                            package = "Matrix")),
         "CsparseMatrix")
## with dimnames(.) - to see that they propagate to L, U :
dimnames(A) <- list(paste0("r", seq_len(nrow(A))),
                    paste0("C", seq_len(ncol(A))))
str(luA <- lu(A)) # p is a 0-based permutation of the rows
Formal class 'sparseLU' [package "Matrix"] with 5 slots
  ..@ L  :Formal class 'dtCMatrix' [package "Matrix"] with 7 slots
  .. .. ..@ i       : int [1:141] 0 1 18 22 25 26 1 18 22 25 ...
  .. .. ..@ p       : int [1:31] 0 6 11 13 16 19 22 26 31 35 ...
  .. .. ..@ Dim     : int [1:2] 30 30
  .. .. ..@ Dimnames:List of 2
  .. .. .. ..$ : chr [1:30] "r2" "r12" "r30" "r28" ...
  .. .. .. ..$ : NULL
  .. .. ..@ x       : num [1:141] 1 -0.993819 -0.004979 -0.000132 0.000132 ...
  .. .. ..@ uplo    : chr "L"
  .. .. ..@ diag    : chr "N"
  ..@ U  :Formal class 'dtCMatrix' [package "Matrix"] with 7 slots
  .. .. ..@ i       : int [1:206] 0 0 1 2 2 3 2 4 2 3 ...
  .. .. ..@ p       : int [1:31] 0 1 3 4 6 8 12 15 20 25 ...
  .. .. ..@ Dim     : int [1:2] 30 30
  .. .. ..@ Dimnames:List of 2
  .. .. .. ..$ : NULL
  .. .. .. ..$ : chr [1:30] "C1" "C2" "C30" "C28" ...
  .. .. ..@ x       : num [1:206] -7178502 -24613411 -18311731 -6399179 715 ...
  .. .. ..@ uplo    : chr "U"
  .. .. ..@ diag    : chr "N"
  ..@ p  : int [1:30] 1 11 29 27 19 28 9 18 17 7 ...
  ..@ q  : int [1:30] 0 1 29 27 19 28 9 18 17 8 ...
  ..@ Dim: int [1:2] 30 30
                  # q is a 0-based permutation of the columns
xA <- expand(luA)
## which is simply doing
stopifnot(identical(xA$ L, luA@L),
          identical(xA$ U, luA@U),
          identical(xA$ P, as(luA@p +1L, "pMatrix")),
          identical(xA$ Q, as(luA@q +1L, "pMatrix")))

P.LUQ <- with(xA, t(P) %*% L %*% U %*% Q)
stopifnot(all.equal(unname(A), unname(P.LUQ), tolerance = 1e-12))

## permute rows and columns of original matrix
pA <- A[luA@p + 1L, luA@q + 1L]
PAQ. <- with(xA, P %*% A %*% t(Q))
stopifnot(all.equal(unname(pA), unname(PAQ.), tolerance = 1e-12))

pLU <- drop0(luA@L %*% luA@U) # L %*% U -- dropping extra zeros
stopifnot(all.equal(pA, pLU, tolerance = 1e-12))

[Package Matrix version 1.5-3 Index]