Apostila de MAT 271

Prévia do material em texto

UFV
Valentín Mendoza
 Cálculo Numérico
Apresentação
Estas notas de aula de Cálculo Numérico foram escritas no semestre 2013-II quando
ministrei esta disciplina para os estudantes da Universidade Federal de Viçosa. Em
princı́pio foram apresentados em forma de 6 apostilas, cada uma correspondendo a um
capı́tulo do conteúdo. Agora serão apresentadas juntas acrescentando os exercı́cios
de cada tópico.
O objetivo ao escrevê-las é proporcionar ao aluno os principais conceitos da análise
numérica que lhe permitam resolver os problemas que se lhe apresentarão no decorrer
de sua atividade de trabalho. Tivemos especial interesse em apresentar a dedução dos
métodos, seguindo as ideias gerais nas quais está baseado cada uma das técnicas
aqui mostradas.
Não está demais dizer que toda a tecnologia de nosso tempo está muito ligada ao
Cálculo Numérico. Sem ele seria impossı́vel calcular as probabilidades da distribuição
normal, resolver equações de difusão, trabalhar com dinâmica de fluidos, governar o
movimento de robots, etc. Pois existem modelos de fenômenos fı́sicos, econômicos,
biológicos, etc, que são analı́ticamente insolúveis, ou cuja solução é tão complicada
como para utilizá-la para fins práticos. Nesses casos o Cálculo Numérico constitui-se a
única ferramenta para abordá-los e para obter soluções aceitáveis. Por isso, se o aluno
sente que a disciplina é um pouco entediante quando utiliza a calculadora para realizar
contas sobre contas num determinado exercı́cio, deve ser consciente que isso é um
pequeno preço a pagar pelo conhecimento e a práxis que adquirirá e que, na hora das
aplicações, serão inestimáveis.
Fico devendo a parte computacional sem a qual a potência do Cálculo Numérico
não é apreciada.
Acredito que este trabalho contenha alguns erros, mas estes seriam muitos mais
se não fosse pelas correições que alguns dos alunos fizeram ao texto inicial. Gostaria
de agradecer a Daniela Abrantes Leal, Beatryz Cardoso Mendes, Gabriel de Rezende
Coelho, Franco Luiz Alves e Marcus Vinicius Miranda, que me indicaram alguns erros
tipográficos ou de conteúdo. Mas a pesar do dedicado esforço deles, escaparam-se
algumas coisas que podem ser melhoradas. Por isso, se o leitor tiver comentários,
crı́ticas e sugestões pode encaminhá-las ao email valentin@ufv.br. Serei grato a
quem comunicar os erros de qualquer natureza encontrados no texto.
Viçosa, março de 2014.
J. Valentı́n Mendoza M.
2
Sumário
1 Métodos numéricos para a solução de equações 5
1.1 Introdução . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1.2 Preliminares . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.3 O Método de Bisseção . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.3.1 Estimativa do número de iterações . . . . . . . . . . . . . . . . . . . 8
1.4 Método da Falsa Posição . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
1.4.1 Erro Absoluto no método da Falsa Posição . . . . . . . . . . . . . . 10
1.4.2 Erro Relativo no método da Falsa Posição . . . . . . . . . . . . . . . 10
1.5 O Método das Aproximações Sucessivas . . . . . . . . . . . . . . . . . . . . 11
1.6 O Método de Newton . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
1.6.1 Convergência do Método de Newton . . . . . . . . . . . . . . . . . 15
1.7 Método da Secante . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
1.8 Exercı́cios . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
2 Interpolação Polinomial 22
2.1 Introdução . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22
2.2 Preliminares . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23
2.3 A Fórmula de Lagrange . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
2.4 A Fórmula de Newton . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
2.5 Exercı́cios . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
3 Resolução Numérica de Sistemas Lineares 34
3.1 Introdução . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
3.2 Preliminares . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
3.2.1 Sistemas Triangulares Inferiores . . . . . . . . . . . . . . . . . . . . 34
3.2.2 Sistemas Triangulares Superiores . . . . . . . . . . . . . . . . . . . . 35
3.3 Método Direto: Decomposição LU . . . . . . . . . . . . . . . . . . . . . . . 36
3.4 Método Iterativo: O método de Gauss-Seidel . . . . . . . . . . . . . . . . . 41
3.4.1 O Caso Geral . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43
3.4.2 Critérios de Convergência . . . . . . . . . . . . . . . . . . . . . . . . 44
3.4.3 Critério de Parada . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44
3.5 Exercı́cios . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
4 Integração Numérica 49
4.1 Introdução . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
4.2 Preliminares . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
4.2.1 Estudo do Erro . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
4.3 Regra dos Trapézios. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
4.3.1 Limitante superior do erro na Regra do Trapézio . . . . . . . . . . . 53
4.3.2 Regra do Trapézio Generalizada . . . . . . . . . . . . . . . . . . . . 54
4.4 Regra 13 de Simpson . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56
4.4.1 Limitante superior para o Erro E2 na regra 1/3 de Simpson . . . . . 57
4.4.2 Regra 1/3 de Simpson Generalizada . . . . . . . . . . . . . . . . . . 60
4.5 Regra 38 de Simpson . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
4.5.1 Limitante superior para o erro na regra 3/8 de Simpson . . . . . . . 63
4.5.2 Regra 3/8 de Simpson Generalizada . . . . . . . . . . . . . . . . . . 64
4.6 Exercı́cios . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66
3
5 Resolução Numérica de Equações Diferenciais Ordinárias 69
5.1 Introdução . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
5.2 Preliminares . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
5.3 O Método de Euler . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
5.4 Métodos de Runge-Kutta . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
5.4.1 Métodos de Runge-Kutta de Ordem 2 . . . . . . . . . . . . . . . . . 75
5.4.2 Métodos de Runge Kutta de Ordem 3 . . . . . . . . . . . . . . . . . 78
5.4.3 O Método de Runge-Kutta de Ordem 4 . . . . . . . . . . . . . . . . 81
5.5 Sistemas de Equações Diferenciais Ordinárias . . . . . . . . . . . . . . . . . 82
5.5.1 Aplicação a equações diferenciais de segunda ordem . . . . . . . . 84
5.6 Exercı́cios . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87
6 Resolução Numérica de Equações Diferenciais Parciais 90
6.1 Introdução . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90
6.2 Preliminares . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90
6.2.1 Equações Diferenciais Parciais . . . . . . . . . . . . . . . . . . . . . 90
6.2.2 Aproximações da primeira e da segunda derivada . . . . . . . . . . 91
6.3 Equações Parabólicas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92
6.3.1 O modelo de Equação Parabólica . . . . . . . . . . . . . . . . . . . . 92
6.3.2 Discretização . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 93
6.3.3 Método Explı́cito . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 94
6.3.4 O Método Implı́cito . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98
6.3.5 O Método de Crank-Nicolson . . . . . . . . . . . . . . . . . . . . . .101
6.3.6 Condições de Fronteira de Neumann ou com Derivadas . . . . . . 105
6.4 Equações Elı́pticas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108
6.4.1 O modelo de Equação Elı́ptica . . . . . . . . . . . . . . . . . . . . . 108
6.5 Equações não lineares . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112
6.6 Exercı́cios . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 113
4
1 Métodos numéricos para a solução de equações
1.1 Introdução
Muitos problemas em diferentes áreas de conhecimento reduzem-se à solução da equação:
f (x) = 0 (1.1)
em que f : [a, b]→ R é uma função contı́nua. Um valor x̄ que satisfaz f (x̄) = 0 é chamado
uma raiz ou um zero de f (x).
É conhecido que se f (x) for um polinômio de grau menor ou igual a 4 existe uma
fórmula algébrica para determinar todas as suas raizes; mas se f (x) for um polinômio de
grau maior ou igual a 5 ou uma função envolvendo funções transcendentais, não existe
um método geral que permita encontrar exatamente as raizes da equação f (x) = 0.
Exemplo 1.1. A figura abaixo mostra o sol na origem e a terra no ponto (1, 0). Existem
5 pontos de equilibrio, no plano determinado pela órbita da terra ao redor do sol, nos
quais um satélite permanece imóvil em relação à terra. Estes pontos são chamados pontos
de Lagrange em honor ao matemático francês Joseph Louis Lagrange quem os descobriu
estudando o problema restrito dos três corpos. Estes pontos são obtidos quando as forças
(gravitacionais, centrı́petas ) que agem no satélite se contrabalancearem. Em dezembro
de 2005 o satélite de pesquisa solar SOHO (Solar and Heliospheric Observatory) foi
colocado no ponto L1. Este ponto L1 é uma boa posição para monitorar o sol desde que
as correntes de partı́culas do sol, como o vento solar, atingem L1 aproximadamente uma
hora antes de atingir a terra.
Se m1 é a massa do sol, m2 é a massa da terra e r =
m2
m1+m2
, então a coordenada x de L1
é a única raiz da equação de grau 5:
p(x) = x5 − (2 + r)x4 + (1 + 2r)x3 − (1 − r)x2 + 2(1 − r)x + r − 1 = 0
Usando o valor aproximado de r ≈ 3.04042 × 10−6, encontre a localização do ponto
L1.
Sol Terra
Figura 1: Pontos de Lagrange no sistema Sol-Terra.
5
Muitos problemas de ciência e tecnologia, como o problema anterior, conduzem à
resolução de equações não lineares para as quais é difı́cil ou impossı́vel encontrar
uma solução exata.
A solução do problema anterior só é possı́vel utilizando um método numérico. No
que segue apresentaremos alguns métodos que nos ajudarão a resolver equações deste
tipo.
1.2 Preliminares
Um teorema que utilizaremos é o seguinte.
TEOREMA 1.1 (Teorema do Valor Intermediário (TVI)). Seja f (x) é uma função contı́nua
num intervalo [a, b]. Se f (a) f (b) < 0 então existe pelo menos uma raiz x̄ de f (x) = 0 entre a e b.
Em palavras, podemos garantir ao menos uma solução da equação f (x) = 0, desde
que as imagens dos extremos do intervalo [a, b] sejam de signais diferentes. Isto é, se
f (a) > 0 e f (b) < 0, ou se f (a) < 0 e f (b) > 0.
Exemplo 1.2. No exemplo 1.1, a função é
f (x) = p(x) = x5 − (2 + r)x4 + (1 + 2r)x3 − (1 − r)x2 + 2(1 − r)x + r − 1.
Observe que f (0) = r−1 < 0 e f (1) = r > 0 então, pelo TVI, existe uma raiz x̄ no intervalo
[0, 1].
Outro resultado que precisaremos é o seguinte teorema:
TEOREMA 1.2 (Teorema do Valor Médio (TVM)). Seja f (x) uma função contı́nua em [a, b]
e diferenciável no intervalo (a, b). Então existe um ponto c ∈ (a, b) tal que
f (b) − f (a) = f ′(c)(b − a).
O TVI é a ferramenta principal do primeiro método numérico.
1.3 O Método de Bisseção
Suponha f (x) contı́nua e que existe um intervalo [a, b] tal que f (a) f (b) < 0. Sem perda
de generalidade, podemos supor f (a) > 0 e f (b) < 0. Pelo TVI, existe uma raiz x̄ de
f (x) = 0 entre a e b, que pode se qualquer ponto do intervalo [a, b]. Assim, temos pouca
informação sobre o valor da raiz. A ideia do método de Bisseção é diminuir o tamanho
do intervalo [a, b] com o intuito de conseguir um intervalo menor no qual podemos
encontrar a raiz.
Como não temos ideia da ubicação da raiz, democraticamente, escolhemos como
primeira aproximação de x̄, o ponto médio do intervalo:
x1 =
a + b
2
6
Figura 2: O Método da Bisseção.
Se por acaso obtemos que f (x1) = 0, então o problema ficaria resolvido pois neste caso
x̄ = x1. Se f (x1) , 0 então x̄ encontra-se no interior de [a, x1] ou no interior de [x1, b].
Para determinar em qual destes intervalos está a raiz, utilizamos o TVI:
(i) Se f (a) f (x1) < 0 então a raiz encontra-se em [a, x1],
(ii) Se f (x1) f (b) < 0 então a raiz encontra-se em [x1, b].
No caso (i), fazemos b = x1 e no caso (ii), fazemos a = x1. Em qualquer caso, temos
diminuido o tamanho do intervalo à metade. Com este novo intervalo repetimos o
processo o obtemos uma segunda aproximação x2 com a qual, se não é a raiz, diminuimos
o intervalo novamente na metade. Dependendo do valor de f (x2) teremos a = x2 ou
b = x2.
Continuando este processo obtemos uma sequência de aproximações xn+1, para
n = 0, 1, 2, ... que pertencem a intervalos In = [an, bn] com a0 = a e b0 = b de forma que
f (an) f (bn) < 0 e o comprimento de cada Ii+1 é a metade do comprimento de Ii. Assim
bn − an =
bn−1 − an−1
2
= ... =
b0 − a0
2n
. (1.2)
e
xn+1 =
an + bn
2
.
O seguinte Teorema mostra a convergência das aproximações xn à raiz x̄.
TEOREMA 1.3. A sequência de aproximações xn obtidas pelo Método da Bisseção converge para
a raiz x̄.
Demonstração. Pela equação (1.2) temos que
lim
n→∞
|bn − an| = limn→∞
b0 − a0
2n
= 0.
Por tanto limn→∞ bn = limn→∞ an = L onde L é um número real. Como an < xn+1 <
7
bn temos que limn→∞ xn = L. Provaremos que L = x̄. Observe que f (an) f (bn) < 0
então
lim
n→∞
f (an) f (bn) = f (L) f (L) ≤ 0,
o que implica que ( f (L))2 ≤ 0. Mas ( f (L))2 ≥ 0. Assim a única possibilidade é que
f (L)2 = 0, ou f (L) = 0. Dessa forma L = x̄. �
1.3.1 Estimativa do número de iterações
Desejamos determinar quantas iterações são necessárias como mı́nimo, para encontrar
um intervalo [an, bn], que contenha a aproximação xn+1, tal que
|xn+1 − x̄| < �,
em que � é o erro absoluto ao aproximar x̄. Pela equação (1.2)
|xn+1 − x̄| ≤ xn+1 − an =
bn − an
2
=
b0 − a0
2n+1
Logo para obter a precisão desejada, é suficiente considerar um n tal que
b0 − a0
2n+1
< �
Assim,
n >
Log(b0 − a0) − Log(�)
Log(2)
− 1.
Observação 1.1. Ao resolver um problema numericamente, não esqueça de colocar a
calculadora no modo radianos.
Exemplo 1.3. Aplique o método de Bisseção para resolver a equação x3 + cos(x) = 0
obtendo os extremos do intervalo inicial a e b pelo TVI. Encontre um resultado xn tal que
|xn − x̄| < 0.01.
Solução:
A função a avaliar é f (x) = x3 + cos(x). Note que f (−1) = −0.45969 e f (0) = 1, assim
podemos considerar a = −1 e b = 0. Observe que precisaremos de
n >
Log(0 + 1) − log(0.01)
Log2
− 1 = 5.64
8
iterações para obter o erro desejado. Assim, é suficiente n = 6. Na seguinte tabela
organizaremos os dados:
n a b xn+1 f (a) f (b) f (xn+1)
0 -1 0 -0.5 -0.45969 1 0.75258
1 -1 -0.5 -0.75 -0.45969 0.75258 0.309813
2 -1 -0.75 -0.875 -0.45969 0.309813 -0.02892
3 -0.875 -0.75 0 -0.8125 -0.02892 0.309813 0.1513
4 -0.875 -0.8125 -0.84375 -0.02892 0.1513 0.75258
5 -0.875 -0.84375 -.859375 -0.02892 0.75258 0.018240
6 -0.875 -0.859375 -0.8671875 -0.02892 0.018240 -0.00516
Em cada passo fomos utilizando o TVI para escolher os valores de a e b. Repare que
|x7 − x̄| ≤
b6 − a6
2
= 0.00781 < 0.01.
Assim uma aproximação da raiz é x7 = −0.8671875.
1.4 Método da Falsa Posição
O método da Falsa Posição é uma variação do método da Bisseção. Neste método mu-
damos a forma de escolher a aproximação xn+1. Na Biseção aproximavamos utilizando
o ponto médio do intervalo. Na Falsa posição aproximaremos a raiz utilizando a retaL
que passa pelos pontos (a, f (a)) e (b, f (b)). Ou seja, f (x) é aproximada por L e assim x1,
a aproximação de x̄, será o ponto onde a ordenada y da reta L é zero.
Figura 3: O Método da Falsa Posição.
A equação da reta L é:
y − f (a)
x − a
=
f (b) − f (a)
b − a
.
Se y = 0, então x = x1. Logo
0 − f (a)
x1 − a
=
f (b) − f (a)
b − a
.
9
Assim
x1 =
a f (b) − b f (a)
f (b) − f (a)
Tendo x1, podemos usar o mesmo critério da Bisseção para escolher o novo intervalo
[a, b], isto é:
(i) Se f (a) f (x1) < 0 então a raiz encontra-se em [a, x1],
(ii) Se f (x1) f (b) < 0 então a raiz encontra-se em [x1, b].
No caso (i) definimos b = x1, e no caso (ii) definimos a = x1. Logo, continuamos o
processo até obter a precisão desejada, obtendo intervalos [an, bn] com a0 = a e b0 = b e
xn+1 =
an f (bn) − bn f (an)
f (bn) − f (an)
. (1.3)
O Método da Bisseção utiliza só a informação dos valores a e b para encontrar xn+1,
entanto que a Falsa Posição utiliza a informação de a, b, f (a) e f (b).
1.4.1 Erro Absoluto no método da Falsa Posição
No método da Falsa Posição pode acontecer que an ou bn permanece fixo a partir de um
certo n, então não podemos garantir que bn − an tende para 0. Por isso utilizaremos o
erro absoluto |xn+1 − xn| < � como critério de parada, onde � é a precisão dada.
Por tanto, o método para quando atingimos a desigualdade |xn+1 − xn| < �.
1.4.2 Erro Relativo no método da Falsa Posição
Em alguns exemplos, o critério de parada será estabelecido pelo erro relativo. Se desejar-
mos um erro relativo menor que �, o método para uma vez que atingimos a desigualdade∣∣∣∣xn+1 − xnxn+1
∣∣∣∣ < �. (1.4)
Exemplo 1.4. Aplique o método da Falsa Posição para resolver a equação x3 +cos(x) = 0
com a = −1 e b = 0. Encontre um resultado xn+1 tal que |xn+1 − xn| < 0.01.
Solução:
Como f (x) = x3 + cos(x) e f (−1) = −0.45969 e f (0) = 1, consideramos a = −1 e b = 0.
O método para se atingimos |xn+1 − xn| < 0.01. Na seguinte tabela organizaremos os
dados:
10
n a b xn+1 f (a) f (b) f (xn+1)
0 -1 0 -0.68507 -0.45969 1 0.45285
1 -1 -0.68507 -0.84135 -0.45969 0.45285 0.07089
2 -1 -0.84135 -0.86254 -0.45969 0.07089 0.00880
3 -1 -0.86254 -0.86512 -0.45969 0.00880 0.00106
Em cada passo, uma vez calculado xn+1, utilizamos o TVI para escolher o intervalo
[a, b]. Note que
|x4 − x3| ≤ 0.00262 < 0.01.
Assim uma aproximação da raiz é x4 = −0.86512.
Observação 1.2. Em todos os métodos que seguem aplicamos um dos dois critérios
de parada:
(i) Erro absoluto: se |xn+1 − xn| < �;
(ii) Erro relativo: se |xn+1−xnxn+1 | < �.
1.5 O Método das Aproximações Sucessivas
Seja f (x) : [a, b] → R uma função contı́nua. O método de Aproximações Sucessivas ou
Método do Ponto Fixo consiste em escrever a equação f (x) = 0 na forma
x = φ(x) (1.5)
em que φ(x) é uma função contı́nua. Desse modo, achar uma raiz x̄ que satisfaça f (x̄) = 0
é equivalente a achar um x̄ tal que x̄ = φ(x̄).
Exemplo 1.5. Encontre uma função φ(x) tal que a equação x = φ(x) seja equivalente à
equação x3 − 13 = 0.
Solução:
Neste exemplo f (x) = x3 − 13. Existem várias maneiras de resolver o problema:
(i) Escreva x3 − 13 = 0 na forma x2 = 13x e então isole x na esquerda x =
√
13
x . Por
tanto φ(x) =
√
13
x .
(ii) Escreva x3 − 13 = 0 na forma 3x3 = 2x3 + 13. Divida entre 3x2 para obter
11
x = 2x
3+13
3x2 . Neste caso φ(x) =
2x3+13
3x2 .
Do exemplo anterior é claro que existem muitas, na verdade infinitas, maneiras de definir
uma função φ(x). Se temos escolhido uma destas funções φ e um ponto x0 ∈ [a, b],
construa a sequência
xn+1 = φ(xn) (1.6)
Se por acaso a sequência {xn} convergir para um número x∗, isto é, limn→∞ xn = x∗, então
temos que
lim
n→∞
xn+1 = limn→∞φ(xn).
Como φ é contı́nua,
lim
n→∞
xn+1 = φ( limn→∞ xn),
o que implica x∗ = φ(x∗). Pela definição da φ, segue que f (x∗) = 0, ou seja, x̄ = x∗.
A análise anterior mostra que podemos encontrar x̄ usando a sequência {xn} da
equação (1.6), desde que o limite limn→∞ xn exista.
O seguinte Teorema fornece condições para que o limite limn→∞ xn exista.
TEOREMA 1.4. Seja φ(x) uma função contı́nua e diferenciável num intervalo [a, b] no qual
vale
|φ′(x)| < 1. (1.7)
tal que a equação f (x) = 0 seja equivalente à equação x = φ(x). Então, dado x0 ∈ [a, b], a
sequência de iterações xn+1 = φ(xn) converge para a raiz x̄ de f (x) = 0.
Demonstração. Seja x̄ satisfazendo x̄ = φ(x̄) ou f (x̄) = 0. Pelo Teorema do Valor
Médio, para cada n existe um cn ∈ (x̄, xn) tal que
|φ(xn) − φ(x̄)| = |φ′(cn)||xn − x̄|.
Como |φ′(x)| < 1 para todo x ∈ [a, b] temos que |φ′(cn)| ≤M, para um certo 0 < M < 1.
Pela definição xn+1 = φ(xn), isto implica que
|xn+1 − x̄| ≤M|xn − x̄|.
Continuando este processo temos:
|xn+1 − x̄| ≤M|xn − x̄| ≤M2|xn−1 − x̄| ≤ . . . ≤Mn+2|x0 − x̄|.
12
Por tanto
lim
n→∞
|xn+1 − x̄| ≤ limn→∞M
n+2
|x0 − x̄| = 0
pois M < 1. Então,
lim
n→∞
|xn+1 − x̄| = 0
o que implica limn→∞ xn+1 = x̄. �
Observação 1.3. Pelo Teorema 1.4 devemos escolher uma φ que satisfaça:
|φ′(x)| < 1. (1.8)
Exemplo 1.6. A função f (x) = xe−x − e−3 possui exatamente uma raiz x̄ em [0.01, 1].
Encontre esta raiz pelo método de Aproximações Sucessivas, utilizando um erro absoluto
menor que 0.01.
Solução:
Escreva a equação xe−x − e−3 = 0 na forma xe−x = e−3, então divida entre e−x para
obter x = ex−3. Assim φ(x) = ex−3 e φ′(x) = ex−3. Como 0.01 ≤ x ≤ 1 temos
−2.99 ≤ x − 3 ≤ −2 e então
e−2.99 ≤ ex−3 ≤ e−2,
pois ex−3 é uma função crescente. Logo, |φ′(x)| ≤ e−2 ≈ 0.1353 < 1. Então φ(x)
satizfaz as condições do Teorema 1.4. Considere o ponto x0 = 0.5, então
x1 = φ(x0) = φ(0.5) = e0.5−3 = 0.08208.
Continuamos calculando os valores pela fórmula xn+1 = φ(xn) na seguinte tabela:
n xn
0 0.5
1 0.08208
2 0.05404
3 0.052551
Observe que |x3 − x2| = 0.00148 < 0.01. A solução aproximada é x3 = 0.052551.
13
1.6 O Método de Newton
A ideia para encontrar uma raiz da equação f (x) = 0 em todos os métodos numéricos
é começar com uma (ou mais) aproximação inicial x0 e então determinar a seguinte
aproximação seguindo uma regra dada. No método de Newton usa-se a reta tangente
para achar a segunda aproximação. Isto é feito da seguinte maneira.
Seja f : [a, b] → R uma função contı́nua em [a, b] e diferenciável em (a, b), e seja
x0 ∈ (a, b) uma aproximação inicial da raiz x̄. A reta tangente no ponto (x0, f (x0)) usa-se
para determinar a seguinte aproximação. Assim a nova aproximação será então o valor
de x1 no qual a ordenada da reta tangente é zero.
Reta tangente
Figura 4: O método de Newton.
A equação da reta tangente no ponto (x0, f (x0)) é:
y − f (x0)
x − x0
= f ′(x0).
Substituindo as condições x = x1 e y = 0 temos:
0 − f (x0)
x1 − x0
= f ′(x0).
Isolando x1, obtemos a fórmula do método de Newton no ponto x0:
x1 = x0 −
f (x0)
f ′(x0)
.
Uma vez calculado x1 podemos fazer o mesmo processo para achar outra aproximação
x2 a partir de x1. É claro que a fórmula para x2 é:
x2 = x2 −
f (x1)
f ′(x1)
.
14
Em forma geral, a aproximação (n + 1)-ésima é dada por
xn+1 = xn −
f (xn)
f ′(xn)
. (1.9)
1.6.1 Convergência do Método de Newton
O seguinte Teorema garante a convergência do Método de Newton.
TEOREMA 1.5. Sejam f (x), f ′(x) e f ′′(x) contı́nuas num intervalo I = [a, b] que contém a raiz
x̄ e suponha que f ′(x̄) , 0. Então existe um intervalo I∗ ⊂ I, contendo a raiz x̄ tal que se x0 ∈ I∗,
então a sequência {xn} gerada pela fórmula
xn+1 = xn −
f (xn)
f ′(xn)
convergirá para a raiz.
Demonstração. Observe que o método de Newton é um caso particular do
método das Aproximações sucessivas com φ(x) = x − f (x)f ′(x) pois xn+1 = φ(xn).
Assim, pelo Teorema 1.4, o método convergirá se |φ′(x)| < 1. Note que
φ′(x) =
f (x) f ′′(x)
[ f ′(x)]2
Como f (x̄) = 0, temos φ′(x̄) = 0. Pela continuidade de f , f ′ e f ′′, existe um
intervalo I∗ ⊂ I contendo x̄, no qual |φ′(x)| < 1 para x ∈ I∗. �
Observação1.4. O Teorema anterior diz que se começarmos com um ponto x0 próximo
da raiz x̄, então o Método de Newton convergirá para x̄.
Exemplo 1.7. Considere a equação não linear x2− ex = 0. Use o método de Newton para
resolvê-la considerando uma aproximação inicial de x0 = −2 e um erro absoluto menor
que 0.001.
Solução:
Dos dados do problema f (x) = x2 − ex e f ′(x) = 2x − ex. Usando a fórmula
xn+1 = xn −
f (xn)
f ′(xn)
,
15
escreveremos as iterações na tabela abaixo:
n xn f (xn) f ′(xn)
0 -2 3.86466 -4.13533
1 -1.06545 0.79061 -2.47547
2 -0.74607 0.08239 -1.96636
3 -0.70417 0.00133 -1.90285
4 -0.70347 0.000004 -1.90180
Repare que |x4 − x3| = 0.0007 < 0.001. Assim um valor aproximado da raiz é
x4 = −0.70347.
1.7 Método da Secante
O Método de Newton precisa da avaliação da derivada f ′(x) nos pontos xn que, depen-
dendo da função f (x), pode ser muito custosa. O Método da Secante surgue como uma
alternativa ao método de Newton no qual não é necessário avaliar a derivada da função.
A ideia é substituir a expressão f ′(xn) na fórmula de Newton por uma aproximação
dela. Suponha que conhecemos duas aproximações x0 e x1 da raiz x̄. Pelo método de
Newton
x2 = x1 −
f (x1)
f ′(x1)
.
Se não querermos derivar f , podemos aproximar a derivada f ′(x1) pela inclinação da
reta que passa pelos pontos (x0, f (x0)) e (x1, f (x1)):
f ′(x1) ≈
f (x1) − f (x0)
x1 − x0
.
Dessa forma:
x2 = x1 −
f (x1)[ f (x1)− f (x0)
x1−x0
]
Simplificando
x2 =
x0 f (x1) − x1 f (x0)
f (x1) − f (x0)
.
Fazendo a mesma análise para x1 e x2, obteremos uma fórmula para x3
x3 =
x1 f (x2) − x2 f (x1)
f (x2) − f (x1)
.
Em geral, a aproximação xn+1 depende de xn−1 e xn:
xn+1 =
xn−1 f (xn) − xn f (xn−1)
f (xn) − f (xn−1)
. (1.10)
16
Figura 5: O Método da Secante.
Observação 1.5. A fórmula do Método da Secante e do método da Falsa Posição é
a mesma. A diferença entre estes dois métodos é que no Método da Secante não é
necessário que os valores iniciais tenham signais opostos.
Exemplo 1.8. Utilize o Método da Secante para determinar uma raiz da equação x3 −
9x + 3 = 0, em [0, 1] com erro relativo menor que 10−3.
Solução:
Neste exemplo x0 = 0, x1 = 1 e f (x) = x3 − 9x + 3. Na tabela abaixo estão os valores
de xn e f (xn). Em cada passo xn+1 é dado pela fórmula
xn+1 =
xn−1 f (xn) − xn f (xn−1)
f (xn) − f (xn−1)
.
n xn−1 xn xn+1 f (xn−1) f (xn) f (xn+1)
1 0 1 0.375 3 -5 -0.32226
2 1 0.375 0.33194 -5 -0.32226 0.04911
3 0.375 0.33194 0.33763 -0.32226 0.04911 -0.00018
4 0.33194 0.33763 0.33760 0.04911 -0.00018 -0.00007
O erro relativo para
∣∣∣∣ x5−x4x5 ∣∣∣∣ = 0.00088 é menor que o erro procurado e assim uma
aproximação da raiz é x5 = 0.33760.
17
1.8 Exercı́cios
Exercı́cio 1.1. Localize graficamentemente as raizes da equação ln(x) + 5x − 6 − x2 = 0.
Separe cada uma delas em intervalos de comprimento máximo unitário.
Exercı́cio 1.2. Localize graficamente as raı́zes das equações:
(a) 4 cos(x) − e2x = 0; (b)
x
2
− tan(x) = 0;
(c) 1 − x ln(x) = 0; (d) 2x − 3x = 0.
Exercı́cio 1.3. Considere a equação f (x) = xex − 1.
(a) Mostre que a equação possui uma raiz no intervalo [0.5, 1].
(b) Qual é o número mı́nimos de iterações necessárias para se obter uma aproximação
da raiz com duas casas decimais corretas, usando o Método da Bisseção?
(c) Se (xk), k = 0, 1, 2, 3, .. é a sequência de aproximações para a solução da equação
f (x) = 0 através do Método da Biseção, determine o valor aproximado de |x4 − x3|.
Exercı́cio 1.4. Aplique o Método de Bisseção e da Falsa Posição para calcular a raiz positiva
de x3 − 15 = 0 com erro absoluto � < 0, 01, partindo do intervalo [2.00, 3.00].
Exercı́cio 1.5. Aplique o Método de Bisseção para resolver:
(a) ex − x − 3x2 = 0; (b) x3 + cos(x) = 0
obtendo os extremos do intervalo inicial a e b graficamente. Encontre um resultado xn
tal que |x̄ − xn| < 0, 01.
Exercı́cio 1.6. Calcular a raiz de 2x3 − cos(x + 1)− 3 = 0, pertencente ao intervalo [−1, 2],
com precisão de 0.01 usando o Método de Bisseção com no máximo 10 passos. Verifique
quantos passos no mı́nimo são necessários para ter uma precisão de 10−8.
Exercı́cio 1.7. Dado o polinómio p(x) = x3 + 2x2 + x − 2, faça o que se pede:
(a) Determine o intervalo onde todos os zeros do polinômio devem estar.
(b) Encontre uma aproximação para uma das raı́zes de p(x) = 0, com precisão de sete
casas decimais, usando o Método de Newton com no máximo 10 passos e x0 = 1.
Exercı́cio 1.8. Faça o exercı́cio anterior usando o polinômio p(x) = x3 + 4x2 + x − 2.
Exercı́cio 1.9. Considere a função f (x) = x3 − x − 1. Resolva a equação f (x) = 0 pelo
Método de aproximações sucessivas usando a função de iteração φ(x) = 1x +
1
x2 e x0 = 1.
Justifique seus resultados.
18
Exercı́cio 1.10. As raı́zes de f (x) = ln(x) − x + 2 podem ser determinadas usando o
processo iterativo na forma xk+1 = φ(xk), i = 1, 2, .... Encontre as raı́zes considerando:
(a) φ(x) = 2 + ln(x); (b) φ(x) = ex−2.
Exercı́cio 1.11. A função f (x) = xe−x − e−3, x ∈ R, possui exatamente duas raı́zes reais:
α1 ∈ [0.01, 1] e α2 ∈ [4, 5]. Considere as funções φ1(x) = ex−3 e φ2(x) = ln x + 3.
(a) φ1 pode ser usada para aproximar α1 pelo método das aproximações sucessivas
com garantia de convergência? E φ2?. Justifique.
(b) φ1 pode ser usada para aproximar α2 pelo método das aproximações sucessivas
com garantia de convergência? E φ2?. Justifique.
Exercı́cio 1.12. Mostre que x3−2x2−5− sin(x) = 0 tem apenas uma raiz real e determine
seu valor correto até 5 casas decimais usando o Método de Newton, com no máximo 6
iterações. Considere x0 = 2.1.
Exercı́cio 1.13. Mostre que a fórmula para determinar a raiz cúbica de Q,
xn+1 =
1
3
(2xn +
Q
x2n
),n = 0, 1, ....
é um caso especial do Método de Newton. Aplique o método para calcular a raiz cúbica
de 2 com precisão de 10−2 usando o erro relativo. Considere x0 = 1.
Exercı́cio 1.14. Calcular as duas raı́zes de sin(x) − ex − 2x2 + 10 = 0 usando o Método de
Newton, com a precisão de |xn+1 − xn| ≤ 10−5 e no máximo 10 passos.
Exercı́cio 1.15. Localize os zeros do polinômio p(x) = x5 − 109 x
3 + 521 . Usando o Método de
Newton, encontre as raı́zes não nulas de p(x) = 0 com precisão de 6 casas decimais.
Exercı́cio 1.16. Utilizando o Método de Newton e chamando de (xk), k = 0, 1, 2, ... a
sequência de aproximações para a solução da equação tan(x) − 2x = 0 no intervalo
[1, 1.5], com x0 = 1.3, encontre um valor aproximado de |x5 − x4| × 107.
Exercı́cio 1.17. Considere a equação não-linear x2 − ex = 0. Use o Método de Newton para
resolvê-la considerando uma aproximação inicial x0 = −2 e precisão menor que � = 10−1.
Use arredondamento e 4 casas decimais após a vı́rgula.
Exercı́cio 1.18. Utilizando o Método de Newton e chamando de (xk), k = 0, 1, 2, ... a
sequência de aproximações para a solução da equação x3 + x − λ = 0 no intervalo
[1, 2], com x0 = 1.3 e λ = 3, encontre um valor aproximado de |x5 − x4| × 107.
Exercı́cio 1.19. Além da solução obvia x = 0, a equação
x sin(x) + 3x cos(x) − 2x = 0
19
tem duas soluções no intervalo [−1, 2], sendo uma positiva e a outra negativa. Sejam
(Sk) e (Pk), k = 1, 2, ... as sequências de aproximações para as raı́zes de f (x) obtidas pelo
Método da Secante e da Falsa Posição, respectivamente. Considerando S1 = P1 = −1 e
S2 = P2 = 1.5, encontre o valor aproximado de |S5 − P5| × 100.
Exercı́cio 1.20. Utilize o Método da falsa posição para determinar a raiz da equação x log(x)−
1 = 0, em [2, 3] com � = 10−4.
Exercı́cio 1.21. Utilize o Método da falsa posição para determinar a raiz da equação x3 −
9x + 3 = 0, em [0, 1] com � = 10−4.
Exercı́cio 1.22. Utilize o Método da secante para determinar a raiz da equação x3−9x+3 = 0,
em [0, 1] com � = 10−4.
Exercı́cio 1.23. Considere um pórtico em L invertido como na figura onde K é o coefici-
ente de mola torsional. Encontre o valor do ángulode equilibrio θ quando aplicamos
uma força P no extremo do pórtico. Resolva o problema com K = 1, L = 2 e P = 2.
Figura 6: Pórtico em L.
Exercı́cio 1.24. Use o Método de Newton para determinar o ponto de mı́nimo global do
polinômio abaixo, se existir tal ponto.
p(x) = x(x − 1)(x + 1)(x − 2) (1.11)
Atenção: Primeiro verifique se existe ou não um ponto de mı́nimo global (isto é, um
ponto α ∈ R tal que p(α) ≤ p(x), ∀x ∈ R). Depois isole os possı́veis candidatos (se
existirem) e só então calcule, pelo método pedido, o ponto de mı́nimo global.
Exercı́cio 1.25. Verifique que a sequência
xn+1 = (
p − 1
p
)xn +
A
pxp−1n
(1.12)
pode ser utilizada para calcular p
√
A. Depois determine 4
√
9 com erro menor a 10−6.
20
Exercı́cio 1.26. Na figura, o comprimento da corda AB é de 4 cm e o comprimento do
arco AB é de 5 cm. Encontre o ângulo central θ, em radianos, correto até a quarta casa
decimal.
Figura 7: Corda e Arco.
Exercı́cio 1.27. Um tanque tem a forma de um cilindro reto, com raio igual a 1 e compri-
mento l. Sua lateral (circular) é transparente e através dela podemos observar o nı́vel do
liquido no cilindro (”deitado”). A porcentagem do liquido no fluido pode ser obtida em
função do ângulo (veja a figura). Por exemplo, o cilindro está cheio quando θ = π e pela
metade para θ = π2 . Calcule com erro menor de 10
−3 o valor de θ para o qual o cilin-
dro tem um quarto de seu volume cheio, através do método de Aproximações sucessivas.
Justifique todos os passos de forma a garantir a convergência.
Figura 8: Lateral do cilindro.
Exercı́cio 1.28. Um fábrica possui material para confecção de uma lata cilı́ndrica com
área total de superfı́cie de 800 cm2. Deseja-se uma lata com esta área de superfı́cie e
volume de 1 litro. Determine através de um Método de Aproximações Sucessivas qual o raio
da lata e sua correspondente altura. Efetue 3 iterações do método e estime o número
de iterações necessárias para um erro menor que 10−5. Justifique a convergência do
método. (Obs: deseja-se a lata mais alta...)
Exercı́cio 1.29. No exemplo (1.1), utilizando o Método de Newton , determine a posição
do Satélite SOHO. Utilize um erro absoluto menor que 10−10!.
21
2 Interpolação Polinomial
2.1 Introdução
Num sentido amplo, interpolar significa obter informação sobre um evento ou fenômeno
a partir de dados conhecidos. Em nosso caso, queremos obter informação sobre uma
certa função f (x) para a qual só conhecemos o seu valor yi = f (xi), i = 0, ...,n em n + 1
pontos denotados
x0 < x1 < .... < xn.
Como não conhecemos a expressão da f , podemos tentar encontrar uma função g(x) que
satisfaça as mesmas condições da f (x), ou seja,
g(xi) = f (xi) = yi, para todo i = 0, ...,n.
Isto é, desejamos encontrar uma função g(x) que coincida com f (x) pelo menos nos
pontos xi. A ideia é que g(x) possa substituir à f (x).
Figura 9: g(x) coincide com f (x) nos pontos xi, i = 0, 1, ...,n.
Existem várias justificativas para procurar uma função g(x) que substitua à função
f (x):
(a) Pode ser que não conhecemos f (x) exceto nos pontos xi, então se encontrar-
mos uma função g(x) que coincide com ela nestes pontos, podemos supor g
aproxima f , ou que g comporta-se como f .
(b) Outra razão é quando conhecemos f (x), mas a sua expressão é difı́cil de
trabalhar, por exemplo, de derivar ou integrar. Então se g coincide com f pelo
menos nos pontos xi, podemos fazer as contas com g ao invés de com f .
22
Um exemplo no qual podemos usar interpolação é o seguinte.
Exemplo 2.1. A tabela abaixo mostra a população (em milhões) do Brasil de 1960 a 2010.
Ano 1960 1970 1980 1990 2000 2010
População 72.7759 96.0604 121.7404 149.6483 174.5049 195.2102
Seja f (x) a população do Brasil no ano x. O problema a ser resolvido é o seguinte:
Estime a população nos anos 1955, 1975 e 2012.
DEFINIÇÃO 2.1. A interpolação denomina-se Interpolação Polinomial se a função g(x)
é um polinômio.
Note que um polinômio está definido em toda a reta real R e, além disso, é fácil de
derivar e integrar.
No que segue, a interpolação será sempre polinomial.
2.2 Preliminares
A pergunta natural então é: É possı́vel encontrar um polinômio que interpola f (x) nos
pontos x0, x1, ..., xn? ou seja, dados (n + 1)-pontos distintos x0, x1, ..., xn e suas respectivas
imagens yi = f (xi), i = 0, 1, ...,n, existirá um polinômio g(x) = P(x) tal que P(xi) = yi?
O Teorema abaixo fornece uma resposta à pergunta anterior.
TEOREMA 2.1. Seja f (x) definida em x0, x1, ..., xn (n + 1) pontos distintos de um intervalo
[a, b], então existe um único polinômio P(x) de grau menor ou igual a n tal que P(xi) = f (xi) = yi,
para todo i = 0, 1, ...,n.
Demonstração. Considere um polinômio de grau n da forma
P(x) = a0 + a1x + a2x2 + ... + anxn.
Sob a hipótese P(xi) = yi temos:
a0 + a1x0 + ... + anxn0 = y0
a0 + a1x1 + ... + anxn1 = y1
...
a0 + a1xn + ... + anxnn = yn
(2.1)
Isto é, se o sistema linear anterior tem solução então podemos encontrar um P(x) que
23
interpola f (x) tal que os coeficientes de P(x) são precisamente os valores da solução
do sistema. Logo, é suficiente provar que o sistema tem solução única.
De fato, a matriz do sistema é:
A =

