2. Sistemas de equações lineares

Jaime E. Villate. Métodos Numéricos com Maxima,
Universidade do Porto, Portugal, 2015

  1. Sistemas de equações lineares
    1. Representação matricial
    2. Método de eliminação de Gauss
    3. Método de Gauss-Jordan

2.1. Representação matricial

Um sistema linear de ordem n é um sistema com n variáveis x 1 , x 2 , …, x n , que satisfazem n equações lineares:

(2.1)
a 11 x 1 + a 12 x 2 +· · · + a 1 n x n = b 1 a 21 x 1 + a 21 x 2 +· · · + a 2 n x n = b 2 .. . a n 1 x 1 + a n 2 x 2 +· · · + a n n x n = b n

Se as equações são linearmente independentes, existe sempre uma única solução. Claro está que basta dar a lista de equações e variáveis aos comandos solve ou linsolve do Maxima para obter a solução, mas o objetivo deste capítulo não é obter a resposta mas sim mostrar, passo a passo, como funcionam alguns métodos de resolução.

Chama-se matriz aumentada do sistema à matriz:

(2.2)
a 11 a 12 ·· · a 1 n b 1 a 21 a 22 ·· · a 2 n b 2 .. . . . . . . . . . . .. . a n 1 a n 2 ·· · a n n b n

2.2. Método de eliminação de Gauss

A qualquer equação no sistema pode somar-se outra das equações, multiplicada por uma constante k , e a solução do sistema não se altera. Na matriz aumentada, essa operação corresponde a substituir uma linha L i por L i + k L j . Essa propriedade pode ser usada para obter uma matriz triangular superior, em que a i j é zero, se i < j :

(2.3)
a 11 a 12 ·· · a 1 n b 1 0 a 22 ·· · a 2 n b 2 .. . . . . . . . . . . .. . 0 0 ·· · a n n b n

que corresponde a outro sistema de equações com a mesma solução do sistema original (as constantes a i j e b j não são as mesmas a i j e b j da matriz inicial).

A última equação nesse sistema permite obter facilmente o valor da variável x n

(2.4)
a n n x n = b n x n = b n a n n

Este resultado pode ser substituído na penúltima equação para obter o valor de x n 1 ,

(2.5)
a n 1 , n 1 x n 1 + a n 1 , n x n = b n 1 x n = b n 1 a n 1 , n 1 a n 1 , n b n a n 1 , n 1 a n n

e assim sucessivamente, até se obter o valor de todas as variáveis.

Exemplo 2.1

Resolva o sistema de equações:

8 x + 4 . 5 z 12 . 4 = 3 . 2 y 2 x 0 . 8 y = 7 . 6 + 6 . 1 z x + 3 y + 10 . 2 z = 8 . 4

usando o método de eliminação de Gauss.

Resolução. No Maxima, convém guardar as 3 equações em 3 variáveis

(%i1) eq1: 8*x + 4.5*z - 12.4 = 3.2*y$
(%i2) eq2: 2*x - 0.8*y = 7.6 + 6.1*z$
(%i3) eq3: x + 3*y + 10.2*z = 8.4$

Este sistema pode ser escrito da mesma forma do sistema geral 2.1 para obter a matriz aumentada. A função augcoefmatrix realiza esse procedimento dando como resultado a matriz aumentada. É necessário dar dois argumentos a augcoefmatrix, a lista das equações e a lista das variáveis, na ordem que se quiser. Neste caso, a matriz aumentada pode ser guardada na variável m com o seguinte comando

(%i4) m: float (augcoefmatrix ([eq1, eq2, eq3], [x, y, z]));
(%o4)    8 . 0 3 . 24 . 5 12 . 42 . 0 0 . 8 6 . 1 7 . 61 . 03 . 010 . 2 8 . 4

optou-se por usar a função float para evitar que a matriz seja armazenada em números racionais.

Há que ter em conta que os números na última coluna não são b i mas sim b i (ou seja, no sistema 2.1 todas as equações foram igualadas a zero).

Antes de começar com o processo de eliminação pode fixar-se um número reduzido de algarismos significativos, por exemplo 4, para os resultados que serão apresentados a seguir, de modo a evitar que as matrizes tenham muitos algarismos:

(%i5) fpprintprec: 4$

Como está a ser usado o formato de vírgula flutuante de dupla precisão (8 bytes), internamente os cálculos continuarão a ser feitos com 16 algarismos significativos, mas todos os resultados apresentados no ecrã serão arredondados para 4 algarismos significativos.

