Método da matriz tridiagonal

Cálculo Numérico - Versão GNU Octave

4.5 Método da matriz tridiagonal

O método da matriz tridiagonal ou algoritmo de Thomas⁸ ou ainda TDMA (do inglês tridiagonal matrix algorithm) é o caso particular da eliminação gaussiana aplicada a matrizes tridiagonais.

Uma matriz tridiagonal é uma matriz quadrada cujos únicos elementos não nulos estão na diagonal principal e nas diagonais imediatamente acima e abaixo da principal. Um sistema tridiagonal é um sistema de equações lineares cuja matriz associada é tridiagonal, conforme a seguir: $\begin{array}{rcl} [\begin{matrix} b_{1} & c_{1} \\ a_{2} & b_{2} & c_{2} \\ a_{3} & b_{3} & ⋱ \\ ⋱ & ⋱ & c_{n - 1} \\ a_{n} & b_{n} \end{matrix}] [\begin{matrix} x_{1} \\ x_{2} \\ x_{3} \\ ⋮ \\ x_{n} \end{matrix}] = [\begin{matrix} d_{1} \\ d_{2} \\ d_{3} \\ ⋮ \\ d_{n} \end{matrix}] . & (4.91) \end{array}$

Observamos que não é necessário armazenar todos os $n^{2}$ elementos da matriz em memória, sendo suficiente armazenar os vetores $a_{n}$ , $b_{n}$ e $c_{n}$ . Por conveniência, a partir daqui, definiremos os elementos inexistentes na matriz $a_{1}$ e $c_{n}$ como zero:

a_{1} = c_{n} = 0 .

(4.92)

