1.2 Métodos Iterativos

Em revisão

A ideia básica é resolver o sistema linear

por um método de projeção. Mais especificamente, busca-se uma solução aproximada $x_m\in {ℝ}^n$ no subespaço afim $x_0+{𝒦}_m$ de dimensão $m$, impondo-se a condição de Petrov¹³¹³endnote: ¹³Georgi Iwanowitsch Petrov, 1912 - 1987, engenheiro soviético. Fonte: Wikipedia.-Galerkin¹⁴¹⁴endnote: ¹⁴Boris Galerkin, 1871 - 1945, engenheiro e matemático soviético. Fonte: Wikipédia.

onde ${ℒ}_m$ também é um subespaço de dimensão $m$. Quando ${𝒦}_m$ é um subespaço de Krylov, i.e.

temos o método de subespaço de Krylov. Aqui, temos o resíduo

sendo $x_0$ uma aproximação inicial para a solução do sistema. Notamos que com isso, temos que a aproximação calculada é tal que

onde $q_{m-1}$ é um dado polinômio de grau $m-1$. No caso particular de $x_0=0$, temos

Diferentes versões deste método são obtidas pelas escolhas do subespaço ${ℒ}_m$ e formas de precondicionamento do sistema.

Em revisão

O GMRES é um método de subespaço de Krylov assumindo ${ℒ}_m=A{𝒦}_m$, com

onde $v_1=r_0/\Vert r_0\Vert$ é o vetor unitário do resíduo $r_0=b-A{x}_0$ para uma dada aproximação inicial $x_0$ da solução do sistema $Ax=b$.

Vamos derivar o método observando que qualquer vetor $x$ em $x_0+{𝒦}_m$ pode ser escrito como segue

onde, $V_m=[v_1,\mathrm{\dots },v_m]$ é a matriz $n\times m$ cujas colunas formam uma base ortogonal $\{v_1,\mathrm{\dots },v_m\}$ de ${𝒦}_m$ e $y\in R^m$. Aqui, $V_m$ é computada usando-se o seguinte método de Arnoldi¹⁵¹⁵endnote: ¹⁵Walter Edwin Arnoldi, 1917 - 1995, engenheiro americano estadunidense. Fonte: Wikipédia.- Gram¹⁶¹⁶endnote: ¹⁶Jørgen Pedersen Gram, 1850 - 1916, matemático dinamarquês. Fonte: Wikipédia.-Schmidt¹⁷¹⁷endnote: ¹⁷Erhard Schmidt, 1876 - 1959, matemático alemão. Fonte: Wikipédia. modificado [9, Subseção 6.3]:

1. 1.

   Dado $v_1$ de norma 1

2. 2.

   Para $j=1,\mathrm{\dots },m$:

   1. (a)

      $w_j\leftarrow A{v}_j$

   2. (b)

      Para $i=1,\mathrm{\dots },j$:

      1. i.

         $h_{i,j}\leftarrow (w_j,v_i)$

      2. ii.

         $w_j\leftarrow w_j-h_{i,j}{v}_i$

   3. (c)

      $h_{j+1,j}\leftarrow \Vert w_j\Vert$

   4. (d)

      Se $h_{j+1,j}\leftarrow 0$, então pare.

   5. (e)

      $v_{j+1}\leftarrow w_j/h_{j+1,j}$

Seja, então, ${\overline{H}}_m={[h_{i,j}]}_{i,j=1}^{m+1,m}$ a matriz de Hessenberg¹⁸¹⁸endnote: ¹⁸Karl Adolf Hessenberg, 1904 - 1959, engenheiro e matemático alemão. Fonte: Wikipédia. cujas entradas não nulas são computadas pelo algoritmo acima (Passos 2(a)i-ii). Pode-se mostrar que [9, Proposição 6.5]

onde, $\beta =\Vert r_0\Vert$.

A aproximação GMRES $x_m$ é então computada como

Observamos que este último é um pequeno problema de minimização, sendo que requer a solução de um sistema $(m+1)\times m$ de mínimos quadrados, sendo $m$ normalmente pequeno.

Em resumo, a solução GMRES $x_m$ é computada seguindo os seguintes passos:

1. 1.

   Escolhemos uma aproximação inicial $x_0$ para a solução de $Ax=b$.

2. 2.

   Computamos o resíduo $r_0=b-A{x}_0$.

3. 3.

   Computamos o vetor unitário $v_1=r_0/\Vert r_0\Vert$.