1 x0 . . . xn0
1 x1 . . . xn1
...
...
. . .
...
1 xn . . . xnn

e, de álgebra linear, podemos provar que o determinante do sistema é:
Det(A) =
n∏
i> j
(xi − x j) (2.2)
Como xi , x j, se i , j, Det(A) , 0. Por tanto o sistema tem solução única. Dessa
forma o polinômio P(x) que interpola f (x) é único. �
DEFINIÇÃO 2.2. O polinômio P(x) do Teorema 2.1 é chamado Polinômio Interpolante
de f (x) nos pontos x0, x1, ..., xn.
Pelo Teorema 2.1, podemos obter os coeficientes do polinômio interpolante resol-
vendo o sistema linear (2.1). Mas, se n é grande, a resolução de um sistema linear
n× n pode ser muito custosa em relação ao tempo computacional utilizado. Assim,
é necessário utilizarmos outros métodos para determinar P(x) como a Fórmula de
Lagrande e a Fórmula de Newton.
2.3 A Fórmula de Lagrange
Para definir a fórmula de Lagrange do polinômio interpolante P(x) precisamos dos
polinômios de Lagrange que são definidos abaixo.
DEFINIÇÃO 2.3. Sejam x0, x1, ..., xn, (n + 1) pontos distintos. Seja k ∈ {0, ..,n}. O k-ésimo
polinômio de Lagrange Lk(x) é definido por:
Lk(x) =
∏
i,k
( x − xi
xk − xi
)
=
(x − x0)(x − x1) · · · (x − xk−1)(x − xk+1) · · · (x − xn)
(xk − x0)(xk − x1) · · · (xk − xk−1)(xk − xk+1) · · · (xk − xn)
Os polinômios de Lagrange tem a seguinte propriedade:
24
(A) Lk(xk) = 1 e Lk(xi) = 0, se i , k. De fato, substituindo x = xk, o numerador e o
denominados são iguais e assim Lk(xk) = 1. Por outro lado, substituindo x = xi
com i , k, vemos que o numerador anula-se no fator (x − xi), e assim Lk(xi) = 0.
DEFINIÇÃO 2.4. Sejam yi = f (xi). Então o polinômio de Lagrange de f (x) é dado por
P(x) = y1L1(x) + y2L2(x) + · · · + ynLn(x).
Pela propriedade (A), a avaliação de P(x) em x = xi é:
P(xi) = y1L1(xi) + y2L2(xi) + · · · + ynLn(xi) = yiLi(xi) = yi = f (xi),
para todo i ∈ {0, 1, ...,n}. Então P(x) interpola f (x), e pela unicidade do Teorema 2.1, P(x)
é o polinômio interpolante de f (x) nos pontos x0, x1, · · · , xn
Por tanto, a fórmula de Lagrange:
P(x) = y1L1(x) + y2L2(x) + · · · + ynLn(x).
fornece uma maneira direta de calcular o polinômio interpolante.
Exemplo 2.2. A seguinte tabela mostra a população (em milhões) mundial de 1970 a
2000.
Ano 1970 1980 1990 2000
População 3710 4450 5280 6080
Encontre um polinômio interpolador de grau 3 e estime a população nos anos 1985
e 1995.
Solução:
Dos dados, x0 = 1970, x1 = 1980, x2 = 1990 e x3 = 2000; y0 = 3710, y1 = 4450,
y2 = 5280 e y3 = 6080. Definimos os polinômios de Lagrange:
L0(x) =
(x − x1)(x − x2)(x − x3)
(x0 − x1)(x0 − x2)(x0 − x3)
=
(x − 1980)(x − 1990)(x − 2000)
(1970 − 1980)(1970 − 1990)(1970 −2000)
L1(x) =
(x − x0)(x − x2)(x − x3)
(x1 − x0)(x1 − x2)(x1 − x3)
=
(x − 1970)(x − 1990)(x − 2000)
(1980 − 1970)(1980 − 1990)(1980 − 2000)
L2(x) =
(x − x0)(x − x1)(x − x3)
(x2 − x0)(x2 − x1)(x2 − x3)
=
(x − 1970)(x − 1980)(x − 2000)
(1990 − 1970)(1990 − 1980)(1990 − 2000)
L3(x) =
(x − x0)(x − x1)(x − x2)
(x3 − x0)(x3 − x1)(x3 − x2)
=
(x − 1970)(x − 1980)(x − 1990)
(2000 − 1970)(2000 − 1980)(2000 − 1990)
25
O polinômio interpolante é:
P(x) = 3710L0(x) + 4450L1(x) + 5280L2(x) + 6080L3(x).
A população estimada no ano 1985 será
P(1985) = 3710L0(1985) + 4450L1(1985) + 5280L2(1985) + 6080L3(1985).
P(1985) = 3710(−0.0625) + 4450(0.5625) + 5280(0.5625) + 6080(−0.0625) = 4861.25.
Então a população estimada no ano 1985 é de 4861.25 milhões.
Dado real: Em 1985 a população mundial atingiu os 5000 milhões de pessoas.
2.4 A Fórmula de Newton
A Fórmula de Newton para determinar o polinômio interpolante de f (x) é uma fórmula
mais eficiente que a de Lagrange, pois precisa de menos custo computacional. A ideia
de Newton foi escrever o polinômio interpolante P(x) na forma
P(x) = b0 + b1(x − x0) + b2(x − x0)(x − x1) + · · · + bn(x − x0)(x − x1) · · · (x − xn−1), (2.3)
em que os coeficientes b0, b1, · · · , bn devem ser determinados. Está fórmula precisa da
definição das diferenças divididas.
Seja f uma função definida em (n + 1)-pontos x0, x1, · · · , xn nos quais yi = f (xi).
DEFINIÇÃO 2.5. As diferenças divididas de ordem 0, são os valores da própria função
f (x), ou seja,
f [xi] := f (xi), i = 0, 1, 2, ...,n
As diferenças divididas de ordem k, para 0 ≤ k ≤ n, são definidas utilizando as
diferenças de ordem k − 1:
f [xi, xi+1, · · · , xi+k] =
f [xi+1, · · · , xi+k] − f [xi, · · · , xi+k−1]
xi+k − xi
, para 0 ≤ i ≤ n − k
Por exemplo se n = 3, as diferenças divididas de ordem 0 de f para os valores
x0, x1, x2, x3 são:
f [x0] = f (x0), f [x1] = f (x1), f [x2] = f (x2), f [x3] = f (x3),
as de ordem 1:
f [x0, x1] =
f [x1] − f [x0]
x1 − x0
, f [x1, x2] =
f [x2] − f [x1]
x2 − x1
, f [x2, x3] =
f [x3] − f [x2]
x3 − x2
,
26
as de ordem 2:
f [x0, x1, x2] =
f [x1, x2] − f [x0, x1]
x2 − x0
, f [x1, x2, x3] =
f [x2, x3] − f [x1, x2]
x3 − x1
,
e as de ordem 3:
f [x0, x1, x2, x3] =
f [x1, x2, x3] − f [x0, x1, x2]
x3 − x0
.
Podemos organizar as diferenças divididas numa tabela, de modo que as diferenças
divididas de ordem 1 são calculadas a partir das diferenças de ordem 0, as de ordem 2 a
partir das de ordem 1, e assim sucessivamente.
xi yi ordem 0 ordem 1 ordem 2 ordem 3
x0 y0 f [x0]
f [x0, x1]
x1 y1 f [x1] f [x0, x1, x2]
f [x1, x2] f [x0, x1, x2, x3]
x2 y2 f [x2] f [x1, x2, x3]
f [x2, x3]
x3 y3 f [x3]
Observação 2.1. Uma propriedade fácil de verificar das diferenças divididas é que
estas não se alteram se permutamos dois valores xi e x j, por exemplo:
f [x1, x2] = f [x2, x1],
f [x1, x2, x3] = f [x2, x1, x3] = f [x2, x3, x1]
Como P(x) é o polinômio interpolante:
P(x0) = f (x0),P(x1) = f (x1), · · · ,P(xn) = f (xn).
Pela igualdade (2.3) temos que P(x0) = b0. Então
b0 = P(x0) = f (x0) = f [x0].
Assim foi provado que o coeficiente b0 deve ser igual à diferença dividida f [x0]. Da
mesma forma, da igualdade (2.3), P(x1) = b0 + b1(x1 − x0), ou, equivalentemente,
b1 =
P(x1) − b0
x1 − x0
=
f (x1) − f (x0)
x1 − x0
=
f [x1] − f [x0]
x1 − x0
= f [x0, x1].
Isto é, b1 é a primeira diferença de ordem 1. Continuando com a mesma ideia, como
f (x0) = P(x2) = b0 + b1(x2 − x0) + b2(x2 − x0)(x2 − x1), temos que
f [x2] = f [x0] + f [x0, x1](x2 − x0) + b2(x2 − x0)(x2 − x1)
27
Efetuando algumas operações
f [x2] − f [x0]
x2 − x0
= f [x0, x1] + b2(x2 − x1)
Ou
b2 =
f [x0, x2] − f [x0, x1]
x2 − x1
Pela Observação (2.1), temos
b2 =
f [x0, x2] − f [x1, x0]
x2 − x1
= f [x1, x0, x2] = f [x0, x1, x2]
Ou seja, b2 é a segunda diferença dividida.
Repetindo o processo para P(xi), obteremos que os coeficientes do Polinômio Inter-
polante na forma de Newton (2.3) são dados pelas diferenças divididas:
bi = f [x0, x1, · · · , xi], para i = 0, 1, 2, ...,n.
Observação 2.2. Resumindo, o polinômio interpolante na forma de Newton é:
P(x) = f [x0] + f [x0, x1](x − x0) + f [x0, x1, x2](x − x0)(x − x1) + · · ·
+ f [x0, · · · , xn](x − x0)(x − x1) · · · (x − xn−1).
Note que para formar o polinômio devemos tomar as diferenças divididas da parte
superior da tabela.
Exemplo 2.3. Um cabo sob ação do seu próprio peso está suspenso entre dois pontos
distantes 24 metros.
12 metros
Figura 10: Um cabo sob ação do seu próprio peso.
A parte mais baixa do cabo ficou a 12 metros de altura. Foram medidas diferentes
alturas em vários pontos as quais estão mostradas na tabela abaixo
28
xi yi xi yi
-2 12.16 2 12.16
-7 14.10 7 14.10
-9 15.53 9 15.53
-12 18.51 12 18.51
Estime a altura a uma distância de 5 metros do centro.
Solução:
Pela simetria do problema em relação ao eixo y, é suficiente considerar os dados
x0 = 0, x1 = 2, x2 = 7, x3 = 9 e x4 = 12.
xi yi ordem 0 ordem 1 ordem 2 ordem 3 ordem 4
0 12.00 12.00
0.08
2 12.16 12.16 0.044
0.388 0.0003
7 14.10 14.10 0.0467 0.000049
0.715 0.00089
9 15.53 15.53 0.0556
0.993
12 18.51 18.51
Pela fórmula de Newton
P(x) = 12.00 + 0.08(x − 0) + 0.044(x − 0)(x − 2) + 0.0003(x − 0)(x − 2)(x − 7)+ (2.4)
+0.000049(x − 0)(x − 2)(x − 7)(x − 9).
Substituindo o valor de x = 5 obtemos: P(5) = 13.0618.
29
2.5 Exercı́cios
Exercı́cio 2.1. Use uma cúbica para determinar uma aproximação para a única raiz
positiva da equação 4 cos(x) − ex = 0.
Exercı́cio 2.2. Considere a função f (x) definida nos pontos, conforme tabela. Determine
o polinômio interpolador, usando a Fórmula de Lagrange, e estime f (0.8).
xi 0 0.5 1.0
f (xi) 1.3 2.5 0.9
Exercı́cio 2.3. Considere a função f (x) = 3+x1+x definida nos pontos conforme tabela.
Determine o polinômio interpolador, usando a Fórmula de Lagrange, e estime f (0.25).
xi 0.1 0.2 0.4
f (xi) 2.82 2.67 2.43
Exercı́cio 2.4. Considere a tabela:
xi 1 3 4 5
f (xi) 0 6 24 60
(a) Determine o polinômio, na forma de Lagrange, sobre todos os pontos.
(b) Calcule f (3.5).
Exercı́cio 2.5. Construir o polinômio de interpolação, na forma de Lagrange, para a função
y = sin(πx), escolhendo os pontos: x0 = 0, x1 = 16 e x2 =
1
2 .
Exercı́cio 2.6. A integral elı́ptica é definida por: K(k) =
∫ π
2
0
dx
(1−k2 sin2 x)1/2
. Por uma tabela
de valores desta integral, encontramos:
K(1) = 1.5708; K(2) = 1.5719; K(3) = 1.5739.
Determine K(2.5), usando o polinômio de interpolação, na forma de Lagrange, sobre todos
os pontos.
Exercı́cio 2.7. Dada a função tabelada f (x) definida pelos pontos
f (0) = 0, f (1) = 0.5, f (1.5) = 0.4, f (2.5) = 0.285, f (3.0) = 0.25.
(a) Determinar o polinômio de interpolação usando a Fórmula de Newton sobre dois
pontos (interpolação linear).
30
(b) Determinar o polinômio de interpolação usando a Fórmula de Newton sobre três
pontos (interpolação quadrática).
(c) Calcular f (0.5) usando os items a) e b).
Exercı́cio 2.8. Seja f (x) uma função definida pelos pontos f (−1) = 0, f (0) = 1 e f (2) = −1.
Usando a forma de Newton, encontre o polinômio p(x) que interpola a função f .
Exercı́cio 2.9. Considere a função f (x) = ln(x) para a qual temos
f (1) = 0, f (2) = 0.6931, f (3) = 1.0986, f (4) = 1.3863.
Usando o Método de Newton encontre uma aproximação de ln(3.7).
Exercı́cio 2.10. Seja f (x) dadas pelos seguintes valores
f (0.2) = 0.16, f (0.34) = 0.22, f (0.4) = 0.27, f (0.52) = 0.29, f (0.6) = 0.32, f (0.72) = 0.37.
Obter f (0.47) usando o polinômio de grau 2.
Exercı́cio 2.11. Considere a função f (x) =
√
x definida nos pontos conforme tabela.
Determinar o valor aproximado de
√
1.12 usando o polinômio de Interpolação de Newton
sobre três pontos.
xi 1.0 1.1 1.15 1.25 1.3
f (xi) 1 1.048 1.072 1.118 1.14
Exercı́cio 2.12. Sabendo-se que a equação x4 + 6x2 − 1 = 0 tem uma raiz em [0, 1],
determinar o valor aproximado dessa raiz usando polinômio de Interpolaçãode Newton
sobre três pontos.
Exercı́cio 2.13. A função y =
∫
∞
x
e−t
t dt é dada pela tabela:
xi 0.01 0.02 0.03 0.04 0.05 0.06
f (xi) 4.0379 3.3547 2.9591 2.6813 2.4679 2.2953
Através da Fórmula de Newton, calcule y para x = 0.0378 usando parábola e uma
cúbica.
Exercı́cio 2.14. Dada a tabela abaixo, calcule e3.1 usando um polinômio de interpolação
sobre três pontos:
xi 2.4 2.6 2.8 3.0 3.2 3.4 3.6 3.8
exi 11.02 13.46 16.44 20.08 24.53 29.96 36.59 44.70
31
Exercı́cio 2.15. Construa a tabela de diferenças divididas para a função f (x) com os
dados:
f (0) = 0, f (0.5) = −2.241, f (1.0) = −1.65, f (1.5) = −0.594, f (2) = 1.34, f (2.5) = 4.564
Estime o valor de f (1.23).
Exercı́cio 2.16. A raiz de uma função contı́nua pode ser aproximada pela raiz de seu
polinômio. Usando um polinômio de grau dois, encontre x̄ tal que f (x̄) = 0. Use 5 casas
decimais. Os valores exatos de f , para os valores dados de x, são:
f (1) = −0.757, f (2) = 0.141, f (3) = 0.842, f (4) = 0.909.
Exercı́cio 2.17. Uma maneira de se calcular a derivada de uma função em um ponto
x0, quando não se conhece a expressão da mesma, é usar uma tabela para formar um
polinômio que aproxime a função, derivar então esse polinômio e avaliar sua derivada
em x = x0. Dados os pontos abaixo, calcule f ′(0.50) usando um polinômio interpolador de
grau 2:
f (0.40) = 1.51, f (0.45) = 1.49, f (0.50) = 1.47, f (0.55) = 1.44, f (0.60) = 1.42.
Exercı́cio 2.18. Engenheiros e programadores de computadores usam interpolação po-
linomial para fazer gráficos por computador. Por exemplo, suponha que deseja-se
construir uma curva que passa pelos pontos A = (2, 3), B = (4, 2), C = (5, 0) e D = (3,−2)
como na figura. Para auxiliar estes professionais, usando Interpolação de Lagrange, en-
contre uma curva que passa por estes 4 pontos. Faça as hipótesis que ache necessário
para a solução do problema.
Figura 11: Desenho no computador.
Exercı́cio 2.19. A seguinte tabela mostra a população do Brasil de 1960 a 2010. Encontre
um polinômio interpolador de grau 5 e estime a população nos anos 1955, 1975 e 2012.
32
A população em 2012 foi aproximadamente 198.656 milhões. Quão precisas você acha
que são as aproximações nos anos 1955 e 1975?
Ano 1960 1970 1980 1990 2000 2010
População( em milhões) 72.7759 96.0604 121.7404 149.6483 174.5049 195.2102
33
3 Resolução Numérica de Sistemas Lineares
3.1 Introdução
Nesta seção estudaremos técnicas para resolver sistemas lineares:
a11x1 + a12x2 + ... + a1nxn = b1
a21x1 + a22x2 + ... + a2nxn = b2
...
an1x1 + an2xn2 + ... + annxn = bn
(3.1)
Estes sistemas podem ser escritos na forma
Ax = b (3.2)
em que A = (ai j) é uma matriz n × n e b é um vetor coluna n × 1. Por álgebra linear
sabemos que se Det(A) , 0 então o sistema tem solução única. Além disso, podemos
resolvê-lo usando a redução de Gauss (ou eliminação Gaussiana), isto é, convertendo o
sistema a um triangular superior.
Aqui estudaremos dois métodos com os quais podemos resolver um ou mais sistemas
de equações lineares: O Método da Decomposição LU e o Método de Gauss-Seidel.
3.2 Preliminares
Existem dois tipos de sistemas lineares que são fáceis de resolver: Os sistemas triangu-
lares Inferiores e os sistemas triangulares Superiores.
Usaremos a letra L (de lower=inferior em inglês) para uma matriz triangular inferior
e U (de upper=superior em inglês) para uma matriz triangular superior.
3.2.1 Sistemas Triangulares Inferiores
São aqueles sistemas Lx = b na qual todos os termos que estão acima da diagonal de
L = (li j) são 0:
l11x1 = b1
l21x1 +l22x2 = b2
...
...
. . . =
...
ln1x1 +ln2x2 · · · +lnnxn = bn
(3.3)
34
Exemplo 3.1. Resolva o sistema
3x1 = 1
x1 +5x2 = 1
2x1 −x2 +2x3 = 3
(3.4)
Solução:
Observe que da primeira equação x1 = 13 . Substituindo na segunda equação
x2 = (1 − x1)/5 = (1 −
1
3
)/5 =
2
15
.
Agora, na terceira equação
x3 = (3 − 2x1 + x2)/2 = (4 −
2
3
+
2
15
)/2 =
26
15
.
Do Exemplo 3.1, podemos fazer a seguinte observação.
Observação 3.1. A solução geral de um sistema triangular inferior como (3.3), é dada
pelas fórmulas:
x1 =
b1
l11
, e xi =
bi −
∑i−1
j=1 li jx j
lii
, para todo i = 2, ...,n
3.2.2 Sistemas Triangulares Superiores
São aqueles sistemas Ux = b na qual todos os termos que estão abaixo da diagonal de
U = (ui j) são 0:
u11x1 +u12x2 . . . +u1nxn = b1
+u22x2 . . . +u2nxn = b2
. . .
... =
...
unnxn = bn
(3.5)
Exemplo 3.2. Resolva o sistema
2x1 +3x2 −3x3 = 4
4x2 −x3 = −3
2x3 = 3
(3.6)
35
Solução:
Observe que da terceira equação x3 = 32 . Substituindo na segunda equação
x2 = (−3 + x3)/4 = (−3 +
3
2
)/4 = −
3
8
.
Agora, isolando x1 na primeira equação
x1 = (4 + 3x3 − 3x2)/2 = (4 +
9
2
+
9
8
)/2 =
59
16
.
Do Exemplo 3.2, podemos fazer a seguinte observação.
Observação 3.2. A solução geral de um sistema triangular superior como (3.5), é dada
pelas fórmulas:
xn =
bn
unn
, e xi =
bi −
∑n
j=i+1 ui jx j
uii
, para todo i = 1, ...,n − 1. (3.7)
3.3 Método Direto: Decomposição LU
Quando resolvemos vários sistemas lineares
Ax = b1,Ax = b2, ...,Ax = bp,
por eliminação Gaussiana, em que todos os sistemas tem a mesma matriz A, teremos que
repetir o método da eliminação para cada vetor coluna b. Isto traz como consequência o
uso de uma grande quantidade de memória para armazenar os dados de cada sistema.
Uma maneira de evitar esse consumo de memória é utilizando o Método da De-
composição LU, que é mais eficiente pois permite calcular a decomposição sem usar os
vetores coluna b.
Antes de descrever o método, apresentaremos um exemplo.
Exemplo 3.3. Resolva o sistema 2 0
−3 3
  3 −40 −6
  x1x2
 =  5
