lu                  package:Matrix                  R Documentation

_T_r_i_a_n_g_u_l_a_r _D_e_c_o_m_p_o_s_i_t_i_o_n _o_f _a _S_q_u_a_r_e _M_a_t_r_i_x

_D_e_s_c_r_i_p_t_i_o_n:

     Computes triangular decompositions of square matrices.

_U_s_a_g_e:

     lu(x, ...)
     ## S4 method for signature 'dgeMatrix':
     lu(x, warnSing = TRUE, ...)
     ## S4 method for signature 'dgCMatrix':
     lu(x, errSing = TRUE, ...)

_A_r_g_u_m_e_n_t_s:

       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.

_D_e_t_a_i_l_s:

     '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.

_V_a_l_u_e:

     An object of class '"LU"', i.e., '"denseLU"' or '"sparseLU"', see
     'sparseLU'; this is a representation of a triangular decomposition
     of 'x'.

_R_e_f_e_r_e_n_c_e_s:

     Golub, G., and Van Loan, C. F. (1989). _Matrix Computations,_ 2nd
     edition, Johns Hopkins, Baltimore.

     Tim Davis (2005) <URL:
     http://www.cise.ufl.edu/research/sparse/CSparse/>

     Timothy A. Davis (2006) _Direct Methods for Sparse Linear
     Systems_, SIAM Series Fundamentals of Algorithms.

_S_e_e _A_l_s_o:

     Class definitions 'LU' and 'sparseLU' and function 'expand'; 'qr',
     'chol'.

_E_x_a_m_p_l_e_s:

     ##--- 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

