Gaussian elimination figures prominently in the matrixDemos package, and is implemented directly in the GaussianElimination function. Gaussian elimination is also used by matrixInverse to compute the inverse of a nonsingular square matrix; by Ginv to compute a generalized inverse of a singular square matrix or a rectangular matrix; and by RREF to compute the reduced row-echelon form of a matrix.

GaussianElimination has one required argument, the matrix A, which is transformed to reduced row-echelon form; for example:

> (A <- matrix(c(1,2,3,4,5,6,7,8,10), 3, 3)) # a nonsingular matrix
     [,1] [,2] [,3]
[1,]    1    4    7
[2,]    2    5    8
[3,]    3    6   10
> GaussianElimination(A)
     [,1] [,2] [,3]
[1,]    1    0    0
[2,]    0    1    0
[3,]    0    0    1
> (B <- matrix(1:9, 3, 3)) # a singular matrix
     [,1] [,2] [,3]
[1,]    1    4    7
[2,]    2    5    8
[3,]    3    6    9
> GaussianElimination(B)
     [,1] [,2] [,3]
[1,]    1    0   -1
[2,]    0    1    2
[3,]    0    0    0

In these applications, GaussianElimination is equivalent to RREF.

A second optional argument to GaussianElimination is B, which may be a matrix or a vector, and which is transformed in the same manner as A. For example,

> A # nonsingular
     [,1] [,2] [,3]
[1,]    1    4    7
[2,]    2    5    8
[3,]    3    6   10
> GaussianElimination(A, diag(3)) # find inverse of A
     [,1] [,2] [,3]       [,4]       [,5] [,6]
[1,]    1    0    0 -0.6666667 -0.6666667    1
[2,]    0    1    0 -1.3333333  3.6666667   -2
[3,]    0    0    1  1.0000000 -2.0000000    1
> b <- 1:3
> GaussianElimination(A, b)  # solve the matrix equation Ax = b
     [,1] [,2] [,3] [,4]
[1,]    1    0    0    1
[2,]    0    1    0    0
[3,]    0    0    1    0

For the last example, the solution is \(x_1 = 1, x_2 = x_3 = 0\).

Other arguments to GaussianElimination:

Some further examples:

> matrixInverse(A, fractions=TRUE)
     [,1] [,2] [,3]
[1,] -2/3 -2/3    1
[2,] -4/3 11/3   -2
[3,]    1   -2    1
> RREF(A, verbose=TRUE)
     [,1] [,2]     [,3]
[1,]    1    2 3.333333
[2,]    0    1 1.333333
[3,]    0    2 3.666667
     [,1] [,2]       [,3]
[1,]    1    0 -0.3333333
[2,]    0    1  1.8333333
[3,]    0    0 -0.5000000
     [,1] [,2] [,3]
[1,]    1    0    0
[2,]    0    1    0
[3,]    0    0    1
     [,1] [,2] [,3]
[1,]    1    0    0
[2,]    0    1    0
[3,]    0    0    1
> Ginv(B, fractions=TRUE)  # a generalized inverse of singular B
     [,1]  [,2]  [,3] 
[1,]  -3/4     0  7/12
[2,]     0     0     0
[3,]   1/4     0 -1/12
> B %*% Ginv(B) %*% B   # check
     [,1] [,2] [,3]
[1,]    1    4    7
[2,]    2    5    8
[3,]    3    6    9
> (C <- matrix(1:12, 3, 4)) # rectangular matrix
     [,1] [,2] [,3] [,4]
[1,]    1    4    7   10
[2,]    2    5    8   11
[3,]    3    6    9   12
> Ginv(C, fractions=TRUE)
     [,1]  [,2]  [,3] 
[1,]  -2/3     0   5/9
[2,]     0     0     0
[3,]     0     0     0
[4,]   1/6     0 -1/18
> round(C %*% Ginv(C) %*% C, 5)   # check
     [,1] [,2] [,3] [,4]
[1,]    1    4    7   10
[2,]    2    5    8   11
[3,]    3    6    9   12