Sistemas Lineares e métodos iterativos

Prévia do material em texto

Cálculo Numérico – Profa. Simone Tomazzoni 
 2014 
 
Página | 1 
 
Sistemas Lineares – Métodos Iterativos 
 
 Os métodos diretos para resolver um SELAS, como a eliminação de Gauss, a fatoração LU e outros 
tornam-se proibitivos em termos de tempo de computador e armazenamento se a matriz A, dos coeficientes, 
é bastante grande. Há situações práticas de engenharia, tais como a discretização de equações diferenciais 
parciais, onde a ordem da matriz pode chegar facilmente a várias centenas de milhares. Exemplos disso são 
a computação gráfica e a termodinâmica. Além disso, para maioria dos problemas de grande porte a matriz 
dos coeficientes é uma matriz esparsa, o que significa que tem muitos zeros, e a esparsidade perde-se nos 
procedimentos de triangulação. Para tais problemas é aconselhável usar uma classe de métodos, chamados 
métodos iterativos, que nunca alteram a esparsidade da matriz e que requerem armazenamento apenas de 
alguns vetores de n componentes por vez. 
 A ideia básica por detrás dos métodos iterativos vê de um teorema de ponto fixo, que, traduzido para 
um SELAS Ax = b, significa que, escrevendo o mesmo de uma maneira equivalente, 
x = Bx + d 
e iniciando com uma solução aproximada x
(0)
, gera-se uma sequência de aproximações 
 ( )kx
, definida 
iterativamente por 
(1) (0)
(2) (1)
( 1) ( ) , 1,2,3...k k k
 
 
  
x Bx d
x Bx d
x Bx d
 
que, sob certas condições, converge para a solução exata, ignorados os arredondamentos. 
 Para resolver um SELAS iterativamente, então precisamos saber como escrevê-lo na forma acima e 
como escolher x
(0)
 para que a sequência convirja, ou sob que condições o processo iterativo converge, 
qualquer que seja a escolha de x
(0)
. 
 
Método de Jacobi
1
 
 
Um SELAS Ax = b, ou 
11 1 12 2 1 1
21 1 22 2 2 2
1 1 2 2
n n
n n
n n nn n n
a x a x a x b
a x a x a x b
a x a x a x b
   
    


    
 
 
1 Carl Gustav Jacobi apresentou o método iterativo que recebe seu nome em 1845. 
Cálculo Numérico – Profa. Simone Tomazzoni 
 2014 
 
Página | 2 
 
 
pode ser escrito, supondo aii  0 para i =1, 2, ..., n, 
1 1 12 2 1
11
2 2 21 1 2
22
1 1 ( 1) 1
1
( ... )
1
( ... )
1
( ... )
n n
n n
n n n n n n
nn
x b a x a x
a
x b a x a x
a
x b a x a x
a
 

   


   




   

 
Matricialmente, 
112
11 11 1
23 221 11
1 1
22 22 22 2
2 2
22
( 1)
( 1)( 1)
( 1)1 2
0
0
0
n
n
n n
n n
n n n
nnn nn n
nn nn nn
aa
a a b
a aa a
x x
a a a b
x x
a
a
x x
a b
aaa a
a a a

 

 
  
 
 
      
    
     
    
        
 
 
   
 
x x
B
 
 
 
 
 
 
 
 
 
  
d
 
o que está na forma x = Bx + d. 
Exemplo 1. (manual) Resolver o SELAS dado por 
1 2 3
1 2 3
1 2 3
10 2 7
5 8
2 3 10 6
x x x
x x x
x x x
  

   
   
 empregando 4 
iterações do método de Jacobi, com x
(0)
 = [0.7, -1.6, 0.6]. 
1 2 3
2 1 3
3 1 2
1
(7 2 )
10
1
( 8 )
5
1
(6 2 3 )
10
x x x
x x x
x x x

  


   


  

 
 
 
Cálculo Numérico – Profa. Simone Tomazzoni 
 2014 
 
Página | 3 
 
 
k 0 1 2 3 4 
x1 0,7 
x2 -1,6 
x3 0,6 
 
Critério de parada 
 Como já vimos anteriormente, a solução exata quase nunca é atingida num processo iterativo. Assim, 
necessitamos ter algum critério para terminar o processo infinito, adotando como solução aproximada o 
último termo calculado da sequência iterativa. 
 Quando há convergência, x
(k+1)
 é uma aproximação melhor da solução x do que x
(k)
. Então um 
critério natural de parada é interromper as iterações se (erro relativo) 
( 1) ( )
( )
k k
k

 

