Dear Doug and Simon, The U matrix in the LU decomposition is upper triangular, but (as I understand it) it's not in row-echelon because the leading entries aren't 1 and the columns with leading entries aren't fully swept.
I'm certainly not making any claims, other than pedagogical, for the REF. In fact, the function that I sent originated as a demonstration for a class. By the way, I think that setting the tolerance in the function to Machine$double.eps (which is what I did) is too small: Gabor Grothendieck sent me an example in which my function failed with the default tolerance. Regards, John > -----Original Message----- > From: Douglas Bates [mailto:[EMAIL PROTECTED] On > Behalf Of Douglas Bates > Sent: Friday, March 05, 2004 8:26 AM > To: Simon Fear > Cc: John Fox; Spencer Graves; [EMAIL PROTECTED] > Subject: Re: [R] row-echelon form (was no subject) > > "Simon Fear" <[EMAIL PROTECTED]> writes: > > > I think one needs an LU decomposition rather than QR. > > However, I couldn't find anything off the shelf to do an LU, other > > than learning that determinant() now uses LU instead of QR > or SVD, so > > the code to do it must be in there for those that want it. > > I was about to write that the version of the Matrix package > for r-devel (to be R-1.9.0) has an LU decomposition but then > I checked and found that, although the underlying C code is > there, I haven't written an accessor function. In any case > the LU decomposition is computed and stored whenever it is > needed, say for calculating a condition number or for solving > a linear system. > > > library(Matrix) # version 0.7-3 for R-1.9.0 > > mm = Matrix(rnorm(9), nc = 3) > > mm > [,1] [,2] [,3] > [1,] -0.4186151 -0.3959071 -0.1717103 > [2,] -0.5334526 -2.2817253 -1.5398915 > [3,] -0.5993042 -1.3396313 0.2062784 > > rcond(mm) > [1] 0.04505896 > > str(mm) > list() > - attr(*, "x")= num [1:9] -0.419 -0.533 -0.599 -0.396 -2.282 ... > - attr(*, "Dim")= int [1:2] 3 3 > - attr(*, "rcond")= Named num 0.0451 > ..- attr(*, "names")= chr "O" > - attr(*, "factorization")=List of 1 > ..$ LU: list() > .. ..- attr(*, "x")= num [1:9] -0.599 0.890 0.699 -1.340 > -1.089 ... > .. ..- attr(*, "Dim")= int [1:2] 3 3 > .. ..- attr(*, "pivot")= int [1:3] 3 2 3 > .. ..- attr(*, "class")= atomic [1:1] LU > .. .. ..- attr(*, "package")= chr "Matrix" > - attr(*, "class")= atomic [1:1] geMatrix > ..- attr(*, "package")= chr "Matrix" > > expand([EMAIL PROTECTED]) > $L > [,1] [,2] [,3] > [1,] 1.0000000 0.0000000 0 > [2,] 0.8901200 1.0000000 0 > [3,] 0.6985019 -0.4955766 1 > > $U > [,1] [,2] [,3] > [1,] -0.5993042 -1.339631 0.2062784 > [2,] 0.0000000 -1.089293 -1.7235040 > [3,] 0.0000000 0.000000 -1.1699244 > > As I understand it the LU decomposition is much more widely > used in numerical linear algebra than is the row-reduced > echelon form. ______________________________________________ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
