Uma matriz simétrica A pode ser decomposta na forma A=L D L^{T}, onde L é triangular inferior com l_{j j}=1 e D diagonal.
a) Descreva um procedimento para determinar L e D.
b) Descreva um procedimento para resolver o sistema linear A x=b utilizando essa decomposição.
Passo 1
O enunciado dá uma matriz simétrica A e afirma que ela pode ser escrita como
A = L D L^{T},
onde L é triangular inferior com l_{j j}=1 para todo j e D é diagonal. O exercício pede duas coisas:
(a) descrever um procedimento para determinar L e D (ou seja, um algoritmo prático para construir essa decomposição a partir de A);
(b) descrever um procedimento para resolver o sistema linear A x=b usando essa decomposição.
Informações importantes dadas: A é simétrica; L tem uns na diagonal; D é diagonal. A resolução começa a seguir.
Passo 2
Entendendo a multiplicação L D L^{T} elemento a elemento. Se L=(\ell_{i s}) e D=\operatorname{diag}(d_{1},\dots,d_{n}), então o produto dá, para quaisquer índices i,j,
Explicação didática: o termo L D é uma matriz cujas colunas s são as colunas de L multiplicadas por d_{s}. Ao multiplicar por L^{T} (que tem entradas \ell_{j s} transpostas), o elemento (i,j) é a soma dos produtos entre a linha i de L D e a linha j de L (ou coluna j de L^{T}), o que resulta na soma acima. Note que o somatório vai até \min(i,j) porque L é triangular inferior (entradas acima da diagonal são zero).
Passo 3
Derivação das fórmulas recursivas para obter D e L. Vamos construir L e D coluna a coluna (ou pivô a pivô), para k=1,2,\dots,n.
Para k=1: como não há termos anteriores,
d_{1}=A_{11}.
E para cada i>1,
\ell_{i1}=\frac{A_{i1}}{d_{1}}.
(Recorde: \ell_{11}=1 por definição.)
Para passo genérico k: compare a expressão de A_{k k} e de A_{i k} com a soma acima, separando os termos já conhecidos (com s<k) do novo termo s=k.
Para a diagonal:
Verifique que d_{k}\neq 0. (Se for zero, é necessário pivotamento ou a decomposição sem permutações não existe; para matrizes simétricas definida-positiva isso não ocorre e d_{k}>0.)
Defina \ell_{k k}=1. As entradas \ell_{i j} com j>i permanecem zero.
Comentário didático: o algoritmo é “coluna por coluna”: primeiro encontra-se o primeiro pivô d_{1} e a primeira coluna de L; depois o segundo pivô ajustando pelo que já foi computado, e assim por diante. O custo é da ordem de n^{3}/3 operações aritméticas (mesma classe de custo da fatoração LU sem aproveitamento de simetria), mas com metade da memória em comparação com uma fatoração LU completa porque exploramos a simetria.
Passo 5
Observação sobre estabilidade e degenerescência: se em algum passo d_{k}=0 o algoritmo pára (não dá para dividir). Para matrizes simétricas que não são definidas positivas pode ser preciso usar pivotamento simétrico (permutar linhas/colunas simultaneamente) para continuar. Para matrizes simétricas definidas positivas (caso comum em muitos problemas) todos os d_{k}>0 e o procedimento funciona sem pivotamento.
Passo 6
Procedimento para resolver A x=b usando A=L D L^{T}. A ideia é resolver em três fases simples, aproveitando que L é triangular inferior com uns na diagonal e D é diagonal:
Resolva L\,y=b por substituição progressiva (forward substitution).
Fórmula: para i=1,\dots,n,y_{i}=b_{i}-\sum_{j=1}^{i-1}\ell_{i j}\,y_{j}.(Aqui usamos \ell_{i i}=1.)
Resolva D\,z=y facilmente dividindo por cada elemento diagonal:z_{i}=\frac{y_{i}}{d_{i}},\qquad i=1,\dots,n.
Resolva L^{T}x=z por substituição regressiva (back substitution). Note que L^{T} é triangular superior com (L^{T})_{i j}=\ell_{j i}. Para i=n,n-1,\dots,1:x_{i}=z_{i}-\sum_{j=i+1}^{n}\ell_{j i}\,x_{j}.
Essas três etapas custam apenas O(n^{2}) operações após a fatoração: a fatoração é a parte custosa (O(n^{3})), e cada resolução com o mesmo A e diferentes b reutiliza L,D.
Passo 7
Exemplo numérico curto para fixar a ideia (ilustrativo). Considere a matriz simétrica
(b) Para resolver A x=b com A=L D L^{T} faça: (i) L y=b por forward substitution \;y_{i}=b_{i}-\sum_{j=1}^{i-1}\ell_{i j}y_{j}; (ii) z=D^{-1}y com z_{i}=y_{i}/d_{i}; (iii) L^{T}x=z por back substitution \;x_{i}=z_{i}-\sum_{j=i+1}^{n}\ell_{j i}x_{j}.