−3
 (3.8)
36
Solução:
Repare que o sistema pode-se escrever na forma Ax = b, em que A = LU, sendo L
uma matriz triangular inferior e U uma matriz triangular superior.
L =
 2 0
−3 3
 , U =  3 −40 −6
 , b =  5
−3

Fazendo
 y1y2
 := Ux, concluimos que encontraremos a solução do sistema se
resolvemos os dois sistemas: 2 0
−3 3
  y1y2
 =  5
−3
 , e  3 −40 −6
  x1x2
 =  y1y2

O primeiro sistema tem solução y1 = 53 e y2 = (−3 + 3y1)/3 = (−3 + 5)/3 =
2
3 . Tendo
estes valores, podemos encontrar os valores de x1 e x2 no segundo sistema. Assim,
−6x2 = y2 = 23 então x2 = −
1
9 , e x1 = (y1 + 4x2)/3 = (
5
3 −
4
9 )/3 =
11
27 .
Observação 3.3. O Exemplo (3.3) mostra que se por acaso podemos decompor a
matriz A como produto de duas matrizes triangulares A = LU (em que L triangular
inferior e U triangular superior), então a resolução do sistema Ax = b se reduz a
resolução de dois sistemas: Ly = b que é triangular inferior, eUx = y que é triangular superior.
Isto será o objetivo da decomposição LU.
Assim, a pergunta natural é: Quando uma matriz A pode ser decomposta como produto
de duas matrizes LU, em que L é triangular inferior e U é triangular superior? ou, toda
matriz A tem decomposição LU?
Para responder estas questões, precisamos da seguinte definição.
DEFINIÇÃO 3.1. Denomina-se menores principais de ordem k de uma matriz A = (ai j),
com i, j = 1, · · · ,n a:
∆k = Det(Ak)
em que Ak = (ai j), com i, j = 1, · · · , k, é formada pelas k primeiras linhas e k primeiras
colunas de A.
37
Exemplo 3.4. Encontre os menores principais da matriz:
A =

