Método de Gauss-Seidel (II)

Critérios de convergência

Vamos lá?

Em uma aula anterior vimos o método de Gauss-Seidel, uma variação inteligente do método de Jacobi. Porém não ficou claro quando ele era convergente nem por quê. Nesta aula, vamos usar o resultado da aula geral sobre convergência de métodos iterativos para tentar entender o que acontece no caso do método de Gauss-Seidel.

Para poder usar o resultado de convergência que deduzimos duas aulas atrás é preciso escrever o método de Gauss-Seidel (GS) no formato do método iterativo geral apresentado naquela aula. Vamos relembrar como é a iteração de GS usando um exemplo.

Considere o sistema linear a seguir $$ \left\{ \begin{array}{rcrcrcr} 5 x_1 & - & 2x_2 & + & x_3 & = 6 \\[2mm] x_1 & - & 4x_2 & - &2x_3 & = 3 \\[2mm] 2 x_1 & + & 2x_2 & - &7x_3 & = 8 \end{array} \right. $$ Para construir a iteração, reescreva o sistema como $$ \left\{ \begin{array}{rcrcrcrcrcr} 5 x_1 & & & & & = 6 & - & 2x_2 & - & x_3 \\[2mm] x_1 & - & 4x_2 & & & = 3 & & & - & 2x_3 \\[2mm] 2 x_1 & + & 2x_2 & - &7x_3 & = 8 \end{array} \right. $$ Os valores do vetor $x$, à direita, são tomados do iterando $x^{(k)}$, e os valores do vetor $x$ à esquerda, são definidos como as componentes do iterando $x^{(k+1)}$. Ou seja, $$ \left\{ \begin{array}{rcrcrcrcrcr} 5 x_1^{(k+1)} & & & & & = 6 & - & 2x_2^{(k)} & - & x_3^{(k)} \\[2mm] x_1^{(k+1)} & - & 4x_2^{(k+1)} & & & = 3 & & & - & 2x_3^{(k)} \\[2mm] 2 x_1^{(k+1)} & + & 2x_2^{(k+1)} & - &7x_3^{(k+1)} & = 8 \end{array} \right. $$ Se estamos no passo $k$ do método, como o vetor $x^{(k)}$ já é conhecido, fica claro que $x^{(k+1)}$ é definido como a solução de um sistema linear triangular inferior acima. Se $A=L+D+U$, com $L$ e $U$ as porções triangulares inferior e superior, respectivamente, de $A$ e $D$ a matriz diagonal contendo a diagonal de $A$, então, matricialmente, a iteração de GS se escreve como \begin{equation}\label{sl2} (L+D)x^{(k+1)} = b - Ux^{(k)}. \end{equation} Desta forma, conseguimos encaixar o método de Gauss-Seidel no formato geral que havíamos visto antes, quer seja \begin{equation}\label{gs} x^{(k+1)} = -(L+D)^{-1}Ux^{(k)} + (L+D)^{-1}b = B^{GS}x^{(k)} + c \equiv \phi(x^{(k)}). \end{equation}

O fato de $x^{(k+1)}$ ser dado por \eqref{gs} não significa que a matriz $B^{GS}$ precise ser computada na prática. O melhor é determinar $x^{(k+1)}$ diretamente pela solução do sistema linear \eqref{sl2} que, por ser triangular inferior, é resolvido por sustituição.

Lembrando do resultado de convergência para métodos iterativos estacionários para sistemas lineares, para garantir a convergência do método de GS, basta pedir que $\|B^{GS}\| \lt 1$, para alguma normal matricial.

Para o método de Gauss-Jacobi foi simples relacionar o valor de $\|B^{J}\|_\infty$ diretamente com as entradas da matriz $A$ do sistema linear original. Já no caso do método de GS essa relação é possível, mas bem menos direta. Dessa análise surge o critério de Sassenfeld.

Critério de Sassenfeld: Seja $A$ é uma matriz quadrada de ordem $n$, com os elementos da diagonal todos não nulos. Defina \begin{align*} \beta_1 & = {1\over |a_{11}|} \sum_{j=2}^n |a_{1j}|,\\[2mm] \beta_i & = {1\over |a_{ii}|} \left[\sum_{j=1}^{i-1} \beta_j |a_{1j}| + \sum_{j=i+1}^n |a_{1j}|\right],\quad j =2,\ldots, n. \end{align*} Se $\max_i \{b_i\} \lt 1$, então o método de Gauss-Seidel produz uma sequência convergente para a solução do sistema linear $Ax=b$, seja qual for o vetor $b$ e o iterando inicial $x^{(0)}.$

Este critério sai de analisar a norma infinito de $B^{GS}$ e impor uma condição que garanta que seu valor seja menor que $1$, que sabemos ser suficiente para a convergência, por força do teorema mais geral demonstrado em uma aula passada.

Um corolário deste resultado é que sistemas lineares com matrizes diagonalmente dominantes por linha também têm a convergência assegurada para as sequências geradas pelo método de Gauss-Seidel.

Vamos considerar um exemplo no Octave.

> A = [5 2 -1; 4 6 3; 0 -2 4];
> b = [-10; 4; 6];

A matriz $A$ definida acima não é diagonalmente dominante por linhas (por quê?). Mas será que satisfaz o critério de Sassenfeld? Os coeficientes $\beta_i$ do critério podem ser computados com:

> n = rows(A);
> beta(1,1) = sum(abs(A(1,2:n))) / abs(A(1,1));
> for i=2:n;
    beta(i,1) = (abs(A(i,1:i-1))*beta(1:i-1) + sum(abs(A(i,i+1:n))))/abs(A(i,i));
  endfor
> beta
beta =

   0.60000
   0.90000
   0.45000

Como o valor máximo do vetor $\beta$ é menor que $1$, temos que a matriz $A$ satisfaz o critério de Sassenfeld.

> D = diag(diag(A))    # Diagonal de A
D =

Diagonal Matrix

   5   0   0
   0   6   0
   0   0   4

> L = tril(A)-D;       # Porção triangular inferior de A
> U = triu(A)-D;       # Porção triangular superior de A
> BJ = -inv(D)*(L+U);  # Matriz de iteração do método de Jacobi
> BGS = -inv(L+D)*U;   # Matriz de iteração do método de Gauss-Seidel

Acima, construímos as matrizes de iteração do método de Jacobi e do método de Gauss-Seidel. Estamos fazendo isto apenas para fins didáticos. Como comentamos acima, a construção explícita dessas matrizes não é necessária para a aplicação desses métodos.

> norm(BJ,inf)
ans =  1.1667
> norm(BGS,inf)
ans =  0.90000

No bloco acima, vemos que $\|B_{J}\|_\infty > 1$ e $\|B_{GS}\|_\infty <1$. Logo, não é possível garantir a convergência do método de Jacobi, mas sim temos garantida a convergência do método de Gauss-Seidel. Vamos iterar o método de Gauss-Seidel partindo do vetor nulo.

> x = [0;0;0];
> x1 = inv(D+L) * (b-U*x);    # 1ª iteração de Gauss-Seidel
> norm(x1-x,inf)              # norma do passo
ans =  2.5000
>
> x = x1; x1 = inv(D+L) * (b-U*x);  # 2ª iteração de Gauss-Seidel
> norm(x1-x,inf)                    # norma do passo
ans =  1.0500
>
> x = x1; x1 = inv(D+L) * (b-U*x;)  # 3ª iteração de Gauss-Seidel
> norm(x1-x,inf)                    # norma do passo
ans =  0.31500
>
> x = x1; x1 = inv(D+L) * (b-U*x);  # 4ª iteração de Gauss-Seidel
> norm(x1-x,inf)                    # norma do passo
ans =  0.015750
>
> x1                                # Aproximação corrente
x1 =

  -2.00075
   0.99988
   1.99994

Neste exemplo, após o cálculo de cada novo iterando foi sendo computado $\|x^{(k+1)} - x^{(k)}\|_\infty,$ isto é a norma entre dois iterandos consecutivos, conhecido como passo. Por vim, exibimos o último iterando. A saber, a solução exata deste sistema linear é $x = (-2,1,2)^T.$

1. Escreva um algoritmo para o método de Gauss-Seidel e conte o número de operações realizadas em cada iteração do método (não esqueça de computar as operações necessárias para a verificação do critério de parada). Implemente seu algoritmo em Octave.

2. Prove que se uma matriz é diagonalmente dominante por linhas, então ela também satisfaz o critério de Sassenfeld.