x x
x 
para um ε positivo preestabelecido chamado comumente de tolerância, e escolhido de acordo com a exatidão 
desejada, , e alguma norma de vetores 
.
. Outro critério é se o número de iterações excede o número 
preestabelecido. Teoricamente, na expressão acima, em lugar de x
(k)
 deveríamos colocar a solução exata x, 
para produzir o erro relativo, mas isso não é possível uma vez que x não é conhecido. 
 
Exemplos de normas para o R
n 
Dado um vetor x  Rn, de componentes xi, definimos as normas: 
 
1 21
1
2 2 2 2 2 2 22 22
1 2 1 22
1
1 2
Norma 1: : ...
Norma 2: : ... ...
Norma : : max ...
n
i n
i
n
i n n
i
i n
x x x x
x x x x x x x
x x x



    
        
    


x
x
x
 
 
Cálculo Numérico – Profa. Simone Tomazzoni 
 2014 
 
Página | 4 
 
Exemplo 2. Calcule o erro relativo do resultado do exemplo 1, considerando as normas 1, 2 e . 
 
Exercício 
 Use iterações de Jacobi para resolver o sistema linear 
1 2 3
1 2 3
1 2 3
10 2 27
3 6 2 61,5
5 21,5
x x x
x x x
x x x
  

    
    
 até que o 
erro relativo caixa abaixo de ε = 0,05. Use estimativa inicial x(0) = [0, 0, 0]. 
 
k 0 1 2 3 4 ... n 
x1 
x2 
x3 
( 1) ( )
( )
k k
k



x x
x
 
- 
 
 Matricialmente, a obtenção da matriz de iteração de Jacobi, a matriz B na expressão do ponto fixo, é 
obtida da seguinte forma: a matriz dos coeficientes A é decomposta em três parcelas A = L + D + U, sendo 
 
21
31 32
1 2 ( 1)
0 0 0 0
0 0 0
0:
0
0n n n n
a
a a
a a a 
 
 
 
 
 
 
 
 
L
, 
11
22
33
0 0 0
0 0 0
: 0 0
0
0 0 0 nn
a
a
a
a
 
 
 
 
 
 
  
D
 e 
12 13 1
23 2
( 1)
0
0 0
: 0 0 0
0 0 0 0
n
n
n n
a a a
a a
a 
 
 
 
 
 
 
  
U
 
Assim, o sistema linear Ax = b, na forma da equação do ponto fixo, resulta 
Cálculo Numérico – Profa. Simone Tomazzoni 
 2014 
 
Página | 5 
 
-1 -1
B
Ax = b
(L + D + U)x = b
x = -D (L + U)x + D b
d
 
 Uma forma de efetuar os cálculos para o método de Jacobi no MATLAB, definir a matriz B e o vetor 
d, conforme acima, e calcular, sucessivamente, x
(k+1)
 = Bx
(k)
 + d. Há comandos especiais para obter as 
matrizes L, U e D, bem como sua inversa, conhecida a matriz A 
No MATLAB 
>>% Dada uma matriz A 
>> A=[ 2 3 6; 1 5 7; -2 9 4] 
A = 
 2 3 6 
 1 5 7 
 -2 9 4 
> L=tril(A,-1) 
L = 
 0 0 0 
 1 0 0 
 -2 9 0 
>> U=triu(A,1) 
U = 
 0 3 6 
 0 0 7 
 0 0 0 
>> diag(diag(A)) %cria uma matriz diagonal com os elementos da diagonal de A 
ans = 
 2 0 0 
 0 5 0 
Cálculo Numérico – Profa. Simone Tomazzoni 
 2014 
 
Página | 6 
 
 0 0 4 
>> D=diag(diag(A)) 
D = 
 2 0 0 
 0 5 0 
 0 0 4 
>>inv(D) 
ans = 
 1/2 0 0 
 0 1/5 0 
 0 0 1/4 
>> L+D+U %Obtém-se a matriz A 
ans = 
 2 3 6 
 1 5 7 
 -2 9 4 
Outros comandos que podem ser úteis: 
>> % Dado um vetor x 
>> x = [-1 3 7 -9] 
x = 
 -1 3 7 -9 
>> max(abs(x)) 
ans = 
 9 
>> sum(abs(x)) 
ans = 20 
 
Cálculo Numérico – Profa. Simone Tomazzoni 
 2014 
 
Página | 7 
 
Obs. No método de Jacobi todos os valores de x da iteração (k+1) dependem dos valores de x da iteração (k) 
por isso o método é também chamado de Método dos deslocamentos simultâneos. 
 
