em , ,

Script em R para Gerar Sistemas Lineares e Resolução por Triangularização de Gauss

Nos artigos anteriores (Parte 1 e Parte 2), abordamos inicialmente um sistema gerador de sistemas lineares de qualquer dimensão, seguido do método de resolução por triangularização de Gauss. Neste artigo, vamos transformar os dois scripts apresentados em funções: gera_lineares() e tri_gauss().

Como mencionado previamente no artigo “Gerador de sistemas lineares com Script em R”, pode-se montar um sistema de equações lineares de qualquer dimensão (n). A grande vantagem desse sistema em particular é que se sabe o resultado de antemão. Assim, se a dimensão for n, o vetor solução é x’=[1,2,…,n].

Outra função elaborada a partir do Script GK1 se refere ao método de triangularização de Gauss para resolução de sistemas lineares. O Script GK3 usa as duas funções e mostra um exemplo de aplicação com n=5. O Leitor poderá testar o script com n que desejar. Na primeira parte do Script GK3, apresenta-se a função gera_lineares().

gera_lineares=function(n){
a=c(rep(0,n*n))
for (i in 1:n){
  j=i
  kb=j+(i-1)*n
  a[kb]=2}
for (i in 1:(n-1)){
  j=i+1
  kb=j+(i-1)*n
  a[kb]=-1}
for (i in 2:n){
  j=i-1
  kb=j+(i-1)*n
  a[kb]=-1}
a[n*n]=1
A=matrix(c(a),nrow=n,ncol=n,byrow=TRUE)
b=matrix(c(rep(0,n)),nrow=n,ncol=1)
b[n]=1
#https://stackoverflow.com/questions/8936099/
#returning-multiple-objects-in-an-r-function
setClass(Class="Res",representation(
A="matrix",b="matrix"))
return(new("Res",A=A,b=b))
}

Observe-se a estrutura de uma função, como ela é definida e, sobretudo, como se pode retornar os valores calculados. Nesse caso, retorna-se dois objetos que se torna possível a partir da definição de uma classe chamada “Res”. Esta solução foi encontrada em stackoverflow.com, que é um site de perguntas e respostas, de acordo com o crédito mencionado nas Linhas 19 e 20.  É conveniente retornar os objetos nas suas formas criadas na função gera_lineares().

Sistema Linear e Triangulação de Gauss

No próximo trecho do Script GK3, adiciona-se a função tri_gauss(), que permite resolver um sistema de equações lineares a partir de um data frame com n linhas por (n+1) colunas. Lembrando que o vetor b é pendurado na matriz dos coeficientes A.

Nesse caso, retorna-se o vetor solução e, portanto, a função return() recebe apenas um objeto, ao contrário da função gera_lineares(). Como se pode verificar, trata-se do mesmo algoritmo disponibilizado como Script GK1. Para alterar a dimensão do sistema, basta modificar a chamada da função gera_lineares() na Linha 52. Os objetos gerados são transformados em um data frame denominado Ab, que por sua vez é o argumento de entrada da função tri_gauss().

O Leitor poderá copiar as Linhas 1 a 52 no editor de scripts do R e executar, fazendo a seleção de todas as Linhas (ctrl-a) e depois clicando em “Executar linha ou seleção”.

tri_Gauss=function(A){
n=nrow(A)  #no. de equacoes
n1=n+1  #no. de colunas apos adicao do vetor b
for (k in 1:(n-1)){  #percorre a diagonal principal
  for (i in (k+1):(n)){
    if (abs(A[i,k])>abs(A[k,k])){
      #faz a condensacao pivotal
      xtemp=c(rep(0,n1))
      for (j in k:n1){xtemp[j]=A[k,j]}
      for (j in k:n1){A[k,j]=A[i,j]}
      for (j in k:n1){A[i,j]=xtemp[j]}
    }
    mult=-A[i,k]/A[k,k]
    for (j in k:n1){
      A[i,j]=A[i,j]+mult*A[k,j]}
  }
}
#substituicao reversa
x=c(rep(0,n))
x[n]=A[n,n1]/A[n,n]
for (i in (n-1):1){
  soma=0
  for (j in (i+1):n){
    soma=soma+A[i,j]*x[j]}
  x[i]=(A[i,n1]-soma)/A[i,i]}
return(x)
}
res=gera_lineares(5)
A=res@A; b=res@b
Ab=data.frame(A,b)
Ab
x=tri_Gauss(Ab)
x

Os resultados obtidos da execução do Script GK3 estão apresentados a seguir.

> Ab
  X1 X2 X3 X4 X5 b
1  2 -1  0  0  0 0
2 -1  2 -1  0  0 0
3  0 -1  2 -1  0 0
4  0  0 -1  2 -1 0
5  0  0  0 -1  1 1
> x=tri_Gauss(Ab)
> x
[1] 1 2 3 4 5
>

O aprendizado da linguagem R pode ser feito por meio do uso de uma função de biblioteca, no caso a função solve(), mas também pela programação do algoritmo, como demonstrado neste artigo e nos próximos que virão. Evidentemente, os nossos artigos têm a finalidade de promover a divulgação dos algoritmos usados nas bibliotecas do R, sempre que possível.

Próximos artigos

Na sequência, os próximos artigos irão abordar a resolução de sistemas lineares pelo método interativo de Gauss Seidel.

Observação importante: Esse script foi executado sob a versão 4.0.0 do R.

Escrito por Jorge Kazuo Yamamoto

Prof. Dr. Jorge Kazuo Yamamoto, fundador da Geokrigagem, é geólogo, foi pesquisador do IPT e docente do Instituto de Geociências da USP, onde se aposentou como Professor Titular do Departamento de Geologia Sedimentar e Ambiental. Atualmente, atua como Professor Sênior do Departamento de Engenharia de Minas e de Petróleo – Escola Politécnica – USP. É responsável pela disciplina “Métodos geoestatísticos” na Pós-Graduação do IPT – Investigação do subsolo: Geotecnia e Meio Ambiente. Dedica-se ao ensino de geoestatística, com ênfase no desenvolvimento de algoritmos e pesquisa de novas aplicações, tais como: variância de interpolação, cálculo da variância global de depósitos minerais e correção do efeito de suavização da krigagem. Ultimamente, seu interesse está voltado para o ensino e divulgação da linguagem R.

Deixe uma resposta

O seu endereço de e-mail não será publicado. Campos obrigatórios são marcados com *

Método de Triangularização de Gauss – Parte 2

Método Iterativo de Gauss Seidel para Resolução de Sistemas Lineares – PARTE 1