# [R] Solving Matrices

Fri Apr 16 18:09:33 CEST 2004

```"Richard" writes:

[ ... ]
>
> So here are three different functions all saying much the same thing
> about x: it is numerically close to a matrix which does NOT span the
> whole of R**4.

Theorem 4 (in my text):

Let A be an m x n matrix. Then the following statements are
logically equivalent. That is, for a particular A, either they
are all true statements or they are all false:

a) For each b in R**m, the equation Ax = b has a solution.
b) The columns of A span R**m.
c) A has a pivot position in every row.

The back of the book is in agreement with you but I can manipulate the
matrix to produce (c) above which should mean that (b) is correct:

> x <- matrix(data=c(7,-5,6,-7,2,-3,10,9,-5,4,-2,2,8,-9,7,15),nrow=4)
> rowEchelonForm(x)
[,1] [,2] [,3] [,4]
[1,]    1    0    0    0
[2,]    0    1    0    0
[3,]    0    0    1    0
[4,]    0    0    0    1

Likewise, I can produce the above using the swap(), scale(), gauss()
and bgauss() functions that the author of the text originally wrote
for maple, Mathematica, etc. but from all I can tell, the row echelon
form of matrix x meets the criteria of theorem 4 above, yet the columns
of x do not span R**4.

Using the gauss() function:

[,1]      [,2]       [,3]          [,4]
[1,]    7  2.000000 -5.0000000  8.000000e+00
[2,]    0 -1.571429  0.4285714 -3.285714e+00
[3,]    0  0.000000  4.5454545 -1.718182e+01
[4,]    0  0.000000  0.0000000 -5.035972e-15

Again, this has a pivot column in each row although [4,4] seems a bit
suspicious.

Elizabeth

:

gauss <- function(A, r, v = c()) {
m <- dim(A)
n <- dim(A)
if(length(r) > 1 | r < 1 | r > m)
stop("The second entry in gauss(...) must be the pivot row index.")
if(length(v) == 0) {
if(r+1 > m)
v <- m
else
v <- (r+1):m
}
for(j in v) {
if(j == r)
stop("Pivot ow cannot change itself.")
}

col <- 1
while(A[r, col] == 0 && col < n)
col <- col + 1
if(col == n && A[r, col] == 0)
stop("Row r has only zeroes and cannot be a pivot row.")
for(j in v) {
multiplier <- -A[j, col]/A[r, col]
A[j,] <- A[j,] + multiplier * A[r,]
}
A
}

```