1 −2 3
−2 3 −5
−2 2 4

Solução:
Pela definição ∆1 = Det(a11) = 1, ∆2 = Det
 1 −2
−2 3
 = −1.
O Teorema abaixo fornece condições para a existência de uma decomposição LU.
TEOREMA 3.1. Considere A = (ai, j), i, j = 1, · · · ,n. Se os menores principais ∆k são distintos
de 0, para k = 1, · · · ,n − 1, então A se decompõe, de maneira única, no produto LU em que
L = (li j) é uma matriz triangular inferior com lii = 1, i = 1, · ·· ,n e U = (ui j) é uma matriz
triangular superior.
Demonstração. Provaremos o Teorema por indução. Se n = 1, é claro que A =
(a11) pode ser decomposta, escrevendo A = LU = (1)(aii).
Suponha que, se A é de ordem k − 1, então pode ser decomposta na forma
Ak−1 = Lk−1Uk−1
em que Lk−1 é triangular inferior e Uk−1 é triangular superior.
Devemos provar que existe uma decomposição para matrizes de ordem k. Seja
A de ordem k, então escrevemos
A =
 Ak−1 sr akk

em que r é um vetor linha e s é um vetor coluna. Dado que Ak−1 é uma matriz de
ordem k− 1 que, pela hipótese, pode ser decomposta na forma Ak−1 = Lk−1Uk−1.
Defina as matrizes
L =
 Lk−1 0rU−1k−1 1
 , e U =  Uk−1 L−1k−1s0 akk − rU−1k−1L−1k−1s
 (3.9)