Exemplo. 
Use 6 iterações de Jacobi para resolver o sistema linear 
1 2 3
1 2 3
1 2 3
8 8
7 2 4
2 9 12
x x x
x x x
x x x
  

   
   
 . Use estimativa 
inicial x
(0)
 = [0, 0, 0]. Calcule o erro relativo na última iteração. 
k 0 1 2 3 4 5 6 
x1 
x2 
x3 
 
Método de Gauss-Seidel 
 
 No método de Jacobi, para calcular as componentes de x
(k+1)
, são usadas todas as componentes de 
x
(k)
, entretanto observemos que, para calcular xi
(k+1)
, poderíamos ter utilizado xi
(k+1)
, i = 1, 2, 3,...,i-1, já 
disponíveis. 
( 1) ( ) ( ) ( )
1 1 12 2 13 3 1
11
( 1) ( 1) ( ) ( )
2 2 21 1 23 3 2
22
( 1) ( 1) ( 1) ( )
1 1 2 2 ( 1) 1
1
( ... )
1
( ... )
1
( ... )
k k k k
n n
k k k k
n n
k k k k
n n n n n n n
nn
x b a x a x a x
a
x b a x a x a x
a
x b a x a x a x
a

 
  
 

   


    




   

 
 
 Repetindo, a ideia é usar cada componente logo que for disponível, para o cálculo da componente 
seguinte. A modificação expressa nas equações acima é conhecida como Método de Gauss-Seidel. 
 
 
Cálculo Numérico – Profa. Simone Tomazzoni 
 2014 
 
Página | 8 
 
Exemplo. 
 Use 4 iterações de Gauss-Seidel para resolver o sistema linear 
1 2 3
1 2 3
1 2 3
8 8
7 2 4
2 9 12
x x x
x x x
x x x
  

   
   
. Use 
estimativa inicial x
(0)
 = [0, 0, 0]. Calcule o erro relativo na última iteração. 
 
k 0 1 2 3 4 
x1 
x2 
x3 
 
 Convém observar que, no método de Jacobi, x
(k+1)
 é totalmente determinada usando-se as 
componentes de x
(k)
. No método de Gauss-Seidel x
(k+1)
 é determinada usando-se as componentes de x
(k)
 e de 
x
(k+1)
 já determinadas, com a vantagem de não exigir o armazenamento simultâneo dos dois vetores x
(k)
 e 
x
(k+1)
 em cada passo e em geral convergir mais rapidamente que o método de Jacobi, como se constata no 
exemplo dado e na análise de sua convergência que será feita adiante. 
 
 Matricialmente, para obtenção da matriz de iteração Gauss-Seidel, a matriz B na expressão do ponto 
fixo, é obtida da seguinte forma: a matriz dos coeficientes A é decomposta em três parcelas A = L + D + U, 
sendo 
 
21
31 32
1 2 ( 1)
0 0 0 0
0 0 0
0:
0
0n n n n
a
a a
a a a 
 
 
 
 
 
 
 
 
L
, 
11
22
33
0 0 0
0 0 0
: 0 0
0
0 0 0 nn
a
a
a
a
 
 
 
 
 
 
  
D
 e 
12 13 1
23 2
( 1)
0
0 0
: 0 0 0
0 0 0 0
n
n
n n
a a a
a a
a 
 
 
 
 
 
 
  
U
 
Assim, o sistema linear Ax = b, na forma da equação do ponto fixo, resulta 
Cálculo Numérico – Profa. Simone Tomazzoni 
 2014 
 
Página | 9 
 
( 1) ( )( ) ( )k k  -1 -1
B
Ax = b
(L + D + U)x = b
x = - D L U x + D L b
d
 
 Uma forma de efetuar os cálculos para o método de Gauss-Seidel no MATLAB, definir a matriz B e 
o vetor d, conforme acima, e calcular, sucessivamente, x
(k+1)
 = Bx
(k)
 + d. 
 
Exemplo. 
Use 7 iterações de Gauss-Seidel para resolver o sistema linear 
1 2 3
1 2 3
1 2 3
5 7
5 7
5 7
x x x
x x x
x x x
  

  
   
. Use 
estimativa inicial x
(0)
 = [0, 0, 0]. Calcule o erro relativo na última iteração. 
 
k 0 1 2 3 4 5 6 7 
x1 
x2 
x3 
 
Convergência 
 É muitas vezes difícil acertar uma boa estimativa inicial x
