Redução Gauss-Jordan

Universidade Federal da Paraíba

Programa de Pós-Graduação em Economia

#------------------------------------------------------------
# Função para redução Gauss-Jordan
# Hilton M B Ramalho - PPGE/UFPB
# Adaptado de Jhon Fox
# Introduction to the R Statistical Computing Environment
#-------------------------------------------------------------

# Definir função
reducaoGJ <- function(A,x=NULL) {
  
  
  # computar número de linhas e colunas da matriz A
  n <- nrow(A)
  m <- ncol(A)
  d <- ncol(A[,-n])
  
  # definir linha e coluna inicial
  i <- j <- 1
  
  # looping sobre linhas e colunas de A
  while (i <= n && j <= m){
    
    # iniciando looping nas colunas  
    while (j <= m){
      
      Colatual <- A[,j]
      Colatual[1:n < i] <- 0
      
      # na coluna atual achar a posição do pivô máximo
      pospivo <- which.max(abs(Colatual))
      pivo <- Colatual[pospivo]
      
      # checar se o pivô é zero
      # caso verdadeiro, passar para a próxima coluna
      if (pivo == 0) { 
        j <- j + 1
        next
      }
      
      # se a linha do pivô for maior que a linha atual - trocar as linhas
      if (pospivo > i) A[c(i, pospivo),] <- A[c(pospivo, i),] 
      
      # multiplicar a linha atual pelo inverso do pivô - primeiro elemento da linha = 1
      A[i,] <- A[i,]/pivo 
      linha <- A[i,]
      A <- A - outer(A[,j], linha) # fazer combinação usando matriz produto da linha pivo e coluna j
      A[i,] <- linha  # restaurar linha atual
      j <- j + 1
      break
    }
    i <- i + 1
  }
  # colocar linhas com zeros abaixo
  zeros <- which(apply(A[,1:m], 1, function(x) all(x == 0)))
  if (length(zeros) > 0){
    zeroLinha <- A[zeros,]
    A <- A[-zeros,]
    A <- rbind(A, zeroLinha)
    rownames(A) <- NULL
  }
  
  # Apenas para o caso de sistema nxn 
  if(!is.null(x) & n==d) {
    cat("\n Solução do sistema:\n")
    for(i in 1:n) cat(x[i],"=",A[i,m],"\n")
  }
  
  cat("\n Matriz reduzida - Gauss-Jordan\n")
  return(A)
}

#----------------------------------------------------------




# .: Exemplo 1
# Considere o sistema linear
#  m = 3, n = 3
# -2x1 + 4x2 + 6x3 = 2
#  x1 -   x2       = 0
# -x1 +   x2 +  x3 = 2
# Definir vetores Colunas -
# coeficientes do sistema linear
c1 <- c(-2,4,6)
c2 <- c(1,-1,0)
c3 <- c(-1,1,1)
b <- c(2,0,2)
# Definir matriz de variáveis 
# do sistema
x <- c("x1","x2","x3")

# Definir matriz A(3,4). 
# Matriz ampliada do sistema (A|b)
A <- cbind(c1,c2,c3,b)

# Matriz A
A

##      c1 c2 c3 b
## [1,] -2  1 -1 2
## [2,]  4 -1  1 0
## [3,]  6  0  1 2

# Redução GJ com solução do sistema
# Observar função reducaoGJ na rotina 
# Módulo 1 - Redução Gauss-Jordan.R
reducaoGJ(A,x)

## 
##  Solução do sistema:
## x1 = 1 
## x2 = 0 
## x3 = -4 
## 
##  Matriz reduzida - Gauss-Jordan

##      c1 c2 c3  b
## [1,]  1  0  0  1
## [2,]  0  1  0  0
## [3,]  0  0  1 -4

# Apenas redução GJ
reducaoGJ(A)

## 
##  Matriz reduzida - Gauss-Jordan

##      c1 c2 c3  b
## [1,]  1  0  0  1
## [2,]  0  1  0  0
## [3,]  0  0  1 -4