Note que L é triangular inferior e U é triangular superior. Além disso:
LU =
 Lk−1 0rU−1k−1 1
  Uk−1 L−1k−1s0 akk − rU−1k−1L−1k−1s
 =  Lk−1Uk−1 sr akk
 = A
Por tanto A se decompõe como um produto LU em que L e U são dadas pelas
fórmulas (3.9). �
38
Observação 3.4. Da prova do Teorema (3.1) podemos concluir que:
(a) Se a matriz é de ordem n×n, é suficiente verificar que os menores determinantes
são distintos do 0 até a ordem n − 1.
(b) A condição que usaremos para a matriz L é que lii = 1,∀i = 1, · · · ,n.
Repare que a hipótese lii = 1,∀i = 1, ...,n, simplifica os cálculos para encontrar os
coeficientes das matrizes L e U. Quando usamos está hipótese, o método é chamado
Método de Doolittle. Mas podia-se considerar outras hipóteses. Por exemplo, outros
dois métodos de decomposição LU são:
(i) Método de Crout : Quando uii = 1,∀i = 1, ...,n, e
(ii) Método de Choleski: Quando uii = lii,∀i = 1, ...,n.
Exemplo 3.5. Resolva o sistema linear
1 −3 1
2 −8 8
−6 3 −15