Para reduzir a matriz à forma triangular, começa-se por substituir a segunda linha, L 2 , por L 2 - (1/4)  L 1 , onde L 1 representa toda a primeira linha da matriz. O objetivo é fazer com que a 21 fique igual a zero (note-se que 1/4 foi calculado dividindo a 21 por a 11 ).

No Maxima, a função rowop( m , i , j , k ) produz uma nova matriz, em que a linha L i da matriz m é substituída por L i - k L j e as restantes linhas ficam iguais. Assim sendo, a substituição L 2 - (1/4)  L 1 L 2 pode ser feita com o seguinte comando:

(%i6) m: rowop (m, 2, 1, 1/4);
(%o6)    8 . 0 3 . 24 . 5 12 . 40 . 00 . 0 7 . 225 4 . 51 . 03 . 010 . 2 8 . 4

A seguir, para eliminar a 31 = 1.0 usando novamente a 11 = 8.0, faz-se a substituição L 3 - k L 1 L 3 , onde k = a 11 / a 31 = 1/8. Ou seja,

(%i7) m: rowop (m, 3, 1, 1/8);
(%o7)    8 . 0 3 . 24 . 5 12 . 40 . 00 . 0 7 . 225 4 . 50 . 03 . 49 . 637 6 . 85

Neste momento deveria usar-se a 22 para eliminar a 23 . No entanto, isso não é possível porque a 22 é nula ( k = a 23 / a 22 não existe). Em casos como este, troca-se a linha em questão (neste caso a segunda) com uma linha posterior (neste caso apenas pode ser com a terceira). Essa troca de linhas equivale a alterar a ordem das equações, sem que por isso se altere a solução. Usa-se a função rowswap para trocar a ordem de duas linhas:

(%i8) m: rowswap (m, 2, 3);
(%o8)    8 . 0 3 . 24 . 5 12 . 40 . 03 . 49 . 637 6 . 850 . 00 . 0 7 . 225 4 . 5

ficando assim uma matriz diagonal superior, que permite passar à segunda fase do método.

Para resolver as equações associadas a cada linha da matriz é conveniente definir uma lista com as variáveis, na mesma ordem em que foram usadas para encontrar a matriz aumentada, seguidas por 1:

(%i9) vars: [x, y, z, 1]$

Multiplicando a terceira linha da matriz m pela lista de variáveis, obtém-se a última equação, na qual só aparece a variável z ; resolve-se essa equação para encontrar o valor de z , que será guardado para ser usado nos passos seguintes

(%i10) sol: float (solve (m[3].vars));
(%o10)    [ z = -0.6228 ]

Multiplicando a segunda linha da matriz pela lista de variáveis e substituindo o valor que já foi encontrado para z , obtém-se uma equação que depende unicamente de y . A resolução dessa equação conduz ao valor de y e é conveniente acrescentar esse resultado na mesma lista em que foi armazenado o valor de z

(%i11) sol: append (float (solve ( subst (sol, m[2].vars))), sol);
(%o11)    [ y = 3.78, z = -0.6228 ]

Finalmente, repete-se o comando anterior, mas agora multiplicando pela primeira linha da matriz, para determinar o valor de x

(%i12) sol: append (float (solve ( subst (sol, m[1].vars))), sol);
(%o12)    [ x = 3.412, y = 3.78, z = -0.6228 ]

2.3. Método de Gauss-Jordan

O método de eliminação de Gauss tem a vantagem de não alterar o valor do determinante da matriz com coeficientes a i j e, como tal, o determinante pode ser calculado facilmente multiplicando os valores a i i na diagonal da matriz triangular obtida.

O método de Gauss-Jordan consiste em transformar a matriz aumentada até os coeficientes a i j corresponder à matriz identidade, ou seja, iguais a 1 se i = j ou zero caso contrário. A grande vantagem é que os valores finais de b i são os valores das variáveis, sem ser preciso realizar mais contas. Em contrapartida, o determinante da matriz já não é o produto dos elementos na diagonal da matriz reduzida.

Para conseguir transformar a matriz na matriz identidade, começa-se por dividir a primeira equação por a 11 , para que a 11 fique igual a 1, e a seguir usa-se essa primeira linha para eliminar o primeiro elemento em todas as outras linhas. A seguir divide-se a segunda linha por a 22 , para que este fique igual a 1 e usa-se essa segunda linha para eliminar o segundo elemento em todas as outras linhas, incluindo a primeira.

