# R functions to calculate row-reduced echelon form of an m x n matrix A.
# Written by S. Spiriti, October 2004 

swap <- function(A, row1, row2)
{
     temp <- A[row1, ]
     A[row1, ] <- A[row2, ]
     A[row2, ] <- temp
     return(A)
}

mult <- function (A, const, row1)
{
     A[row1, ] <- const*A[row1, ]
     return(A)
}

mult.add <- function(A, const, row1, row2)
{
     A[row2, ] <- const*A[row1, ] + A[row2, ]
     return(A)
}

pivot <- function(A, row, col, m, n)
{
      pcol <- A[row:m, col]
      row.index <- which.max(abs(pcol)) - 1 + row
      A <- swap(A, row, row.index)
      return(A)
}

RREF <- function(A, m, n)
{
    i <- 1
    j <- 1
    while ((i <= m)  & (j <= n))
    {
        A <- pivot(A, i, j, m, n)
        if (A[i, j] == 0)
            j <- j + 1
        else
        {
            piv.entry <- A[i, j]
            A <- mult(A, 1/piv.entry, i)
            for (k in 1:m)
            {
                if (k != i)
                   A <- mult.add(A, -A[k, j], i, k)
                A <- zapsmall(A)
            }
            print(A)
            i <- i + 1
            j <- j + 1
        }
     }
     print("Row-reduced echelon form of A: ")
     print(A)
}