4. 4.

   Usamos o método de Arnoldi-Gram-Schmidt modificado para calculamos uma base ortogonal $V_m$ de ${𝒦}_m$ e a matriz de Hessenberg ${\overline{H}}_m$ associada.

5. 5.

   Computamos $y_m={\min }_y\Vert \beta e_1-{\widehat{H}}_{m}y\Vert$.

6. 6.

   Computamos $x_m=x_0+V_{m}y$.

Em revisão

Em revisão

A ideia é reformular o sistema $Ax=b$ como um problema de minimização. Vamos começar definindo o funcional

O vetor $y$ que minimiza $J$ é a solução de $Ax=b$. De fato, denotando $x$ a solução de $Ax=b$, temos

O último termo é independente de $y$ e, portanto, $J$ é mínimo quando

é minimizado. Agora, como $A$ é definida positiva²²²²endnote: ²²$x^{T}Ax>0$ para todo $x\ne 0$., o menor valor deste termo ocorre quando $y-x=0$, i.e. $y=x$.

Observamos, também, que o gradiente de $J$ é

i.e., é o oposto do resíduo $r=b-Ay$. Com isso, temos que $y=x$ é a única escolha tal que $\nabla J=0$. Ainda, temos que $\nabla J$ é o vetor que aponta na direção e sentido de maior crescimento de $J$. Isso nos motiva a aplicarmos a seguinte iteração²³²³endnote: ²³Iteração do método do máximo declive.

onde, ${\alpha }_k>0$ é um escalar que regula o tamanho do passo a cada iteração. Lembrando que $-\nabla J=r$, temos que a iteração é equivalente a

Notamos que $x^{(k+1)}$ é um ponto na reta $\{x^{(k)}+\alpha r^{(k)}:\alpha \in ℝ\}$ que tem a mesma direção de $\nabla J\left(x^{(k)}\right)$ e passa pelo ponto $x^{(k)}$. O procedimento de escolher um ${\alpha }^{(k)}$ entre todos os possíveis, é conhecido como pesquisa linear (em inglês, line search).

A cada iteração, queremos escolher ${\alpha }_k$ de forma que $J\left(x^{(k+1)}\right)\le J\left(x^{(k)}\right)$. Isso pode ser garantido fazendo a seguinte escolha²⁴²⁴endnote: ²⁴Chamada de pesquisa linear exata. Qualquer outra escolha para $\alpha$ é conhecida como pesquisa linear não exata.

A fim de resolver este problema de minimização, vamos denotar

Então, observamos que

Agora, usando o fato de $A$ ser simétrica, obtemos

a qual, é uma função quadrática. Seu único mínimo, ocorre quando

Logo, encontramos

Com isso, temos a iteração do Método do Gradiente

Em revisão

Em revisão

O método do gradiente consiste em uma iteração da forma

com $p^{(k)}=r^{(k)}$. Ou seja, a nova aproximação $x^{(k+1)}$ é buscada na direção de $p^{(k)}$. Aqui, a ideia é usar uma melhor direção para buscar a solução.

O método do gradiente conjugado é um método de gradiente que busca encontrar a solução de $Ax=b$ pela computação do mínimo do seguinte funcional²⁷²⁷endnote: ²⁷Compare com o funcional $J$ dado em (1.87).

onde, denota o produto interno padrão e

é o produto interno induzido por $A$, lembrando que $A$ é positiva definida²⁸²⁸endnote: ²⁸Mostre que é de fato um produto interno.. Associada a este produto interno, temos a norma

chamada de norma da energia. O produto interno associado é também conhecido como produto interno da energia. Com isso, definimos que dois vetores $x$ e $y$ são conjugados, quando eles são ortogonais com respeito ao produto interno da energia, i.e. quando

Aqui, a ideia é desenvolver um método iterativo em que o erro a cada passo seja conjugado a todas as direções de busca anteriores. Consulte o desenvolvimento detalhado do método em [11, Seção 7.7].

Em revisão

Em revisão

A ideia é tomar $P$ igual a uma fatoração LU incompleta (ILU, do inglês, Incomplete LU). Incompleta no sentido que entradas de $L$ e de $U$ sejam adequadamente removidas, buscando-se uma boa esparsidade e ao mesmo tempo uma boa aproximação $LU$ para $A$.

Em revisão

O precondicionamento ILU(0) impõe que as matrizes $L$ e $U$ tenham o mesmo padrão de esparsidade da matriz $A$.

Em revisão