Não existe um comando específico do Maxima para dividir uma linha de uma matriz por uma constante; no entanto, a divisão por k pode ser feita com a função rowop( m , i , i , 1 1/ k ), que substitui a linha L i da matriz m , por L i - (1-1/ k )  L i = L i / k .

Exemplo 2.2

Determine a solução do seguinte sistema, pelo método de Gauss-Jordan

2 x 1 1 . 2 x 2 4 x 3 = 5 3 x 1 + 3 x 2 + 5 x 3 = 122 . 5 x 1 + 4 x 2 3 x 3 = 6

Resolução. Para mostrar como é possível obter soluções exatas, na resolução deste exemplo só será usada a função float no fim e os dois números de vírgula flutuante 1.2 e 2.5 na primeira e terceira equações serão escritos como as fracções 12/10 e 25/10.

Em vez de se usar a função augcoefmatrix, a matriz aumentada pode também ser inserida diretamente usando a função matrix

(%i13) m: matrix ([2,-12/10,-4,5],[-3,3,5,12],[25/10,4,-3,-6]);
(%o13)    2 6 5 45 33 5 12 5 24 3 6

Começa-se por dividir a primeira linha pelo seu primeiro elemento (ou seja 2):

(%i14) m: rowop (m, 1, 1, 1 - 1/m[1][1]);
(%o14)    1 3 5 25 2 33 5 12 5 24 3 6

A seguir, subtrai-se a primeira linha, multiplicada pelo primeiro elemento da segunda linha, para eliminar o primeiro elemento de segunda linha:

(%i15) m: rowop (m, 2, 1, m[2][1]);
(%o15)    1 3 5 25 20 6 5 139 25 24 3 6

e faz-se o mesmo com a primeira e a terceira linhas, para eliminar o primeiro elemento da terceira linha:

(%i16) m: rowop (m, 3, 1, m[3][1]);
(%o16)    1 3 5 25 20 6 5 139 20 11 22 49 4

Repete-se o mesmo procedimento para a segunda coluna; ou seja, divide-se a segunda linha pelo seu segundo elemento,

(%i17) m: rowop (m, 2, 2, 1 - 1/m[2][2]);
(%o17)    1 3 5 25 20 1 5 665 40 11 22 49 4

e usa-se essa nova linha para eliminar o segundo elemento na primeira e na terceira linhas:

(%i18) m: rowop (rowop (m, 1, 2, m[1][2]), 3, 2, m[3][2]);
(%o18)    10 5 249 40 1 5 665 40 0 79 12 813 8

Finalmente, repete-se o mesmo procedimento para a terceira coluna:

(%i19) m: rowop (m, 3, 3, 1 - 1/m[3][3]);
(%o19)    10 5 249 40 1 5 665 40 0 1 2439 158
(%i20) m: rowop (rowop (m, 1, 3, m[1][3]), 2, 3, m[2][3]);
(%o20)    10 0 2081 790 1 0 535 1580 0 1 2439 158

Os três números na última coluna são os valores exatos das 3 variáveis. Para aproximar esses números por números de vírgula flutuante e colocá-los numa lista, extrai-se a quarta linha da matriz transposta e aplica-se a função float

(%i21) float (transpose(m)[4]);
(%o21)    [ -26.34, 3.386, -15.44 ]

Note-se que se tivesse sido usada a função augcoefmatrix para construir a matriz, seria necessário mudar os sinais desses 3 valores obtidos.

Este algoritmo pode ser automatizado facilmente. O mesmo resultado obtido com os comandos (%i13) até (%i20) podia ter-se obtido com o seguinte ciclo:

(%i22) m: matrix ([2,-12/10,-4,5],[-3,3,5,12],[25/10,4,-3,-6]);
(%i23) for i:1 thru 3 do (
    m: rowop (m, i, i, 1 - 1/m[ i ][ i ] ),
    for j:1 thru 3 do (
        if (i # j) then m: rowop (m, j, i, m[ j ][ i ] )));
(%o23)    done

onde # é o operador "diferente de".

Se em alguma das iterações m[ i ][ i ] for nulo, 1/m[ i ][ i ] produz um erro. Para evitar esse tipo de erro seria necessário acrescentar instruções para conferir que m[ i ][ i ] não seja nulo e, caso contrário, trocar a linha L i com alguma das linhas L j com j > i .