O algoritmo para a solução do sistema tridiagonal (4.91) pelo algoritmo de Thomas é dada pelas seguintes expressões: $\begin{array}{rcl} c_{i}^{'} & = & \{\begin{matrix} \begin{matrix} \frac{c_{i}}{b_{i}}, & i = 1 \\ \frac{c_{i}}{b_{i} - a_{i} c_{i - 1}^{'}}, & i = 2, 3, \dots, n - 1 \end{matrix} \end{matrix} & (4.93) \\ e \\ d_{i}^{'} & = & \{\begin{matrix} \begin{matrix} \frac{d_{i}}{b_{i}}, & i = 1 \\ \frac{d_{i} - a_{i} d_{i - 1}^{'}}{b_{i} - a_{i} c_{i - 1}^{'}}, & i = 2, 3, \dots, n . \end{matrix} \end{matrix} & (4.94) \end{array}$

Finalmente a solução final é obtida por substituição reversa: $\begin{array}{rcl} x_{n} & = & d_{n}^{'} & (4.95) \\ x_{i} & = & d_{i}^{'} - c_{i}^{'} x_{i + 1}, i = n - 1, n - 2, \dots, 1 . & (4.96) \end{array}$

Teorema 4.5.1. A aplicação da eliminação gaussiana sem pivotamento ao sistema (4.91) produz o algoritmo dado em (4.93) e (4.95).

Demonstração. O primeiro passo consiste em dividir todos os elementos da primeira linha de (4.91) por

b_{1}

\begin{array}{rcl} [\begin{matrix} 1 & c_{1}^{'} \\ a_{2} & b_{2} & c_{2} \\ a_{3} & b_{3} & ⋱ \\ ⋱ & ⋱ & c_{n - 1} \\ a_{n} & b_{n} \end{matrix}] [\begin{matrix} x_{1} \\ x_{2} \\ x_{3} \\ ⋮ \\ x_{n} \end{matrix}] = [\begin{matrix} d_{1}^{'} \\ d_{2} \\ d_{3} \\ ⋮ \\ d_{n} \end{matrix}], & (4.97) \end{array}

onde $c_{1}^{'} = \frac{c_{1}}{b_{1}}$ e $d_{1}^{'} = \frac{d_{1}}{b_{1}} .$

O segundo passo consiste em substituir a segunda linha por ela mesma subtraída da linha 1 multiplicada por $a_{2}$ ( $l_{2} \leftarrow l_{2} - a_{2} l_{1}$ ): $\begin{array}{rcl} [\begin{matrix} 1 & c_{1}^{'} \\ 0 & b_{2} - a_{2} c_{1}^{'} & c_{2} \\ a_{3} & b_{3} & ⋱ \\ ⋱ & ⋱ & c_{n - 1} \\ a_{n} & b_{n} \end{matrix}] [\begin{matrix} x_{1} \\ x_{2} \\ x_{3} \\ ⋮ \\ x_{n} \end{matrix}] = [\begin{matrix} d_{1}^{'} \\ d_{2} - a_{2} d_{1}^{'} \\ d_{3} \\ ⋮ \\ d_{n} \end{matrix}] . & (4.98) \end{array}$

Em seguida, dividimos a segunda linha por $b_{2} - a_{2} c_{1}^{'}$ , a fim de normalizar a diagonal principal: $\begin{array}{rcl} [\begin{matrix} 1 & c_{1}^{'} \\ 0 & 1 & c_{2}^{'} \\ a_{3} & b_{3} & ⋱ \\ ⋱ & ⋱ & c_{n - 1} \\ a_{n} & b_{n} \end{matrix}] [\begin{matrix} x_{1} \\ x_{2} \\ x_{3} \\ ⋮ \\ x_{n} \end{matrix}] = [\begin{matrix} d_{1}^{'} \\ d_{2}^{'} \\ d_{3} \\ ⋮ \\ d_{n} \end{matrix}] . & (4.99) \end{array}$

onde $c_{2}^{'} = \frac{c_{2}}{b_{2} - a_{2} c_{1}^{'}}$ e $d_{2}^{'} = \frac{d_{2} - a_{2} d_{1}^{'}}{b_{2} - a_{2} c_{1}^{'}}$ .

O próximo passo consiste em substituir a terceira linha por ela mesma subtraída da linha 2 multiplicada por $a_{3}$ ( $l_{3} \leftarrow l_{3} - a_{3} l_{2}$ ): $\begin{array}{rcl} [\begin{matrix} 1 & c_{1}^{'} \\ 0 & 1 & c_{2}^{'} \\ 0 & b_{3} - a_{3} c_{2}^{'} & ⋱ \\ ⋱ & ⋱ & c_{n - 1} \\ a_{n} & b_{n} \end{matrix}] [\begin{matrix} x_{1} \\ x_{2} \\ x_{3} \\ ⋮ \\ x_{n} \end{matrix}] = [\begin{matrix} d_{1}^{'} \\ d_{2}^{'} \\ d_{3} - a_{3} d_{2}^{'} \\ ⋮ \\ d_{n} \end{matrix}] . & (4.100) \end{array}$

A fim de normalizar o elemento da diagonal da terceira linha, dividimos toda a linha por $d_{3} - a_{3} d_{2}^{'}$ : $\begin{array}{rcl} [\begin{matrix} 1 & c_{1}^{'} \\ 0 & 1 & c_{2}^{'} \\ 0 & 1 & ⋱ \\ ⋱ & ⋱ & c_{n - 1} \\ a_{n} & b_{n} \end{matrix}] [\begin{matrix} x_{1} \\ x_{2} \\ x_{3} \\ ⋮ \\ x_{n} \end{matrix}] = [\begin{matrix} d_{1}^{'} \\ d_{2}^{'} \\ d_{3}^{'} \\ ⋮ \\ d_{n} \end{matrix}] . & (4.101) \end{array}$

Este procedimento é realizado até que se atinja a última linha e temos o seguinte sistema: $\begin{array}{rcl} [\begin{matrix} 1 & c_{1}^{'} \\ 0 & 1 & c_{2}^{'} \\ 0 & 1 & ⋱ \\ ⋱ & ⋱ & c_{n - 1}^{'} \\ 0 & 1 \end{matrix}] [\begin{matrix} x_{1} \\ x_{2} \\ x_{3} \\ ⋮ \\ x_{n} \end{matrix}] = [\begin{matrix} d_{1}^{'} \\ d_{2}^{'} \\ d_{3}^{'} \\ ⋮ \\ d_{n}^{'} \end{matrix}] . & (4.102) \end{array}$

Neste estágio, podemos encontrar os $x_{n}$ através de substituição reversa, isto é: a última linha diz

x_{n} = d_{n}^{'} .

(4.103)

A penúltima linha diz

x_{n - 1} + c_{n - 1}^{'} x_{n} = d_{n - 1}^{'} \Rightarrow x_{n - 1} = d_{n - 1}^{'} - c_{n - 1}^{'} x_{n} .

(4.104)

Esse mesmo procedimento aplicado à linha $i = 1, \dots n - 1$ , nos dá

x_{i} = d_{i}^{'} - c_{i}^{'} x_{i + 1} .

(4.105)

$■$

Exemplo 4.5.1. Considere a resolução do seguinte sistema tridiagonal pelo algoritmo de Thomas: $\begin{array}{rcl} [\begin{matrix} 2 & 1 & 0 & 0 & 0 \\ 1 & 2 & 1 & 0 & 0 \\ 0 & 1 & 2 & 1 & 0 \\ 0 & 0 & 1 & 2 & 1 \\ 0 & 0 & 0 & 1 & 2 \end{matrix}] [\begin{matrix} x_{1} \\ x_{2} \\ x_{3} \\ x_{4} \\ x_{5} \end{matrix}] = [\begin{matrix} 4 \\ 4 \\ 0 \\ 0 \\ 2 \end{matrix}] . & (4.106) \end{array}$

Primeiramente identificamos os vetores $a$ , $b$ , $c$ e $d$ : $\begin{array}{rcl} a & = & (0, 1, 1, 1, 1) & (4.107) \\ b & = & (2, 2, 2, 2, 2) & (4.108) \\ c & = & (1, 1, 1, 1, 0) & (4.109) \\ d & = & (4, 4, 0, 0, 2) & (4.110) \end{array}$

Agora, calculamos os vetores $c^{'}$ e $d^{'}$ : $\begin{array}{rcl} c_{1}^{'} & = & \frac{c_{1}}{b_{1}} = \frac{1}{2} & (4.111) \\ c_{2}^{'} & = & \frac{c_{2}}{b_{2} - a_{2} c_{1}^{'}} = \frac{1}{2 - 1 \cdot \frac{1}{2}} = \frac{2}{3} & (4.112) \\ c_{3}^{'} & = & \frac{c_{3}}{b_{3} - a_{3} c_{2}^{'}} = \frac{1}{2 - 1 \cdot \frac{2}{3}} = \frac{3}{4} & (4.113) \\ c_{4}^{'} & = & \frac{c_{4}}{b_{4} - a_{4} c_{3}^{'}} = \frac{1}{2 - 1 \cdot \frac{3}{4}} = \frac{4}{5} & (4.114) \\ d_{1}^{'} & = & \frac{d_{1}}{b_{1}} = \frac{4}{2} = 2 & (4.115) \\ d_{2}^{'} & = & \frac{d_{2} - a_{2} d_{1}^{'}}{b_{2} - a_{2} c_{1}^{'}} = \frac{4 - 1 \cdot 2}{2 - 1 \cdot \frac{1}{2}} = \frac{4}{3} & (4.116) \\ d_{3}^{'} & = & \frac{d_{3} - a_{3} d_{2}^{'}}{b_{3} - a_{3} c_{2}^{'}} = \frac{0 - 1 \cdot \frac{4}{3}}{2 - 1 \cdot \frac{2}{3}} = - 1 & (4.117) \\ d_{4}^{'} & = & \frac{d_{4} - a_{4} d_{3}^{'}}{b_{4} - a_{4} c_{3}^{'}} = \frac{0 - 1 \cdot (- 1)}{2 - 1 \cdot \frac{3}{4}} = \frac{4}{5} & (4.118) \\ d_{5}^{'} & = & \frac{d_{5} - a_{5} d_{4}^{'}}{b_{5} - a_{5} c_{4}^{'}} = \frac{2 - 1 \cdot \frac{4}{5}}{2 - 1 \cdot \frac{4}{5}} = 1 & (4.119) \end{array}$

Finalmente, calculamos o vetor $x$ : $\begin{array}{rcl} x_{5} & = & d_{5}^{'} = 1 & (4.120) \\ x_{4} & = & d_{4}^{'} - c_{4}^{'} \cdot x_{5} = \frac{4}{5} - \frac{4}{5} \cdot 1 = 0 & (4.121) \\ x_{3} & = & d_{3}^{'} - c_{3}^{'} \cdot x_{4} = - 1 - \frac{3}{4} \cdot 0 = - 1 & (4.122) \\ x_{2} & = & d_{2}^{'} - c_{2}^{'} \cdot x_{3} = \frac{4}{3} - \frac{2}{3} \cdot (- 1) = 2 & (4.123) \\ x_{1} & = & d_{1}^{'} - c_{1}^{'} \cdot x_{2} = 2 - \frac{1}{2} \cdot 2 = 1 & (4.124) \end{array}$

E assim, obtemos o vetor $x = [1, 2, - 1, 0, 1]$ .

Código GNU Octave: Método da matriz tridiagonal

function [x] = TDMA(a,b,c,d)

  n=size(a,1);   %Recupera ordem do sistema.

  cl=zeros(n,1); %Inicializa cl
  dl=zeros(n,1); %Inicializa dl
  x=zeros(n,1);  %Inicializa x

  cl(1)=c(1)/b(1);
  for i=2:n-1
    cl(i)=c(i)/(b(i)-a(i)*cl(i-1));
  end

  dl(1)=d(1)/b(1);
  for i=2:n
    dl(i)=(d(i)-a(i)*dl(i-1))/(b(i)-a(i)*cl(i-1));
  end

  x(n)=dl(n)   ;
  for i=n-1:-1:1
  x(i)=dl(i)-cl(i)*x(i+1);
  end
endfunction

Nesse código, usou-se cl e dl para denotar $c^{'}$ e $d^{'}$ . Observe que se for desnecessário preservar os valores originais dos vetores $c$ e $d$ , eles podem, com economia de memória e simplicidade de código, ser sobrescritos pelos vetores $c^{'}$ e $d^{'}$ , respectivamente. Eis uma nova implementação:

function [x] = TDMA2(a,b,c,d)

  n=size(a,1);   %Recupera ordem do sistema.
  x=zeros(n,1);  %Inicializa x

  c(1)=c(1)/b(1);
  for i=2:n-1
    c(i)=c(i)/(b(i)-a(i)*c(i-1));
  end

  d(1)=d(1)/b(1);
  for i=2:n
    d(i)=(d(i)-a(i)*d(i-1))/(b(i)-a(i)*c(i-1));
  end

  x(n)=d(n)   ;
  for i=n-1:-1:1
  x(i)=d(i)-c(i)*x(i+1);
  end
endfunction

A solução do sistema do Exemplo 4.5.1 pode ser obtida através dos seguintes comandos:

>>a=[0; 1; 1; 1; 1]
>>b=[2; 2; 2; 2; 2]
>>c=[1; 1; 1; 1; 0]
>>d=[4; 4; 0; 0; 2]
>>TDMA(a,b,c,d)

E 4.5.1. Considere o problema linear tridiagonal dado por $\begin{array}{rcl} [\begin{matrix} 5 & 4 & 0 & 0 & 0 & 0 \\ 1 & 3 & 1 & 0 & 0 & 0 \\ 0 & 2 & 4 & 1 & 0 & 0 \\ 0 & 0 & 1 & 2 & 1 & 0 \\ 0 & 0 & 0 & 2 & 3 & 2 \\ 0 & 0 & 0 & 0 & 1 & 2 \end{matrix}] [\begin{matrix} x_{1} \\ x_{2} \\ x_{3} \\ x_{4} \\ x_{5} \\ x_{6} \end{matrix}] = [\begin{matrix} 13 \\ 10 \\ 20 \\ 16 \\ 35 \\ 17 \end{matrix}] . & (4.125) \end{array}$

Identifique os vetores $a$ , $b$ , $c$ e $d$ relativos ao algoritmo da matriz tridiagonal. Depois resolva o sistema usando o computador.

Resposta.

\begin{array}{rcl} a & = & (0, 1, 2, 1, 2, 1) & (4.126) \\ b & = & (5, 3, 4, 2, 3, 2) & (4.127) \\ c & = & (4, 1, 1, 1, 2, 0) & (4.128) \\ d & = & (13, 10, 20, 16, 35, 17) & (4.129) \\ x & = & (1, 2, 3, 4, 5, 6) & (4.130) \end{array}

E 4.5.2. Considere o seguinte sistema de equações lineares: $\begin{array}{rcl} x_{1} - x_{2} & = & 0 \\ - x_{j - 1} + 5 x_{j} - x_{j + 1} & = & cos (j ∕ 10), 2 \leq j \leq 10 \\ x_{11} & = & x_{10} ∕ 2 & (4.131) \end{array}$

Identifique os vetores $a$ , $b$ , $c$ e $d$ relativos ao algoritmo da matriz tridiagonal no sistema linear dado. Depois resolva o sistema usando o computador. Veja também Exercício 4.7.4

Resposta.

\begin{array}{rcl} a & = & (0, - 1, - 1, - 1, - 1, - 1, - 1, - 1, - 1, - 1, - 1 ∕ 2) & (4.132) \\ b & = & (1, 5, 5, 5, 5, 5, 5, 5, 5, 5, 1) & (4.133) \\ c & = & (- 1, - 1, - 1, - 1, - 1, - 1, - 1, - 1, - 1, - 1, 0) & (4.134) \\ d & = & (0, cos (2 ∕ 10), cos (3 ∕ 10), cos (4 ∕ 10), cos (5 ∕ 10), & (4.135) \\ cos (6 ∕ 10), cos (7 ∕ 10), cos (8 ∕ 10), cos (9 ∕ 10), cos (1), 0) & (4.136) \\ x & = & (0, 324295, 0, 324295, 0, 317115, 0, 305943, 0, 291539, & (4.137) \\ 0, 274169, 0, 253971, 0, 230846, 0, 20355, 0, 165301, 0, 082650) & (4.138) \end{array}

Este texto é disponibilizado nos termos da licença Creative Commons Atribuição-CompartilhaIgual 3.0 Não Adaptada (CC-BY-SA 3.0). Página gerada em 19/8/2020 às 17:36:34.

Recursos

Projeto

IME - UFRGS
- Página do IME
- Página da UFRGS

UFRGS - IME - Recursos Educacionais Abertos de Matemática. Contato: reamat@ufrgs.br.

Informe erros ou edite você mesmo!