x1
x2
x3
 =

4
−2
9

Solução:
Os menores principais são ∆1 = 1 e ∆2 = −2. Como ambos são distintos do 0, existe
a decomposição LU.
Sejam
L =

1 0 0
l21 1 0
l31 l32 1
 ,U =

u11 u12 u13
0 u22 u23
0 0 u33

Da equação A = LU, obtemos
1 −3 1
2 −8 8
−6 3 −15
 =

1 0 0
l21 1 0
l31 l32 1


u11 u12 u13
0 u22 u23
0 0 u33

Encontraremos a primeira linha de U multiplicando a primeira linha de L por U:
1.u11 = 1, 1.u12 = −3, e 1.u13 = 1, então u11 = 1,u12 = −3,u13 = 1.
Também, encontramos a primeira coluna de L multiplicando L pela primeira coluna
39
de U:
l21u11 = 2, l31u11 = −6 então l21 = 2, l31 = −6.
Trabalharemos agora com a segunda linha de U. Esta é determinada multiplicando
a segunda linha de L por U:
l21u12 + u22 = −8, l21u13 + u23 = 8, então u22 = −2,u23 = 6.
A segunda coluna de L é encontrada multiplicando L pela segunda coluna de U:
l31u12 + l32u22 = 3, então l32 =
15
2
Por último, determinamos a terceira linha de U multiplicando a terceira linha de L:
l31u13 + l32u23 + u33 = −15, então u33 = −54.
Fazendo y =

y1
y2
y3
 := Ux, obtemos dois sistemas lineares que devem ser resolvidos:
Ly = b, e Ux = y, ou
1 0 0
2 1 0
−6 152 1


y1
y2
y3
 =

4
−2
9
 e

1 −3 1
0 −2 6
0 0 −54


x1
x2
x3
 =

y1
y2
y3

Do primeiro sistema y1 = 4, y2 = (−2 − 2y1) = −10, e y3 = (9 + 6y1 −
15y2
2 ) = 108. Por
fim, com os valores dos yi, encontramos a solução:
x3 = −
108
54
= −2, x2 = (y2 − 6x3)/(−2) = −1, e x1 = (y1 + 3x2 − x3) = 3.
Observação 3.5. Podemos generalizar as fórmulas do exemplo anterior:
ui j = ai j −
i−1∑
k=1
likukj, para i, j = 1, ...,n, i ≤ j.
começando com i = 1, e
li j =
ai j −
∑ j−1
k=1 likukj
u j j
, para i, j = 1, ...,n, i ≥ j.
começando com j = 1.
40
3.4 Método Iterativo: O método de Gauss-Seidel
Assim como nos métodos de resolução de equações não lineares usavamos métodos
iterativos para encontrar uma solução aproximada, nesta seção apresentaremos um
método iterativo para resolver sistemas lineares. Este método é chamado Método de Gauss-
Seidel e permite determinar uma solução aproximada de um sistema linear.
A ideia central por traz dos métodos iterativos é escrever o problema na forma
x = φ(x),
e então determinar aproximações da solução utilizando a sequência x(k+1) = φ(x(k)), a
partir de uma aproximação inicial x(0). Se a sequência {x(k)} é convergente, então conver-
girá para a solução procurada. Aqui seguiremos esta mesma ideia aplicada a sistemas
lineares.
Apresentaremos o método com um exemplo.
Exemplo 3.6. Encontre uma solução aproximada do sistema linear:
6x1 −x2 −x3 = 3
6x1 +9x2 +x3 = 40
−3x1 +x2 +12x3 = 50
Solução:
Seja x = (x1, x2, x3)t. Note que o sistema se pode escrever
x1 −
x2
6 −
x3
6 =
1
2
2x1
3 +x2 +
x3
9 =
40
9
−
x1
4 +
x2
12 +x3 =
25
6
no qual os coeficientes da diagonal são todos iguais a 1. Dessa forma, isolando as
variáveis da diagonal
x1 = 12 +
x2
6 +
x3
6
x2 = −
2x1
3 +
40
9 −
x3
9
x3 =
x1
4 −
x2
12 +
25
6
,
Repare que este sistema tem a forma x = φ(x). Assim, podemos utilizar o método
iterativo
x(k+1)1 =
1
2 +
x(k)2
6 +
x(k)3
6
x(k+1)2 = −
2x(k)1
3 +
40
9 −
x(k)3
9
x(k+1)3 =
x(k)1
4 −
x(k)2
12 +
25
6
, (3.10)
em que estamos utilizando os ı́ndices na parte superior. Suponha que a aproximação
inicial, para k = 0, é x(0) = (x(0)0 , x
(0)
1 , x
(0)
3 )
t = (0, 0, 0)t. Cada aproximação é um vetor
41
coluna. Da primeira equação:
x(1)1 =
1
2
+
x(0)2
6
+
x(0)3
6
=
1
2
+
0
6
+
0
6
=
1
2
.
O ponto central: se já temos uma aproximação x(1)1 de x1, podemos utilizá-la para
determinar a aproximação x(1)2 de x2, substituindo x
(0)
1 por x
(1)
1 na segunda equação
de (3.10):
x(1)2 = −
2x(1)1
3
+
40
9
−
x(0)3
9
= −
2(1/2)
3
+
40
9
−
0
9
=
37
9
= 4.11111
Seguindo esta ideia, podemos substituir x(0)1 por x
(1)
1 e x
(0)
2 por x
(1)
2 na terceira equação
de (3.10):
x(1)3 =
x(1)1
4
−
x(1)2
12
+
25
6
=
(1/2)
4
−
(37/9)
12
+
25
6
= 3.94907.
Logo, o processo iterativo, na iteração k, pode ser modificado ao seguinte:
x(k+1)1 =
1
2
+
x(k)2
6
+
x(k)3
6
x(k+1)2 = −
2x(k+1)1
3
+
40
9
−
x(k)3
9
x(k+1)3 =
x(k+1)1
4
−
x(k+1)2
12
+
25
6
Se k = 1 então
x(2)1 =
1
2
+
x(1)2
6
+
x(1)3
6
=
1
2
+
4.11111
6
+
3.94907
6
= 1.84336
x(2)2 = −
2x(2)1
3
+
40
9
−
x(1)3
9
= −
2(1.84336)
3
+
40
9
−
3.94907
9
= 2.77675
x(2)3 =
x(2)1
4
−
x(2)2
12
+
25
6
=
1.84336
4
−
2.77675
12
+
25
6
= 4.39611.
Para k = 2 obtemos x(3)1 = 1.69547, x
(3)
2 = 2.82567 e x
(3)
3 = 4.35506.
42
3.4.1 O Caso Geral
Tendo como base o Exemplo (3.6), escreveremos as fórmulas para o caso geral. Considere
um sistema de equações lineares:
a11x1 + a12x2 + ... + a1nxn = b1
a21x1 + a22x2 + ... + a2nxn = b2
...
an1x1 + an2xn2 + ... + annxn = bn
(3.11)
que pode ser escrito na forma Ax = b, em que A = (ai j), b = (bi). Suponha aii , 0 para
todo i = 1, ...,n, então podemos dividir cada equação pelo correspondente elemento da
diagonal. Então o sistema converte-se a:
x1 + a∗12x2 + ... + a
∗
1nxn = b
∗
1
a∗21x1 + x2 + ... + a
∗
2nxn = b
∗
2
...
a∗n1x1 + a
∗
n2xn2 + ... + xn = b
∗
n
(3.12)
em que a∗i j =
ai j
aii
para todo i, j = 1, ...,n e b∗i =
bi
aii
, para todo i = 1, ...,n. O sistema
anterior pode-se escrever na forma A∗x = b∗, em que A∗ = (a∗i j). Esta matriz A
∗ pode ser
decomposta na forma
A∗ = L∗ + I + R∗ (3.13)
em que L∗ = (l∗i j) é uma matriz triangular inferior com os elementos da diagonal iguais a
0, I é a matriz identidade e R∗ = (r∗i j) é uma matriz triangular superior com os elementos
da diagonal iguais a 0. Ou seja,
l∗i j =