(0)
. Assim, seria bom ter condições que 
garantam a convergência dos métodos de Gauss-Seidel e Jacobi, seja qual for a escolha da aproximação 
inicial. 
Teorema. Se A é uma matriz estritamente diagonal dominante, então os processos de Jacobi e de Gauss-
Seidel convergem, independente da estimativa inicial x
(0)
. Ou seja, 
1
, para 1,2,...,
n
ii ij
j
j i
a a i n


 
 
 
Cálculo Numérico – Profa. Simone Tomazzoni 
 2014 
 
Página | 10 
 
 A expressão acima indica que o valor absoluto do coeficiente da diagonal em cada equação deve ser 
maior do que a soma dos valores absolutos dos coeficientes na equação. 
 Este critério é suficiente, mas não necessário para a convergência, ou seja, embora o método possa às 
vezes funcionar se o mesmo não for satisfeito, a convergência é garantida se o mesmo for satisfeito. 
 Felizmente, muitos problemas de importância prática na engenharia e na ciência satisfazem tal 
exigência. 
 Uma consequência de um conhecido teorema de Stein-Rosenberg é que, para uma matriz A com 
elementos diagonais positivos e os demais negativos, 
 Os métodos de Jacobi e Gauss-Seidel ambos convergem ou ambos divergem; 
 Quando há convergência, o método de Gauss-Seidel converge mais rapidamente. 
 
 Para o caso geral, infelizmente, não é possível fazer tais afirmações sobre convergência. Mas quando 
ambos os métodos convergem, o de Gauss-Seidel costuma ser preferido devido ao seu menor volume de 
armazenamento exigido na memória do computador. Atualmente, devido à possibilidade de processamento 
paralelo, o método de Jacobi vem sendo amplamente empregado devido à sua estrutura que permite a 
paralelização. 
 
 Um teorema mais geral acerca da convergência dos métodos iterativos de ponto fixo, como Gauss-
Seidel e Jacobi é o seguinte: 
 
Teorema. A sequência definida por x
(k+1)
 = Bx
(k)
 + d converge, seja qual for a estimativa inicial x
(0)
, se e 
somente se B
k 0. 
 
 É possível mostrar que 
( ) 1 1k     B 0 B B
, para alguma norma matricial
2
. 
 
Na expressão acima, (B) é o raio espectral da matriz B. 
 
Comparação dos Métodos de Solução de Sistemas Lineares 
 
Métodos Diretos 
 
1. Processos finitos (convergência para qualquer sistema não singular); 
2. Apresentam problemas com erros de arredondamento; 
3. Utilizam-se técnicas de pivoteamento para amenizar os problemas de arredondamento; 
4. O processo de triangularização destrói a esparsidade da matriz de coeficientes. 
5. Para matrizes cheias a solução requer 
3n
operações sem considerar o pivoteamento. 
 
 
2 Um exemplo de norma matricial é a norma euclidiana 
2
2
1 1
n n
ijE
j i
a
 
 A
 
 
Cálculo Numérico – Profa. Simone Tomazzoni 
 2014 
 
Página | 11 
 
Métodos Iterativos 
 
1. Provavelmente mais eficientes para sistemas de grande porte, principalmente com a utilização de 
computação de alto desempenho (vetorial e paralela); 
2. Tem convergência assegurada apenas sob certas condições; 
3. Conserva a esparsidade da matriz de coeficientes; 
4. Os métodos de Jacobi e Gauss-Seidel 
22n
operações por iteração.Exercícios. 
 
1. O sistema de equações a seguir é projetado para determinar as concentrações c1, c2 e c3 (em g/m
3
) em 
uma série de reatores acoplados em função de uma quantidade de entrada de massa para cada reator ( 
lado direito das equações em g/dia). Resolva esse problema com o método de Gauss-Seidel com 
tolerância ε = 5% para o erro relativo. 
 
1 2 3
1 2 3
1 2 3
15 3 3800
3 18 6 1200
4 12 2350
c c c
c c c
c c c
  

   
   
 
 
2. Dos três conjuntos de equações lineares a seguir, identifique o(s) conjunto(s) que você não poderia 
resolver utilizando Jacobi ou Gauss-Seidel. Mostre, utilizando qualquer número de iterações, que 
necessariamente sua solução não converge. Estabeleça, de modo claro, seu critério de convergência 
(como você sabe que a solução não está convergindo). 
 
 
Conjunto 1 Conjunto 2 Conjunto 3 
8x + 3y + z = 12 
 
x + y + 5z = 7 
 
– x +3 y + 5z = 7 
-6x + 7z = 1 
 
x + 4y – z = 4 
 
–2 x +4 y – 5z = –3 
2x + 4y – z = 5 
 
3x + y – z = 3 
 
2y –z = 1