ai j
aii
sei > j
0 se i ≤ j
e r∗i j =

ai j
aii
sei < j
0 se i ≥ j
O Sistema Ax = b é agora dado por
(L∗ + I + R∗)x = b∗
que é equivalente à equação (L∗ + I)x = −R∗x + b∗.
DEFINIÇÃO 3.2. Ométodo de Gauss-Seidel é dado pelo processo iterativo:
(L∗ + I)x(k+1) = −R∗x(k) + b∗
com condição inicial x(0) = (x(0)1 , x
(0)
2 , x
(0)
3 )
t. Cada aproximação x(k) é um vetor coluna.
43
Acostuma-se escrevê-lo na forma
x(k+1) = −L∗x(k+1) − R∗x(k) + b∗
para expressar o fato que podemos aproveitar as aproximações x(k+1)1 , x
(k+1)
2 , ..., x
(k+1)
i−1
para calcular a seguinte aproximação x(k+1)i , para i = 2, ...,n.
3.4.2 Critérios de Convergência
Só falta fornecer condições sobre as quais a sequência x(k) converge para a solução do
sistema linear Ax = b. Dois critérios são os mais utilizados.
O método de Gauss-Seidel converge se um dos dois critérios for satisfeito:
(a) O critério das linhas, isto é, se
max
1≤i≤n
{βi} < 1
em que
βi =
∑
j,i
|a∗i j| =
i−1∑
j=1
|a∗i j| +
n∑
j=i+1
|a∗i j|,∀i = 1, ...,n
é a soma de todos os termos da linha i da matriz A∗, exceto o termo da diagonal.
(b) O critério de Sassenfeld, isto é, se
max
1≤i≤n
{βi} < 1
em que
βi =
i−1∑
j=1
|a∗i j|β j +
n∑
j=i+1
|a∗i j|,∀i = 1, ...,n.
3.4.3 Critério de Parada
Precisamos de um critério para decidir quando parar de realizar as iterações no método
de Gauss-Seidel. Ao igual que na resolução de equações não lineares, temos dois critérios
de parada.
44
O método de Gauss-Seidel para se atingimos os seguintes erros
• Erro Absoluto. Se
max
1≤i≤n
|x(k+1)i − x
(k)
i | < �. (3.14)
ou seja, se cada componente atinge um erro absoluto menor que �.
• Erro Relativo. Se
max
1≤i≤n
|x(k+1)i − x
(k)
i |
|x(k+1)i |
< �. (3.15)
ou seja, se cada componente atinge um erro relativo menor que �.
Exemplo 3.7. Aplique o método de Gauss-Seidel para resolver:
9x1 −3x2 +x3 = 2
−2x1 +5x2 +3x3 = −1
−x1 2x2 −7x3 = 3
Considerando um erro absoluto de � = 0.04.
Solução:
Dividindo entre os coeficientes da diagonal obtemos:
x1 − 13 x2 +
1
9 x3 =
2
9
−
2
5 x1 +x2 +
3
5 x3 = −
1
5
1
7 x1 −
2
7 x2 +x3 = −
3
7
Repare que o método das linhas não é satisfeito, pois os coeficientes são:
β1 =
1
3
+
1
9
=
4
9
, β2 =
2
5
+
3
5
= 1 e β3 =
1
7
+
2
7
=
3
7
e então max1≤i≤n βi = 1. Assim, o critério das linhas não fornece informação sob a
convergência.
Mas, o critério de Sassenfeld é satisfeito pois:
β1 =
1
3
+
1
9
=
4
9
< 1,
β2 = (
2
5
)(
4
9
) +
3
5
=
7
9
< 1,
β3 = (
1
7
)(
4
9
) + (
2
7
)(
7
9
) =
2
7
< 1.
45
Dessa forma, o método de Gauss-Seidel, dado pelo processo iterativo:
x(k+1)1 =
2
9
+
x(k)2
3
−
x(k)3
9
x(k+1)2 = +
2x(k+1)1
5
−
1
5
−
3x(k)3
5
x(k+1)3 = −
x(k+1)1
7
+
2x(k+1)2
7
−
3
7
com condição inicial x(0) = (0, 0, 0), é convergente. As aproximações estão tabeladas
abaixo.
k x(k)1 x
(k)
2 x
(k)
3 erro 1 erro 2 erro 3
0 0 0 0
1 0.22222 -0.11111 -0.49206
2 0.23985 0.19117 -0.40821
3 0.33130 0.17744 -0.35760 0.09145 0.01373 0.05061
4 0.32110 0.143 -0.37910 0.0102 0.03444 0.0215
Por tanto, uma solução aproximada é (0.32110, 0.143,−0.37910).
46
3.5 Exercı́cios
Exercı́cio 3.1. Para cada um dos sistemas lineares seguintes, obtenha uma solução por
meio de um gráfico, se possı́vel for. Explique os resultados do ponto de vista geométrico.
(a)
x1 + 2x2 = 3
x1 − x2 = 0
(b)
x1 + 2x2 = 3
2x1 + 4x2 = 6
(c)
x1 + 2x2 = 3
2x1 + 4x2 = −6
(d)
2x1 + x2 + x3 = 1
2x1 + 4x2 − x3 = −1
Exercı́cio 3.2. Dê a Fatoração LU de cada matriz e resolva o sistema:
(a)
x1 + x2 + x3 = −1
3x1 + 2x2 − 2x3 = 2
4x1 + 3x2 − 4x3 = 7
(b)
2x1 + 3x2 + x3 − x4 = 6, 9
−x1 + x2 − 4x3 + x4 = −6, 6
x1 + x2 + x3 + x4 = 10, 2
4x1 − 5x2 + x3 − 2x4 = −12, 3
(c)
x + y + 2z + 4t = 7, 12
2x + 5y + z + 2t = 14, 90
x + y + 5z + 6t = 12, 02
4x + 6y + 2z + t = 20, 72
Exercı́cio 3.3. Utilize a Fatoração LU para resolver o sistema linear a seguir:
x1 + x2 + 2x3 + 4x4 = 7, 12
2x1 + 5x2 + x3 + 2x4 = 14, 90
x1 + x2 + 5x3 + 6x4 = 12, 02
4x1 + 6x2 + 2x3 + x4 = 20, 72
Exercı́cio 3.4. Considere o sistema linear
5x1 + 2x2 + x3 = 7
−x1 + 4x2 + 2x3 = 3
2x1 − 3x2 + 10x3 = −1
(a) Verificar a possibilidade de aplicação do método de Gauss-Siedel, usando o Critério
de Sassenfeld.
(b) Se possı́vel, resolvê-lo pelo método do item (a), obtendo o resultado com erro
absoluto < 10−2.
Exercı́cio 3.5. Considere os seguintes sistemas de equações lineares
(I)
−9x1 + 5x2 + 6x3 = 11
2x1 + 3x2 + x3 = 4
−x1 + x2 − 3x3 = −2
(II)
5x1 + 2x2 + x3 = −12
−x1 + 4x2 + 2x3 = 20
2x1 − 3x2 + 10x3 = 3
(a) Resolva os sistemas dados usando o método de Decomposição LU.
(b) Caso haja convergência garantida, resolva os sistemas dados pelo Método Iterativo
de Gauss-Siedel com (x(0)1 , x
(0)
2 , x
(0)
3 ) = (0.1, 0.2, 0.5) e � = 0.01.
47
Exercı́cio 3.6. Dados os sistemas lineares:
(I)
4x1 + 2x2 + 6x3 = 1
4x1 + x2 + 3x3 = 2
−x1 + 5x2 + 3x3 = 3
(II)
2x1 + x2 + 3x3 = 9
−x2 + x3 = 1
x1 + 3x3 = 3
(III)
2x2 + 2x3 = 8
x1 + 3x2 = 6
3x1 + x2 + x3 = 4
mostrar que, reordenando as equações e incógnitas, podemos fazer com que o critério
de Sassenfeld seja satisfeito, mas não o critério das linhas. Resolva-os.
48
4 Integração Numérica
4.1 Introdução
Cada vez que descrevemos um problema utilizando derivadas podemos inferir que a
solução do mesmo envolve integrais do tipo
I( f ) =
∫ b
a
f (x)dx (4.1)
em que f : [a, b]→ R é uma função contı́nua no intervalo [a, b] com valores nos números
reais. Estudando as técnicas de integração do Cálculo, sabemos que existem funções f (x)
para as quais não existe uma integral em termos de funções elementares. Por exemplo,
as seguintes integrais são difı́ceis de calcular∫
sin(x)
x
dx,
∫
e−x
2
dx.
Assim, é necessario o uso de integração numérica.
Em outros casos, quando temos pouca informação sobre a função, por exemplo,
quando conhecemos a função em um número finito de pontos, também é necessária
alguma técnica de integração numérica. Por exemplo no seguinte problema:
Exemplo 4.1. Um radar foi usado para medir a velocidade de um corredor durante os
primeiros 5 segundos de uma corrida (Veja a tabela ). Estime a distância que o corredor
cobriu durante aqueles 5 segundos.
t(s) 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5
v(m/s) 0 4.67 7.34 8.86 9.73 10.22 10.51 10.67 10.76 10.81 10.81
Note que só conhecemos a função v(t) em 10 pontos, mas mesmo assim, desejamos
determinar o valor aproximado da integral∫ 5
0
v(t)dt,
dado que a distância é a integral da velocidade em relação ao tempo.
Numa terceira situação, pode ser que exista uma fórmula para integrar f (x) mas tal
vez está fórmula não seja a maneira mais eficiente de calcular I( f ).
49
4.2 Preliminares
Considere uma função contı́nua f : [a, b]→ R definida num intervalo [a, b].
A ideia básica da integração numérica é obter I( f ) via interpolação. Isto é, primeiro
aproximar f (x) por um polinômio interpolante Pn(x) de grau n em [a, b], um polinômio
sobre n + 1 pontos x0 = a, x1, ..., xn−1, xn = b tais que a separação entre eles é dada por
h =
b − a
n
,
ou xi = x0 + ih, para todo i = 0, 1, ...,n.
Figura 12: Polinômio Pn(x) que interpola f (x).
Então pode-se aproximar a integral de f (x) usando a integral de Pn(x).
Se Pn(x) interpola f (x) nos pontos x0 = a, x1, ..., xn = b, então
I( f ) =
∫ b
a
f (x)dx ≈
∫ b
a
Pn(x)dx.
Assim, aproximar I( f ) reduz-se a determinar a integral de um polinômio, o qual é fácil
de realizar.
4.2.1 Estudo do Erro
Uma vez feita a aproximação, fica ainda o estudo do erro cometido no processo de
aproximar. Este erro pode ser estudado com a equação:
f (x) = Pn(x) + Rn(x) (4.2)
obtida na interpolação de Newton, em que Pn(x) é o polinômio interpolante e Rn(x) é o
resto dado por
Rn(x) = (x − x0)(x − x1)...(x − xn) f [x0, x1, ..., xn, x]
50
em que f [x0, x1, ..., xn, x] é a diferença dividida de ordem (n + 1).
Integrando a igualdade (4.2) entre a e b:
I( f ) =
∫ b
a
Pn(x)dx +
∫ b
a
Rn(x),
a partir do qual deduzimos que o erro cometido En é dado por
En =
∫ b
a
Rn(x).
Repare que pode-se escrever o resto na forma Rn(x) = ψ(x) f [x0, x1, ...,xn, x] em que
ψ(x) = (x − x0)(x − x1)...(x − xn)
Usando o Teorema de Rolle, pode-se mostrar que existe um ξ ∈ [a, b] tal que
f [x0, x1, ..., xn, x] =
f (n+1)(ξ)
(n + 1)!
.
Assim, o erro é dado pela integral
En =
∫ xn
x0
f (n+1)(ξ)
(n + 1)!
ψ(x)dx. (4.3)
O primeiro resultado é o seguinte:
PROPOSIÇÃO 4.1. O erro En é limitado pela seguinte expressão:
|En| ≤
K
(n + 1)!
∫ xn
x0
|ψ(x)|dx, (4.4)
em que
K = max
{
| f (n+1)(x)| : x0 ≤ x ≤ xn
}
.
Demonstração. Pela Equação (4.3), obtemos
|En| ≤ max
{
|
f (n+1)(ξ)
(n + 1)!
|, ξ ∈ [x0, xn]
} ∫ xn
x0
|ψ(x)|dx ≤
K
(n + 1)!
∫ xn
x0
|ψ(x)|dx.
pela definição de K. �
51
Observação 4.1. Pela Proposição 4.1, é suficiente estudar o comportamento de∫ xn
x0
|ψ(x)|dx (4.5)
para determinar um limitante superior do erro En. Aqui vamos trabalhar dois casos:
(a) Caso 1: Se ψ(x) ≥ 0 para todo x ∈ [x0, xn], ou se ψ(x) ≤ 0 para todo x ∈ [x0, xn].
(b) Caso 2: Se ψ(x) satisfaz
∫ xn
x0
ψ(x)dx = 0.
4.3 Regra dos Trapézios.
A regra dos trapézios é a mais simple das regras de integração, ou seja, quando con-
sideramos n = 1. Isto é a interpolação será dada sobre dois pontos x0 = a e x1 = b
e o polinômio interpolante P1(x) é uma reta. Definamos h = b − a e seja y0 = f (x0) e
y1 = f (x1). O polinômio P1(x) pode ser encontrado usando as diferenças divididas:
xi yi ordem 0 ordem 1
x0 y0 y0
y1−y0
h
x1 y1 y1
Figura 13: Aproximação usando um polinômio de grau 1.
Pela fórmula de interpolação de Newton, o polinômio interpolante é
P1(x) = y0 +
y1 − y0
h
(x − x0)
Então a aproximação de I( f ) será
I( f ) =
∫ b
a
f (x)dx ≈
∫ x1
x0
P1(x)dx =
∫ x1
x0
[
y0 +
y1 − y0
h
(x − x0)
]
dx.
52
Desenvolvendo as contas da última integral e fazendo x = x0 + u, obtemos:
I( f ) ≈
∫ h
0
[
y0 +
( y1 − y0
h
)
u
]
du =
h
2
[y0 + y1] =
h
2
[ f (x0) + f (x1)].
Observação 4.2. A fórmula
I( f ) ≈
h
2
[ f (x0) + f (x1)].
é denominada Regra do Trapézio porque a área abaixo do polinômio P1(x) corres-
ponde à área de um trapézio.
4.3.1 Limitante superior do erro na Regra do Trapézio
Note que ψ(x) = (x − x0)(x − x1), então ψ(x) ≤ 0 para todo x ∈ [x0, x1]. Isto implica que
|ψ(x)| = −ψ(x) = (x − x0)(x1 − x) para todo x ∈ [x0, x1]. Pela fórmula (4.4) para o erro E1
obtemos
|E1| ≤
K
(1 + 1)!
∫ x1
x0
(x − x0)(x1 − x)dx,
em que
K = max{| f (2)(x)| : x0 ≤ x ≤ x1}.
Fazendo a mudança x = x0 + u, obtemos
|E1| ≤
K
(2)!
∫ h
0
u(h − u)du =
h3
12
max{| f (2)(x)| : x0 ≤ x ≤ x1}.
Resumindo:
A aproximação de I( f ) dada pela Regra do Trapézio é
I( f ) ≈
h
2
[ f (x0) + f (x1)]. (4.6)
e o limitante superior para o erro E1 é dado por
|E1| ≤
h3
12
max{| f (2)(x)| : x0 ≤ x ≤ x1}. (4.7)
Exemplo 4.2. Usando a Regra dos Trapézios, encontre uma aproximação da integral∫ 2
1
1
x
dx
53
e determine um limitante superior para o erro cometido ao aproximar a integral.
Solução:
Definimos x0 = a = 1 e x1 = b = 2 e então h = x1 − x0 = 1. Assim, f (x0) = f (1) = 1 e
f (x1) = f (2) = 0, 5. Por tanto, a aproximação da integral dada pela fórmula (4.6) é:
I( f ) ≈
1
2
[1 + 0, 5] = 0, 75.
Para determinar o limitante superior do erro precisamos da derivada f (2)(x) = 2x3 .
Assim
max{| f (2)(x)| : x0 ≤ x ≤ x1} = max{
2
x3
: 1 ≤ x ≤ 2} =
2
1
= 2.
Por tanto, pela inequação (4.7),
|E1| ≤
13
12
(2) =
2
12
= 0, 16666.
Logo a integral pertence ao intervalo I( f ) ∈ [0.75 − 0.16666, 0.75 + 0.16666] =
[0.58334, 0.91666].
4.3.2 Regra do Trapézio Generalizada
A Regra dos Trapézios Generalizada consiste na subdivisão do intervalo de integração
[a, b] em n sub-intervalos Ii = [xi−1, xi] para i = 1, ...,n, cada um de comprimento h = b−an ,
definidos pelos pontos
x0 = a, x1 = a + h, x2 = a + 2h, ..., xi = a + ih, ..., xn = b,
para logo aplicar a regra dos trapézios em cada sub-intervalo Ii para i = 1, ...,n. Desse
modo, a integral I( f ) será aproximada pela soma das integrais em cada sub-intervalo
calculadas utilizando a regra do trapézio. A figura abaixo apresenta o caso de 5 sub-
intervalos.
Figura 14: Regra do Trapézio Generalizada.
54
A aproximação no intervalo Ii é dada por
Ai =
h
2
[ f (xi−1) + f (xi)].
Assim, a aproximação será:
I( f ) ≈ A1 + A2 + .... + An
Substituindo os valores:
I( f ) ≈
h
2
[ f (x0) + f (x1)] +
h
2
[ f (x1) + f (x2)] + .... +
h
2
[ f (xn−1) + f (xn)]
Simplificando
I( f ) ≈
h
2
[ f (x0) + 2
n−1∑
i=1
f (xi) + f (xn)] (4.8)
O erro na Regra do Trapézios Generalizada é limitado pela soma dos erros em cada
subintervalo Ii. Denotemos por E(i) o erro no intervalo Ii, e então, pela inequação 4.7,
|E(i)| ≤
h3
12
max{| f (2)(x)| : xi−1 ≤ x ≤ xi},∀i = 1, ...,n.
Então o erro total E da regra generalizada é limitado pela expressão:
|E| ≤ E(1) + E(2) + ... + E(n) ≤ n
h3
12
max{| f (2)(x)| : x0 ≤ x ≤ xn}.
Sabendo que nh = xn − x0, obtemos
|E| ≤
h2
12
(xn − x0) max{| f (2)(x)| : x0 ≤ x ≤ xn}. (4.9)
Exemplo 4.3. Encontre aproximadamente a integral∫ 2
0
√
1 + x2dx
e determine um limitante superior para o erro cometido ao aproximá-la pela Regra do
Trapézio Generalizada com n = 4.
55
Solução:
Dos dados f (x) =
√
1 + x2, h = 2−04 = 0, 5 e x0 = 0, x1 = 0, 5, x2 = 1, x3 = 1, 5, x4 = 2.
Pela fórmula
I( f ) =
∫ 2
0
√
1 + x2dx ≈
h
2
[ f (x0) + 2
3∑
i=1
f (xi) + f (x4)]
deduzimos
I( f ) ≈
0, 5
2
[ f (x0) + 2( f (x1) + f (x2) + f (x3)) + f (x4)].
Substituindo os valores dos f (xi) obtemos:
I( f ) ≈
0, 5
2
[1 + 2(1, 1180 + 1, 4142 + 1, 8027) + 2, 2360] = 2, 9763
Observe que f (2)(x) = 1(1+x2)3/2 e
max{| f (2)(x)| : 0 ≤ x ≤ 2} = max{
1
(1 + x2)3/2
: 0 ≤ x ≤ 2} = 1
Dessa forma
|E| ≤
(0.5)2
12
(2 − 0)(1) = 0, 0416.
4.4 Regra 13 de Simpson
A Regra 1/3 de Simpson ocorre quando o polinômio interpolante Pn(x) tem grau n = 2.
Ou seja, f (x) é aproximada por uma parábola. Definamos
h =
b − a
2
, e , x0 = a, x1 = a + h, x2 = b.
e
y0 = f (x0), y1 = f (x1), y2 = f (x2).
Encontraremos P2(x) usando as diferenças divididas dadas na tabela abaixo:
xi yi ordem 0 ordem 1 ordem 2
x0 y0 y0
y1−y0
h
x1 y1 y1
y2−2y1+y0
2h2y2−y1
h
x2 y2 y2
56
Figura 15: Aproximação usando um polinômio de grau 2
Pela fórmula de Interpolação de Newton, o polinômio P2(x) é dado por:
P2(x) = y0 +
y0 − y1
h
(x − x0) +
y2 − 2y1 + y0
2h2
(x − x0)(x − x1).
Integrando entre x0 e x1 e fazendo a mudança de variável x = x0 + u obtemos∫ x2
x0
P2(x)dx =
∫ 2h
0
[y0 +
y0 − y1
h
u +
y2 − 2y1 + y0
2h2
u(u − h)]du.
Desenvolvendo as contas: ∫ x1
x0
P2(x)dx =
h
3
[y0 + 4y1 + y2].
Logo,
I( f ) ≈
h
3
[ f (x0) + 4 f (x1) + f (x2)]
.
Observação 4.3. A Fórmula∫ b
a
f (x)dx ≈
h
3
[ f (x0) + 4 f (x1) + f (x2)] (4.10)
é chamada regra 1/3 de Simpson.
4.4.1 Limitante superior para o Erro E2 na regra 1/3 de Simpson
Neste caso ψ(x) = (x − x0)(x − x1)(x − x2) e desde que x1 = a + h = a+b2 , segue que∫ x2
x0
ψ(x)dx = 0 (4.11)
57
Procedemos da seguinte maneira para estimar o erro cometido. Acrescentamos um
ponto artificial x3 := x1 ao conjunto de pontos x0, x1, x2 e determinaremos o polinômio
interpolante P∗(x) que é definido pelos pontos
x0, x1, x2, x3.
A tabela de diferenças divididas para estes 4 pontos é
xi yi ordem 0 ordem 1 ordem 2 ordem 3
x0 y0 y0
y1−y0
h
x1 y1 y1
y2−2y1+y0
2h2
y2−y1
h
4y1−3y2−y0
2h3
x2 y2 y2
y1−y2
h2y1−y2
h
x1 y1 y1
Pela fórmula interpolátoria de Newton,
P∗(x) = y0+
y0 − y1
h
(x−x0)+
y2 − 2y1 + y0
2h2
(x−x0)(x−x1)+
4y1 − 3y2 − y0
2h3
(x−x0)(x−x1)(x−x2).
Em outras palavras
P∗(x) = P2(x) +
4y1 − 3y2 − y0
2h3
ψ(x).
Integrando entre x0 e x2 e como
∫ x2
x0
ψ(x)dx = 0 obtemos∫ x2
x0
P∗(x)dx =
∫ x2
x0
P2(x)dx +
4y1 − 3y2 − y0
2h3
∫ x2
x0
ψ(x)dx =
∫ x2
x0
P2(x)dx.
Inferimos que aproximar I( f ) pela integral de P2(x) é igual a aproximá-la pela integral
de P∗(x). Por tanto, o erro cometido usando P2(x) é igual ao erro cometido usando P∗(x).
Assim, é suficiente determinar o erro para P∗(x) que é um polinômio de grau 3 = 2 + 1.
Isto é, devemos realizar a mesma análise da seção