APOSTILA
Cálculo Numérico
Universidade Tecnológica Federal do Paraná
UTFPR
Lauro César Galvão, Dr. e Luiz Fernando Nunes, Dr.
ii
Índices
1
NOÇÕES BÁSICAS SOBRE ERROS ...................................................................................1-1
1.1
1.2
1.2.1
1.2.2
1.3
1.3.1
1.3.2
1.4
1.5
1.5.1
1.5.2
1.5.3
1.6
1.6.1
1.6.2
1.6.3
2
ERROS...............................................................................................................................................1-1
ERROS A BSOLUTOS E RELATIVOS ................................................................................................1-1
Erro Absoluto..................................................................................................................................1-1
Erro Relativo ou Taxa de Erro ....................................................................................................1-2
ERROS DE A RREDONDAMENTO E TRUNCAMENTO.....................................................................1-2
Erro de Arredondamento..............................................................................................................1-2
Erro de Truncamento ....................................................................................................................1-2
A RITMÉTICA DE PONTO FLUTUANTE ...........................................................................................1-3
CONVERSÃO DE BASES ..................................................................................................................1-3
Conversão da Base β para a Decimal (β⇒10) ........................................................................1-3
Conversão da Base Decimal para a β (10⇒β) ........................................................................1-4
Exercícios: Conversão de Bases..................................................................................................1-5
OPERAÇÕES DE PONTOS FLUTUANTES ........................................................................................1-7
Representações...............................................................................................................................1-7
Exercícios........................................................................................................................................1-7
Exercícios complementares..........................................................................................................1-8
ZEROS REAIS DE FUNÇÕES REAIS .............................................................................. 2-11
2.1
2.2
2.3
2.3.1
2.3.2
2.3.3
2.3.4
3
INTRODUÇÃO .................................................................................................................................2-11
FASE I: ISOLAMENTO DAS RAÍZES...............................................................................................2-11
FASE II: REFINAMENTO - CRITÉRIOS DE PARADA ....................................................................2-15
Método da Bissecção (ou Método da Dicotomia) ................................................................. 2-15
Método do Ponto Fixo (ou Método da Iteração Linear ou Método das Aproximações
sucessivas).................................................................................................................................... 2-19
Método de Newton, Newton-Raphson (ou Método das Tangentes).................................... 2-27
Comparação entre os métodos.................................................................................................. 2-30
RESOLUÇÃO DE SISTEMAS DE EQUAÇÕES LINEARES .................................... 3-32
3.1
3.1.1
3.1.2
3.1.3
3.1.4
3.1.5
3.1.6
3.2
3.2.1
3.2.2
3.2.3
3.3
3.3.1
3.3.2
3.3.3
3.3.4
3.3.5
4
INTRODUÇÃO .................................................................................................................................3-32
Forma Algébrica de S n ............................................................................................................... 3-32
Forma Matricial de Sn............................................................................................................... 3-32
Matriz Aumentada ou Matriz Completa do Sistema ............................................................. 3-32
Solução do Sistema ..................................................................................................................... 3-32
Classificação de um Sistema Linear........................................................................................ 3-33
Classificação quanto ao Determinante de A.......................................................................... 3-33
M ÉTODOS DIRETOS.......................................................................................................................3-33
Método de Eliminação de Gauss.............................................................................................. 3-33
Estratégia de Pivoteamento Completo .................................................................................... 3-36
Refinamento de Soluções........................................................................................................... 3-37
M ÉTODOS ITERATIVOS .................................................................................................................3-39
Testes de parada.......................................................................................................................... 3-39
Método de Gauss-Jacobi........................................................................................................... 3-39
Método de Gauss-Seidel............................................................................................................ 3-42
Comparação entre os métodos.................................................................................................. 3-43
Critério de Sassenfeld................................................................................................................ 3-44
INTERPOLAÇÃO.................................................................................................................... 4-47
4.1
4.1.1
4.1.2
4.1.3
4.2
4.2.1
4.3
4.3.1
4.3.2
4.4
INTERPOLAÇÃO POLINOMIAL ......................................................................................................4-47
Existência e Unicidade do Polinômio Interpolador Pn(x).................................................... 4-47
Forma de Lagrange.................................................................................................................... 4-48
Forma de Newton........................................................................................................................ 4-50
ESTUDO DE ERRO NA INTERPOLAÇÃO ........................................................................................4-52
Estimativa para o Erro............................................................................................................... 4-52
INTERPOLAÇÃO INVERSA : CASOS EXISTENTES..........................................................................4-54
Encontrar
x tal que Pn ( x ) .................................................................................................. 4-54
Interpolação inversa................................................................................................................... 4-54
FUNÇÕES SPLINE EM INTERPOLAÇÃO.........................................................................................4-56
iii
4.4.1
4.4.2
4.4.3
5
Função Spline .............................................................................................................................. 4-56
Spline linear interpolante.......................................................................................................... 4-57
Spline cúbica interpolante......................................................................................................... 4-58
AJUSTE DE CURVAS PELO MÉTODO DOS MÍNIMOS QUADRADOS ........... 5-64
5.1
5.2
5.3
5.4
6
INTRODUÇÃO .................................................................................................................................5-64
CASO DISCRETO............................................................................................................................5-65
CASO CONTÍNUO ...........................................................................................................................5-70
FAMÍLIA DE FUNÇÕES NÃO LINEARES NOS PARÂMETROS .....................................................5-72
INTEGRAÇÃO NUMÉRICA ............................................................................................... 6-74
6.1
6.1.1
6.1.2
6.1.3
6.1.4
7
FÓRMULAS DE NEWTON -COTES .................................................................................................6-74
Regra dos Trapézios................................................................................................................... 6-74
Regra dos Trapézios repetida ................................................................................................... 6-76
Regra 1/3 de Simpson................................................................................................................. 6-77
Regra 1/3 de Simpson repetida................................................................................................. 6-80
SOLUÇÃO NUMÉRICA DE EQUAÇÕES DIFERENCIAIS ORDINÁRIAS ....... 7-83
7.1
7.2
7.2.1
7.2.2
7.2.3
7.2.4
7.2.5
INTRODUÇÃO .................................................................................................................................7-83
PROBLEMA DE VALOR INICIAL (PVI) .........................................................................................7-84
Solução numérica de um PVI de primeira ordem ................................................................. 7-84
Método de Euler.......................................................................................................................... 7-84
Métodos de Runge-Kutta............................................................................................................ 7-87
Método de Euler Aprimorado (Método de Runge-Kutta de Segunda Ordem) ................. 7-88
Fórmulas de Runge-Kutta de Quarta Ordem......................................................................... 7-89
iv
Índices de Figuras
[FIG. 1]:
[FIG. 2]:
M ODELAGEM E RESOLUÇÃO DE PROBLEMAS..............................................................................1-1
O GRÁFICO DE UMA FUNÇÃO y = f ( x ) E SEUS ZEROS.........................................................2-11
[FIG. 3]:
EXEMPLO DE UMA FUNÇÃO ESTRITAMENTE CRESCENTE NUM INTERVALO DE a ATÉ b ..2-12
[FIG. 4]:
O GRÁFICO DE
[FIG. 5]:
O S GRÁFICOS DE
[FIG. 6]:
O GRÁFICO DE
f ( x ) = x 3 − 9 x + 3 .......................................................................................2-12
g ( x ) = x 3 E h( x ) = 9 x − 3 . ....................................................................2-13
[FIG. 9]:
f ' ( x ) = 3 x 2 − 9 . ...........................................................................................2-13
GRÁFICO DA FUNÇÃO f ( x ) = x ln x − 3,2 ..........................................................................2-14
GRÁFICO DA FUNÇÃO f ' ( x ) = 1 + ln x .................................................................................2-14
O S GRÁFICOS DE g ( x ) = 5 log x E h( x ) = 2 − 0,4 x . ......................................................2-15
[FIG. 10]:
[FIG. 11]:
[FIG. 12]:
[FIG. 13]:
O S GRÁFICOS DE g ( x ) = x E h( x ) = 5e
. ................................................................2-15
O MÉTODO DA BISSECÇÃO OU DICOTOMIA ................................................................................2-16
O TANQUE DE COMPRIMENTO L . ..............................................................................................2-17
UM EXEMPLO DE UMA FUNÇÃO DE PONTO FIXO.......................................................................2-19
[FIG. 14]:
O S GRÁFICOS DAS FUNÇÕES
y =x
E
φ 2 ( x ) = 6 − x ......................................................2-20
[FIG. 15]:
O S GRÁFICOS DAS FUNÇÕES
y =x
E
φ1 ( x ) = 6 − x 2 . ........................................................2-21
[FIG. 16]:
A SEQÜÊNCIA
[FIG. 7]:
[FIG. 8]:
−x
[FIG. 19]:
{x k } CONVERGE PARA O ZERO α (CONVERGÊNCIA DO TIPO ESCADA).......2-22
A SEQÜÊNCIA {x k } CONVERGE PARA O ZERO α (CONVERGÊNCIA DO TIPO CARACOL)....2-22
A SEQÜÊNCIA {x k } NÃO CONVERGE PARA O ZERO α.............................................................2-22
A SEQÜÊNCIA {x k } NÃO CONVERGE PARA O ZERO α.............................................................2-23
[FIG. 20]:
CASOS EM QUE B É O EXTREMO MAIS PRÓXIMO DE α..............................................................2-24
[FIG. 21]:
O S GRÁFICOS DE
[FIG. 22]:
[FIG. 23]:
INTERPRETAÇÃO GEOMÉTRICA DO MÉTODO DE NEWTON.....................................................2-28
O S GRÁFICOS DAS FUNÇÕES g (x ) = X E h (x) = cos x . .......................................................2-30
[FIG. 24]:
INTERPOLAÇÃO DE
[FIG. 25]:
[FIG. 26]:
INTERPOLAÇÃO POR LAGRANGE .................................................................................................4-50
GRÁFICO DO POLINÔMIO P10 ( x ) INTERPOLANDO f ( x ) .....................................................4-56
[FIG. 27]:
[FIG. 28]:
[FIG. 29]:
[FIG. 30]:
[FIG. 31]:
[FIG. 32]:
[FIG. 33]:
[FIG. 34]:
[FIG. 35]:
[FIG. 36]:
[FIG. 37]:
SPLINE LINEAR INTERPOLANDO 4 PONTOS. ...............................................................................4-57
DOMÍNIO DISCRETO......................................................................................................................5-64
DOMÍNIO CONTÍNUO.....................................................................................................................5-64
O MÉTODO DO MÍNIMOS QUADRADOS........................................................................................5-65
DIAGRAMA DE DISPERSÃO...........................................................................................................5-68
REGRA DOS TRAPÉZIO ..................................................................................................................6-74
REGRA DOS TRAPÉZIOS REPETIDA ..............................................................................................6-76
REGRA 1/3 DE SIMPSON...............................................................................................................6-78
REGRA 1/3 DE SIMPSON REPETIDA .............................................................................................6-80
GRÁFICO DA SOLUÇÃO NUMÉRICA DE UM PVI.........................................................................7-84
GRÁFICO DO MÉTODO DE EULER................................................................................................7-85
[FIG. 17]:
[FIG. 18]:
h (x ) = e x
E
g ( x ) = x 2 − 4 . ...................................................................2-26
f ( x ) PELO POLINÔMIO P ( x ). .............................................................4-47
Cálculo Numérico
Noções básicas sobre erros
1-1
1 Noções básicas sobre Erros
Fenômenos da natureza podem ser descritos através do uso de modelos matemáticos.
MODELAGEM
PROBLEMA
MODELO
MATEMÁTICO
RESOLUÇÃO
SOLUÇÃO
[Fig. 1]: Modelagem e resolução de problemas.
• MODELAGEM: é a fase de obtenção de um modelo matemático que descreve o
comportamento do problema que se quer estudar.
• RESOLUÇÃO: é a fase de obtenção da solução do modelo matemático através da
aplicação de métodos numéricos.
1.1 Erros
Para se obter a solução do problema através do modelo matemático, erros são
cometidos nas fases: MODELAGEM e RESOLUÇÃO.
Exercício 1
Calcular a área da superfície terrestre usando a formulação A =4π r 2 .
Aproximações (ERROS):
MODELAGEM:
Resolução:
RESOLUÇÃO:
OBS. 1:
Características do planeta Terra.
• Características Físicas:
Diâmetro Eq uatorial: 12756Km;
Diâmetro Polar: 12713Km;
Massa: 5,98× 1024 Kg;
Perímetro de Rotação Sideral: 23h 56min 04seg;
Inclinação do Equador Sobre a Órbita: 23o 27’.
• Características Orbitais:
Raio da Órbita, isto é, 1U.A. (unidade astronômica): 149897570Km;
Distância Máxima do Sol: 152100000Km;
Distância Mínima do Sol: 147100000Km;
Período de Revolução Sideral: 365dias 6h 9min 9,5seg;
Velocidade Orbital Média: 29,79Km/seg.
1.2 Erros Absolutos e Relativos
1.2.1 Erro Absoluto
É o módulo da diferença entre um valor exato x de um número e seu valor
aproximado x .
EAx = x − x ,
onde x é o valor exato e x é o valor aproximado.
(Eq.1)
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
1-2
Cálculo Numérico
Noções básicas sobre erros
Geralmente não se conhece o valor exato x . Assim, o que se faz é obter um limitante
superior ( k1 majorante) ou uma estimativa para o módulo do erro absoluto.
EAx ≤ k1 .
(Eq.2)
1.2.2 Erro Relativo ou Taxa de Erro
Erro relativo de x é o módulo do quociente entre o erro absoluto EAx e o valor exato
x ou o valor aproximado x , se x ou x ≠ 0.
(Eq.3)
ERx =
EAx
x−x
EAx x − x
=
ou ERx =
=
.
x
x
x
x
Calcular os erros absoluto e relativo, nos itens a) e b).
a) x =1,5 e x =1,49;
b) y =5,4 e y =5,39.
Exercício 2
Resolução:
1.3 Erros de Arredondamento e Truncamento
1.3.1 Erro de Arredondamento
Arredondar um número na casa d i é desconsiderar as casas d i + j ( j =1,…,∞) de tal
forma que:
d i seja a última casa se d i +1 <5;
d i +1 seja a última casa se d i +1 ≥5.
Exercício 3
Arredondar π na quarta casa decimal, sendo que π=3,1415926535…
Resolução:
1.3.2 Erro de Truncamento
Truncar um número na casa d i é desconsiderar as casas d i + j ( j =1,…,∞).
Exercício 4
Aproximar π truncando na quarta casa decimal, sendo que π=3,1415926535…
Resolução:
Exercício 5
∞
xi
Sabendo-se que e pode ser escrito como e = ∑ , faça a aproximação de
i =0 i!
x
x
e 2 através de um truncamento após quatro termos da somatória.
Resolução:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Noções básicas sobre erros
1-3
1.4 Aritmética de Ponto Flutuante
Um número é representado, internamente, na máquina de calcular ou no computador
através de uma seqüência de impulsos elétricos que indicam dois estados: 0 ou 1, ou seja, os
números são representados na base 2 ou binária.
De maneira geral, um número x é representado na base β por:
(Eq.4)
d d
d
d 
x =±  1 + 22 + 33 +…+ tt  ∗ βexp .
β 
β
β β
Onde:
• d i ⇒ são números inteiros contidos no intervalo 0≤ d i <β; i =1, 2, …, t ;
• exp ⇒ representa o expoente de β e assume valores entre I ≤ exp ≤ S ;
• I , S ⇒ limite inferior e limite superior, respectivamente, para a variação do expoente;
d d
d
d 
•  1 + 22 + 33 +…+ tt  ⇒ é chamada de mantissa e é a parte do número que representa
β 
β
β β
seus dígitos significativos;
• t ⇒ número de dígitos do sistema de representação.
Considerando no sistema de base 10, β=10, represente os seguintes números,
em aritmética de ponto flutuante:
a) 0,34510 ;
b) 31,41510 .
Exercício 6
Resolução:
Os números assim representados estão NORMALIZADOS, isto é, a mantissa é
um número entre 0 e 1.
OBS. 2:
Considerando no sistema binário, β=2, represente o número 1012 em aritmética
de ponto flutuante.
Exercício 7
Resolução:
1.5 Conversão de Bases
1.5.1 Conversão da Base β para a Decimal (β⇒10)
Um número na base β pode ser escrito, na base decimal, como:
m
(Eq.5)
∑ aiβi = am βm + am −1β m−1 +…+ a2β2 + a1 β+ a0 + a−1β−1 + a− 2β−2 +…+ an +1βn +1 + anβn .
i=n
Onde:
• ai ⇒ 0≤ ai <β;
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
• n , m ⇒ números inteiros, com n ≤0 e m ≥0.
Noções básicas sobre erros
1-4
Para a conversão, faz-se a operação entre a mantissa do número normalizado e a base
β
exp
.
Nos exercícios a seguir, faça a conversão da base indicada para a decimal,
determinando o valor da variável x .
Exercício 8
10112 = x10 .
Resolução:
Exercício 9
11,012 = x10 .
Resolução:
Exercício 10
403,125 = x10 .
Resolução:
1.5.2 Conversão da Base Decimal para a β (10⇒β)
Aplica-se um processo para a parte inteira e um outro para a parte fracionária.
• a) PARTE INTEIRA ( N ):
• a.1) N <β
⇒ N10 = Nβ .
• a.2) N ≥β
N
r1
β
q1
r2
β
q2
O
O
O
qn −1
rn
β
qn
⇒
Exercício 11
Até que q n <β
N10 =( q n rn rn −1 … r3 r2 r1 ) β
Converta 5910 para a base 2.
Resolução:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Exercício 12
Noções básicas sobre erros
1-5
Converta 5910 para a base 3.
Resolução:
• b) PARTE FRACIONÁRIA ( F ):
Multiplica-se F por β e toma-se a parte inteira do produto como o primeiro dígito do
número na base β. Repete-se o processo com a parte fracionária do produto tomando sua parte
inteira. Continua-se até que a parte fracionária seja igual a zero.
Nos exercícios a seguir, determinar o valor de x :
Exercício 13
0,187510 = x2 .
Resolução:
Exercício 14
0,610 = x2 .
Resolução:
Exercício 15
13,2510 = x2 .
Resolução:
1.5.3 Exercícios: Conversão de Bases
Transforme para a base que se pede (determine o valor de x ).
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Exercício 16
100101,10012 = x10 .
Noções básicas sobre erros
Resolução:
Exercício 17
19,3867187510 = x4 .
Resolução:
Transforme a medida 35 h 48 min 18 seg para minutos.
DICA: 35:48,1860 = x10 min .
Exercício 18
Resolução:
Transforme 35,805 horas para horas, minutos e segundos.
DICA: 35,80510 = x60 .
Exercício 19
Resolução:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
1-6
Cálculo Numérico
Noções básicas sobre erros
1-7
1.6 Operações de Pontos Flutuantes
1.6.1 Representações
• Precisão dupla: “dobra” a mantissa (2∗ t );
• O zero em ponto flutuante é em geral representado com o menor expoente (exp = I )
possível na máquina;
• Ao converter um número para determinada aritmética de ponto flutuante, emprega-se
sempre o arredondamento;
• Não é possível representar todos os números reais em determinada aritmética de ponto
flutuante (reta furada).
Um exemplo da reta furada é: Considere a aritmética de pontos flutuantes com
parâme tros β=10 e t =3. Tome os números consecutivos 3,57 e 3,58. Existem infinitos
números reais entre 3,57 e 3,58 que não podem ser representados nesta aritmética de pontos
flutuantes. Por exemplo: 3,571 ou 3,57437.
OBS. 3:
1.6.2 Exercícios
Exercício 20
Preencher a tabela a seguir, com base nos parâmetros: t =3, β=10, I =−5, S =5
e −5≤ exp ≤5.
Número
−6,48
0,0002175
3498,3
−0,00000001452
2379441,5
Truncamento
Arredondamento
Deve-se converter os valores para a aritmética de ponto flutuante com 3
algarismos significativos.
OBS. 4:
Nos exercícios seguintes, calcular o valor das expressões utilizando aritmética de
ponto flutuante com 3 algarismos significativos.
Exercício 21
(4,26 + 9,24) + 5,04
Resolução:
Exercício 22
4,26 + (9,24 + 5,04)
Resolução:
Exercício 23
(4210 − 4,99) − 0,02
Resolução:
Exercício 24
4210 − (4,99 + 0,02)
Resolução:
Exercício 25
2
∗(4,0237 − 6,106)
7
Resolução:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
2 ∗ ( 4,0237 − 6,106)
Exercício 26
7
Noções básicas sobre erros
1-8
Resolução:
OBS. 5:
Em aritmética de ponto flutuante não valem as propriedades associativas nem
distributivas.
Exercício 27
a) 42450 +
Sendo β=10, t =4 e exp ∈[−5,5], calcule:
10
∑3;
i =1
10
b)
∑ 3 + 42450.
i =1
Resolução:
1.6.3 Exercícios complementares
Nos exercícios seguintes, converter os números para a base decimal, determinando o
valor da variável x :
Exercício 28
11000112 = x10 .
Resolução:
Exercício 29
11111112 = x10 .
Resolução:
Exercício 30
10101012 = x10 .
Resolução:
.
Exercício 31
101,00112 = x10 .
Resolução:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Exercício 32
Noções básicas sobre erros
1-9
0,01111112 = x10 .
Resolução:
Exercício 33
1,0100112 = x10 .
Resolução:
Nos exercícios seguintes, converter os números para a base binária, determinando o
valor da variável x :
Exercício 34
3710 = x2 .
Resolução:
Exercício 35
234510 = x2 .
Resolução:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Noções básicas sobre erros
Exercício 36
Determine x com 36 dígitos: 0,121710 = x2 .
Resolução:
Exercício 37
Determine x com 8 dígitos: 2,4710 = x2 .
Resolução:
• Logo: 2,4710 = 210 + 0, 4710 = 102 + 0,011110002 = 10, 011110002 .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
1-10
Cálculo Numérico
Zeros reais de funções reais
2-11
2 Zeros reais de funções reais
2.1 Introdução
Dada uma função real f definida e contínua em um intervalo aberto I , chama-se de
zero desta função em I , a todo x ∈ I , tal que f( x ) = 0.
Neste capítulo são apresentados alguns processos iterativos para calcular de forma
aproximada os zeros reais de uma função real f dada. Por um processo iterativo entende-se
um processo que calcula uma seqüência de aproximações x1 , x 2 , x3 ,… da solução desejada. O
cálculo de uma nova aproximação é feito utilizando aproximações anteriores. Dizemos que a
seqüência x1 , x 2 , x3 ,… converge para x , se dado ε>0, ∃ N ∈
(
números naturais), tal
que qualquer que seja n > N , xn − x <ε. Neste caso tem-se que lim xn = x , o que também
n →∞
poderá ser indicado por xn → x . Nos processos iterativos que serão apresentados, a
determinação dos zeros de uma função real de variável real será feita em duas etapas:
• Fase I: Isolar cada zero que se deseja determinar da função f em um intervalo [ a , b ],
sendo que cada intervalo deverá conter um e somente um zero da função f .
• Fase II: Cálculo dos zeros aproximados utilizando um método iterativo, com precisão
prefixada ou não.
2.2 Fase I: Isolamento das raízes
Seja f ( x ) uma função contínua num intervalo [a , b ]. Se f ( a )⋅ f ( b )<0,
então existe pelo menos um zero de f ( x ) entre a e b .
y
y = f (x)
Teorema 1
x
a
b
[Fig. 2]: O gráfico de uma função
y = f ( x ) e seus zeros.
Sob as hipóteses do teorema 1, o zero x =α será definido e único em [ a , b ] se
a derivada f ' ( x ) existir e preservar o sinal dentro do intervalo ] a , b [, isto é se f ' ( x )>0,
∀x ∈] a , b [ ou f ' ( x )<0, ∀x ∈] a , b [. Isto significa dizer que a função f ( x ) é estritamente
crescente ou estritamente decrescente, respectivamente, no intervalo ] a , b [.
OBS. 6:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
2-12
Zeros reais de funções reais
y
y = f (x)
x
a
b
[Fig. 3]: Exemplo de uma função estritamente crescente num intervalo de
a
até
b.
Na pesquisa dos zeros reais de funções reais é muito útil o uso do Teorema 1Erro! A
origem da referência não foi encontrada. (que fornece condições de existência de zeros em
um intervalo), bem como da OBS 6. (que garante a unicidade, isto é, garante que no intervalo
considerado existe um e somente um zero da função f ).
Outro recurso bastante empregado é: a partir da equação f ( x )=0, obter a equação
equivalente g ( x )= h ( x ) e esboçar os gráficos destas funções obtendo os pontos onde as
mesmas se intersectam, pois f (α)=0 ⇔ g (α)= h (α).
Exercício 38
Resolução:
x
f (x)
Isolar os zeros da função f ( x )= x 3 −9 x +3.
Pode-se construir uma tabela de valores para f ( x ) e analisar os sinais:
0
1
2
−4
−3
−2
−1
3
y
y = f (x)
α1
-4 -3 -2 -1
[Fig. 4]: O gráfico de
α3
α2
1
2
x
3
4
f ( x) = x3 − 9x + 3 .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Zeros reais de funções reais
y
α1
h(x )
x
-4 -3 -2 -1
[Fig. 5]: Os gráficos de
g(x)
α2 1
2 α3 3
4
g ( x ) = x 3 e h( x ) = 9 x − 3 .
y
y = f ’( x)
-
3
-4 -3 -2 -1
[Fig. 6]: O gráfico de
Exercício 39
Resolução:
x
3
1
2
3
4
f '( x ) = 3 x 2 − 9 .
Isolar os zeros da função f ( x ) = x ln x − 3,2 .
Pode-se construir uma tabela de valores para f ( x ) e analisar os sinais:
x
1
2
3
4
f ( x)
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
2-13
Cálculo Numérico
Zeros reais de funções reais
y
y = f(x)
0,3
0,2
0,1
0
-0,1
2,6
2,8
3,0
3,2
x
3,4
-0,2
-0,3
-0,4
-0,5
-0,6
-0,7
-0,8
-0,9
-1,0
[Fig. 7]: Gráfico da função
f ( x ) = x ln x − 3,2 .
y
f’ (x )
1
1
[Fig. 8]: Gráfico da função
Exercício 40
Resolução:
x
f ' ( x ) = 1 + ln x .
Isolar os zeros da função f ( x ) = 5 log x − 2 + 0,4 x .
Pode-se construir uma tabela de valores para f ( x ) e analisar os sinais:
x
1
2
3
f ( x)
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
2-14
Cálculo Numérico
Zeros reais de funções reais
y
2-15
g (x)
2
1
h(x)
α2
1
[Fig. 9]: Os gráficos de
Exercício 41
Resolução:
3
x
g ( x ) = 5 log x e h( x ) = 2 − 0,4 x .
Isolar os zeros da função f ( x ) = x − 5e − x .
Pode-se construir uma tabela de valores para f ( x ) e analisar os sinais:
x
0
1
2
3
f ( x)
y
2
g (x)
1
h(x)
1
[Fig. 10]: Os gráficos de
α
g ( x ) = x e h( x ) = 5e− x
2
3
x
.
2.3 Fase II: Refinamento - Critérios de Parada
2.3.1Método da Bissecção (ou Método da Dicotomia)
Este método é normalmente utilizado para diminuir o intervalo que contém o zero da
função, para a aplicação de outro método, pois o esforço computacional cresce
demasiadame nte quando se aumenta a precisão exigida.
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Zeros reais de funções reais 2-16
O processo consiste em dividir o intervalo que contém o zero ao meio e por aplicação
do Teorema 1, aplicado aos subintervalos resultantes, determinar qual deles contém o zero.
 a + b  a + b 
a, 2  ,  2 , b 

 

O processo é repetido para o novo subintervalo até que se obtenha uma precisão
prefixada. Desta forma, em cada iteração o zero da função é aproximado pelo ponto médio de
cada subintervalo que a contém.
y
f (x)
m3
a
m2
m1
b
x
α
[Fig. 11]: O método da bissecção ou dicotomia.
Assim, na figura anterior tem-se:
m1 =
a +b
a + m1
m + m1
, m2 =
, m3 = 2
,…
2
2
2
Desta forma, o maior erro que se pode cometer na:
• 1a iteração ( n =1): é
(b − a )
2
• 2a iteração ( n =2): é
(b − a )
22
• 3a iteração ( n =3): é
(b − a )
23
M
M
• n a iteração:
é
(b − a )
2n
Se o problema exige que o erro cometido seja inferior a um parâmetro ε, determina-se
(b − a )
a quantidade n de iterações encontrando o maior inteiro que satisfaz a inequação:
≤ε
2n
que se resolve da seguinte maneira:
(b − a )
(b − a )
≤ε
⇒
log
≤ log ε ⇒ log( b − a ) − log 2 n ≤ log ε ⇒ log( b − a ) − n log 2 ≤ log ε
n
n
2
2
log( b − a ) − log ε
⇒n ≥
log 2
Exercício 42
Determinar um valor aproximado para
Universidade Tecnológica Federal do Paraná (UTFPR)
5 , com erro inferior a 10−2 .
LAURO / NUNES
Cálculo Numérico
Resolução:
n
1
2
3
4
5
6
7
Portanto
Determinar
a
Zeros reais de funções reais
2-17
5 é equivalente a obter o zero positivo da função f ( x ) = x 2 −5.
x
b
f (a )
f (x)
f (b )
( b − a )/2
5≅
Um tanque de comprimento L tem uma secção transversal no formato de um
semicírculo com raio r (veja a figura). Quando cheio de água até uma distância h do topo, o


h
volume V da água é: V= L ⋅ 0,5 ⋅ π ⋅ r 2 − r 2arcsen  − h ( r 2 − h 2 )  . Supondo que L =10 ft ,
r


3
r=1 ft e V=12,4 ft , encontre a profundidade da água no tanque com precisão de 0,01 ft .
Exercício 43
r
h
[Fig. 12]: O tanque de comprimento
θ
h
L.
Resolução:
Pode-se construir uma tabela de valores para f ( x ) e analisar os sinais:
h
0
1
−1
f (h )
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Zeros reais de funções reais 2-18
Para se confirmar a unicidade deste zero neste intervalo, pode-se utilizar a OBS. 6: , isto
é, calcula-se a derivada f , ( h) de f (h ) para verificar que a mesma preserva o sinal no
intervalo ]0,1[.
n
a
1
2
3
4
5
6
7
Assim, h =
h
b
f (a )
f (h )
f (b)
(b−a)/2
2.3.1.1 Algoritmo do Método da Bissecção
Seja f (x ) uma função contínua em um intervalo [a,b], com f (a ) . f (b) <0 e a raiz de
f (x ) isolada em [ a , b ].
• Dados de Entrada : Pontos extremos a e b do intervalo; precisão ou tolerância (ε) e o
número máximo de iterações (ITMAX).
• Saída: Solução aproximada x ou mensagem de "solução não encontrada" com a precisão
desejada no número máximo de iterações.
PASSO 1
Faça i =1
FA= f (a )
PASSO 2
Enquanto i ≤ ITMAX execute os passos de 3 a 6
PASSO 3
Faça x =
( a + b)
e FX = f ( x )
2
PASSO 4
Se FX = 0 ou
(b − a )
< ε, então
2
Saída ( x ) (Procedimento executado com sucesso)
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
FIM
Zeros reais de funções reais
2-19
PASSO 5
Faça i = i +1
PASSO 6
Se FA·FX > 0 então faça a = x e FA = FX
Caso contrário faça b = x
PASSO 7
Saída (Solução não encontrada com a precisão exigida)
FIM
2.3.2Método do Ponto Fixo (ou Método da Iteração Linear
ou Método das Aproximações sucessivas)
Neste método a seqüência de aproximações do zero α de uma função f (x )
( f (α) = 0 ) é obtida através de uma relação de recorrência da forma:
(Eq.6)
x n+1 = φ ( x n ) , n = 0, 1, 2, …
O ponto x 0 será considerado uma aproximação inicial do zero α da função f (x ) e
φ( x) é uma função que tem α como ponto fixo, isto é, α= φ (α ) .
A primeira pergunta a ser respondida é: dada uma função f (x ) com zero α, como
encontrar uma função φ( x) que tenha α como ponto fixo ? Isto pode ser feito através de uma
série de manipulações algébricas sobre a equação f (x ) =0, transformando-a em uma equação
equivalente da forma x = φ( x) . Nestas transformações deve-se tomar os devidos cuidados
para que φ( x) esteja definida em α e para que α pertença à imagem de φ . Como o zero α é
desconhecido, é necessário determinar um intervalo I que contenha α e que esteja contido
tanto no domínio quanto na imagem de φ . É necessário que o zero α de f (x ) seja único no
intervalo I, caso contrário não será possível discernir qual o zero determinado.
y
y= x
φ (x)
Ponto fixo de φ (x)
(Zero de f (x) )
α
x
[Fig. 13]: Um exemplo de uma função de ponto fixo.
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Zeros reais de funções reais
2-20
Exercício 44
Obter algumas funções de ponto fixo para a função f (x ) = x 2 + x − 6 .
Resolução:
Efetuando diferentes manipulações algébricas sobre a equação f (x ) =0 ou
x 2 + x − 6 =0, pode-se obter diferentes funções de ponto fixo, como por exemplo:
No próximo passo algumas destas funções serão utilizadas na tentativa de gerar
seqüências aproximadoras dos zeros α de f (x ) .
Exercício 45
Aproximar o maior zero da função f (x ) = x 2 + x − 6 , utilizando a função
φ 2 ( x ) = 6 − x , e x 0 =1,5.
Resolução:
Neste caso a fórmula de recorrência x n+1 = φ ( x n ) , n = 0, 1, 2, … será:
x n+1 = φ 2 ( x n ) = 6 − x n , e pode-se construir a seguinte tabela:
n
xn
x
0
1
2
3
4
M
n +1
= φ 2 (x n ) = 6 − x n
M
y
M
y= x
6
φ 2 (x)
x0
[Fig. 14]: Os gráficos das funções
x2 x3 x1
α =2
6
x
y = x e φ2 ( x) = 6 − x .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Exercício 46
Zeros reais de funções reais
2-21
Aproximar o maior zero da função f (x ) = x 2 + x − 6 , utilizando a função
φ1 ( x) = 6 − x 2 , e x 0 =1,5.
Resolução:
Neste caso a fórmula de recorrência x n+1 = φ ( x n ) , n =0, 1, 2, … será:
x n+1 = φ1( x n ) = 6 − xn2 , e pode-se construir a seguinte tabela:
n
xn
0
1
2
3
M
x n+1 = φ1 ( x n ) = 6 − x 2
M
M
y
y= x
6
x2
x0
α =2
x1
x
φ1 (x )
[Fig. 15]: Os gráficos das funções
y = x e φ1 ( x ) = 6 − x 2 .
Assim, os dois exercícios anteriores mostram que dependendo da transformação
x = φ( x) escolhida, a relação de recorrência x n+1 = φ ( x n ) pode ou não fornecer uma
seqüência {x n } convergente. Desta forma, como determinar a priori, quais transformações
fornecerão seqüências convergentes? As figuras que seguem ilustram alguns casos onde
ocorrem convergência e alguns casos onde não ocorrem convergência.
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Zeros reais de funções reais
y
y= x
φ (x)
α x3 x2
[Fig. 16]: A seqüência
x1
x0
x
{x k } converge para o zero α (Convergência do tipo escada).
y
y= x
φ(x)
x1
[Fig. 17]: A seqüência
x3 α x 4 x 2
x0
{x k } converge para o zero α (Convergência do tipo caracol).
y
y= x
φ (x)
α x0
[Fig. 18]: A seqüência
x
x1
x2
x3
x
{x k } não converge para o zero α.
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
2-22
Cálculo Numérico
Zeros reais de funções reais
φ (x)
y
x3
[Fig. 19]: A seqüência
x1 α x0
2-23
y= x
x2
x
{x k } não converge para o zero α.
O Teorema que segue estabelece condições suficientes para garantir a convergência do
processo iterativo.
Como as condições que o teorema que segue são apenas suficientes, dada uma
função φ que não satisfaça estas condições, não se pode garantir que a seqüência gerada
x1 , x 2 , x3 , … diverge.
OBS. 7:
2.3.2.1 Convergência do Método das Aproximações Sucessivas
Seja α um zero de uma função f, isolada em um intervalo I=[a,b], e seja φ
uma função tal que φ(α ) = α . Se:
Teorema 2
i)
φ e φ ' são funções contínuas em I;
ii)
k = max φ' (x ) < 1
iii)
x 0 ∈ I e x n+1 = φ ( x n ) ∈ I , para n = 0, 1, 2, …
x∈I
Então a seqüência {x n } converge para o zero α .
Para se resolver um problema com o método das aproximações sucessivas,
utiliza-se o teorema anterior da seguinte forma: inicialmente determina-se um intervalo I
onde o zero α de f (x ) esteja isolado, e uma função φ que tenha α como ponto fixo.
Analisando φ e φ' , pode-se verificar se as condições i) e ii) do Teorema 2 estão satisfeitas.
Estas condições podem não estar satisfeitas pelo fato do intervalo I ter sido
superdimensionado. Neste caso procura-se por um intervalo I ’ satisfazendo as condições do
teorema. Na demonstração do Teorema 2 , que pode ser vista em HUMES, Ana Flora C., et
al. Noções de Cálculo Numérico. São Paulo: McGraw-Hill, p. 16, 1984, tem-se que as
condições i) e ii) garantem que se x n−1 ∈ I então α − x n < α − xn −1 . Entretanto, isto não
OBS. 8:
implica que x n ∈ I . Uma maneira simples para garantir que x n ∈ I = [a, b] ∀n ≥ 0 é tomar
como valor inicial x0 o extremo de I mais próximo do zero α. Na seqüência, será mostrado
que neste caso x1 = φ( x 0 ) ∈ I : Supondo que a seja o extremo de I mais próximo de α, temse: x1 − α < x0 − α = a − α ≤ b − α , logo x1 ∈ I . A demonstração é análoga para o caso em
que b o extremo de I mais próximo de α.
A condição iii) do Teorema 2 pode ser substituída por: iii’) o zero α é o ponto
médio do intervalo I . Na verdade, se para o intervalo I = [a, b ] , estão satisfeitas as condições
OBS. 9:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Zeros reais de funções reais 2-24
i) e ii) do Teorema 2, e se a estiver mais próximo de α do que de b então, denotando a − α
por r, tem-se que para qualquer x0 ∈ [a, α + r ] a hipótese iii) do teorema é verificada. Mais
ainda, para todo I = [a, b ] nas condições do teorema 2, existe I ’⊂ I tal que qualquer que seja
x0 ∈ I ’ tem-se que x n ∈ I ’, n ≥1.
A determinação do extremo de I = [a, b] mais próximo do zero α pode ser feito
da seguinte maneira: Suponhamos satisfeitas as hipóteses i) e ii) do Teorema 2. Nestas
(a + b )
condições, seja x̂ =
(ponto médio do intervalo I ). Sabe-se que φ( xˆ ) está mais
2
próximo de α do que x̂ . Se x̂ < φ( xˆ ) , então α está entre x̂ e b , ou seja, b é o extremo de I
mais próximo de α . Analogamente, se x̂ > φ( xˆ ) , então a é o extremo de I mais próximo de
α. Se x̂ = φ( xˆ ) , então x̂ é o zero procurado.
OBS. 10:
a
x
b
φ (x)
a
x
α
x
b
α
φ (x) x
[Fig. 20]: Casos em que b é o extremo mais próximo de α.
Sejam dados φ( x ), α e k = max φ' ( x ) satisfazendo as hipóteses do teorema
OBS. 11:
x∈I
k
anterior. Se xn =φ( xn−1 ), então α − xn ≤
xn − xn−1 . Desta forma, obtém-se um limitante
1− k
superior para o erro cometido na n -ésima iteração ( xn ).
Verificar as condições i) e ii) do teorema anterior quando do uso da função
φ 2 ( x ) = 6 − x no exercício anterior.
Exercício 47
Resolução:
Verificação da condição i):
Verificação da condição ii):
Logo,
Verificar as condições i) e ii) do teorema anterior quando do uso da função
φ1 ( x) = 6 − x .
Exercício 48
2
Resolução:
Verificação da condição i):
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Verificação da condição ii):
Zeros reais de funções reais
2-25
Logo,
2.3.2.2 Algoritmo do Método das aproximações sucessivas
Para encontrar uma solução para p = φ( p ) dada um aproximação inicial p 0 .
• Dados de Entrada : Aproximação inicial p 0 , precisão ou tolerância (ε) e o número
máximo de iterações (ITMAX).
• Saída: Solução aproximada p ou mensagem de “solução não encontrada”.
PASSO 1
Faça i = 1
PASSO 2
Enquanto i ≤ ITMAX, execute os passos 3 – 6
PASSO 3
Faça p = φ( p ) (calcular pi )
PASSO 4
Se p − p 0 < ε então
Saída ( p ) (procedimento efetuado com sucesso)
FIM
PASSO 5
Faça i = i + 1
PASSO 6
Faça p 0 = p (atualize p 0 )
PASSO 7
Saída (solução não encontrada após ITMAX iterações)
FIM
OBS. 12:
Outros critérios de parada podem ser utilizados:
•
pn − pn −1 < ε
•
pn − p n−1
<ε
pn
•
f ( pn ) < ε
Encontrar o zero de f (x ) = e x − x 2 + 4 com precisão ε = 10 −6 , utilizando o
método do ponto fixo.
Exercício 49
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Zeros reais de funções reais
Resolução:
Pode-se construir uma tabela de valores para f ( x ) e analisar os sinais:
x
−3
−2
−1
f (x )
5
h(x ) = e x
g (x) = x 2 - 4
y
4
3
2
1
-3 α -2 -1
-1
1
2
3
x
-2
-3
-4
[Fig. 21]: Os gráficos de
h (x ) = e x e g ( x ) = x 2 − 4 .
Procurando uma função de ponto fixo adequada pode-se fazer:
Verificando as hipóteses i) e ii) do Teorema 2:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
2-26
Cálculo Numérico
Zeros reais de funções reais
n
xn
x n+1
2-27
xn+1 − x n
0
1
2
3
Portanto, x =
2.3.3Método de Newton, Newton-Raphson (ou Método das
Tangentes)
Este método é uma particularidade do método das aproximações sucessivas. A idéia é
construir uma função φ( x) para a qual exista um intervalo contendo o zero α , onde
φ ' ( x) < 1 . Esta construção é feita impondo φ' (α ) = 0 . Como φ' ( x) deve ser uma função
contínua, existe sempre uma vizinhança I de α onde max φ' ( α) <1.
x∈I
Obtenção da função φ( x) : A forma mais geral de x = φ( x) equivalente a f (x ) =0 é
dada por:
(Eq.7)
x = x + A( x) f ( x) = φ( x)
onde A( x) é uma função contínua tal que A( α) ≠ 0 . Escolhe-se A( x) de forma que
φ' (α ) = 0 . Derivando-se a (Eq.7), obtém-se φ' ( x) =1+ A( x) f ' ( x ) + A' ( x) f ( x) . Calculando
esta derivada no ponto α , obtém-se: φ ' (α ) =1+ A( α) f ' (α ) . Supondo que f ' (α)?0, para que
1
φ' (α ) = 0 , deve-se ter A(α) = −
. Assim, uma escolha satisfatória para A( x) será
f ' ( α)
portanto:
(Eq.8)
A( x) = −
1
, uma vez que x ≅ α .
f '( x)
Substituindo (Eq.8) em (Eq.7), tem-se:
(Eq.9)
φ( x ) = x −
f ( x)
f '( x )
Assim, o processo iterativo de Newton é definido por:
(Eq.10)
xn+1 = x n −
f ( xn )
, n = 0, 1, 2, …
f '( x n )
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Zeros reais de funções reais
OBS. 13:
A (Eq. 10) é válida mesmo que f ' ( α ) = 0, uma vez que x n ≠ α .
2-28
2.3.3.1 Interpretação Geométrica do Método de Newton
O ponto xn+1 é obtido traçando-se a tangente ao gráfico da função f (x ) no ponto
( xn , f ( x n )) . A intersecção da reta tangente com o eixo das abscissas fornece a nova
aproximação xn+1 . Esta interpretação justifica o nome de método das tangentes.
y
f (xn)
θ
xn +1
α
f (x)
x2
xn
x1 x0 x
[Fig. 22]: Interpretação Geométrica do Método de Newton.
tgθ = f ' ( xn ) =
f (x n )
f ( xn )
⇒ x n+1 = x n −
x n − x n+1
f ' (x n )
2.3.3.2 Convergência do Método de Newton
Teorema 3
que:
i)
ii)
iii)
Seja f : [a , b] → ℜ , duas vezes diferenciável, com f " (x ) contínua. Suponha
f (a ) ⋅ f (b ) < 0
f ' ( x ) ≠ 0, ∀x ∈ [ a, b]
f ' ' ( x) não troca de sinal em [a, b]
Então, a seqüência gerada pelas iterações do método de Newton-Raphson utilizando a
f ( xn )
f (x )
função φ( x ) = x −
que equivale a x n+1 = xn −
converge para o único zero α de
f ' ( xn )
f '( x )
f , isolado em [a, b] , se x0 ∈ [a , b] for escolhido convenientemente.
OBS. 14:
Para se escolher o ponto inicial x 0 , pode-se, por exemplo, fazer x 0 = a se
φ(a ) ∈ [a , b] ou x0 = b caso contrário.
2.3.3.3 Algoritmo do Método de Newton
Para encontrar uma solução para
aproximação inicial p 0 .
f (x ) =0, dada a derivada de
f (x ) e uma
• Dados de Entrada : Aproximação inicial p 0 , precisão ou tolerância (ε) e o número
máximo de iterações (ITMAX).
• Saída: Solução aproximada p ou mensagem de “solução não encontrada”.
PASSO 1
Faça i =1
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
PASSO 2:
Zeros reais de funções reais
Enquanto i ≤ ITMAX, execute os passos 3 – 6
PASSO 3
Faça p = p0 − f ( p 0 ) / f ' ( p0 ) (calcular pi )
PASSO 4
Se p − p 0 < ε então
Saída (p) (procedimento efetuado com sucesso)
FIM
PASSO 5
Faça i = i + 1
PASSO 5
Faça p 0 =p (atualize p 0 )
Passo 7:
Saída (solução não encontrada após ITMAX iterações)
FIM
OBS. 15:
•
•
•
Outros critérios de parada podem ser utilizados:
pn − pn −1 < ε
p n − p n−1
pn
<ε
f ( pn ) < ε
OBS. 16:
O Método de Newton irá falhar se para algum n, f ' ( p n−1 ) = 0.
Exercício 50
Encontrar a solução para a equação x = cos x com precisão ε = 10 −6 .
Resolução:
Pode-se construir uma tabela de valores para f ( x ) e analisar os sinais:
π
x
0
2
f ( x)
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
2-29
Cálculo Numérico
Zeros reais de funções reais
y
2-30
g (x) = x
1
π
α
3π
2
h(x )= cos x
π
2
2π
x
-1
[Fig. 23]: Os gráficos das funções
n
g (x ) = x e h (x) = cos x .
xn
xn+1 − xn
xn +1
0
1
2
Portanto, x =
2.3.4 Comparação entre os métodos
Nos exercícios seguintes, considerando cada método especificado, determine uma
aproximação para o zero da função.
Pelo método da Bissecção, determine uma aproximação para x ∈(1,2) da
Exercício 51
função f ( x )= e − x − cos x com aproximação ε1 = 10 −4 tal que ( b − a )/2< ε1 .
2
Resolução:
n
1
2
3
4
5
6
7
8
9
a
x
b
f (a ) f ( x ) f (b )
Universidade Tecnológica Federal do Paraná (UTFPR)
( b − a )/2
LAURO / NUNES
Cálculo Numérico
10
11
12
13
14
Logo, x =
Zeros reais de funções reais
2-31
Pelo método da Ponto Fixo ou Aproximações Sucessivas, determine uma
Exercício 52
aproximação para x ∈(1,2) da função f ( x )= e − x − cos x com aproximação ε1 = ε 2 =10 −4 tal
que | f ( x n )|< ε1 ou | xn+1 − x n |< ε 2 . Utilize x 0 =1,5.
2
Resolução:
xn
n
0
1
2
3
4
5
xn+1
| xn+1 − x n |
| f ( xn+1 )|
Parada
Logo, x =
Pelo método de Newton-Raphson, determine uma aproximação para x ∈(1,2)
Exercício 53
da função f ( x )= e − x − cos x com aproximação ε1 = ε 2 = 10 −4 tal que | f ( x n )|< ε1 ou
| xn+1 − x n |< ε 2 . Utilize x 0 =1,5.
2
Resolução:
n
0
1
xn
xn+1
| xn+1 − x n |
| f ( x n )|
Parada
Logo, x =
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Resolução de sistemas de equações lineares
3-32
3 Resolução de sistemas de equações
lineares
3.1 Introdução
Vários problemas, como cálculo de estruturas de redes elétricas e solução de equações
diferenciais, recorrem a resolução numérica de um sistema linear S n de n equações com n
incógnitas.
3.1.1 Forma Algébrica de Sn
(Eq.11)
 a11x1
a x

S n =  21 1
 M
 an1 x1
+ a12 x2
+ a 22 x2
M
+ an 2 x2
+ L +
+ L +
+ L +
a1n xn
a2n xn
M
ann x n
= b1
= b2
M
= bn
ou
n
(Eq.12)
S n = ∑ aij x j = bi , i =1, 2, …, n .
j =1
3.1.2 Forma Matricial de Sn
(Eq.13)
A⋅ x=b
 a11 a12 L a1n   x1   b1 
a
    
 21 a22 L a2n  ⋅  x2  = b2  .
 M
M
M  M M

    
an1 an2 L ann   xn  bn 
Onde:
• A ⇒ matriz dos coeficientes;
• x ⇒ vetor das incógnitas (ou vetor solução);
• b ⇒ vetor dos termos independentes.
3.1.3 Matriz Aumentada ou Matriz Completa do Sistema
 a11 a12 L a1n b1 
a

21 a22 L a2 n b2 

B =[ A : b ]=
.
 M
M
M
M


an1 an2 L ann bn 
3.1.4 Solução do Sistema
x =( x1 , x2 , …, xn ) T .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Resolução de sistemas de equações lineares
3-33
3.1.5 Classificação de um Sistema Linear
• COMPATÍVEL: apresenta soluções;
• INCOMPATÍVEL: caso contrário.
3.1.6 Classificação quanto ao Determinante de A
• det A ≠0 (SPD) ⇒ sistema linear possível e determinado (SOLUÇÃO ÚNICA);
• det A =0 (SPI) ou (SI): a matriz A é SINGULAR.
(SPI) ⇒ Sistema possível e indeterminado,
(SI) ⇒ Sistema impossível.
Se bi =0, i =1, 2, …, n , isto é, se b =0, o sistema é dito HOMOGÊNEO. Todo
sistema homogêneo é compatível, pois admite sempre a solução x =0. A solução é chamada
TRIVIAL.
OBS. 17:
3.2 Métodos diretos
São métodos que determinam a solução de um sistema linear com um número finito de
operações.
Definição:
Dois sistemas lineares são equivalentes quando possuem a mesma solução.
3.2.1 Método de Eliminação de Gauss
Com ( n −1) passos, o sistema linear A ⋅ x = b é transformado num sistema triangular
superior equivalente. Tome det A ≠0 como hipótese.
A ⋅ x = b ≈ U ⋅ x = c , o que se resolve por substituição.
[ A : b ] ≈ [U : c ]
 a11 a12 L a1n b1  u11 u12 L u1n
a
 
 21 a22 L a2n b2  ≈  0 u22 L u 2n
 M
M
M
M  M
M
M

 
0 L unn
an1 an2 L ann bn   0
Exercício 54
c1 
c2 
.
M

cn 
2 x1 + 3x2

Resolver o sistema S3 , com S3 = 4 x1 + 4 x2
2 x − 3x
 1
2
−
x3
=
5
− 3 x3 = 3 .
+ x3 = − 1
Resolução:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Resolução de sistemas de equações lineares
3-34
(0 )
• Etapa 1: em B0 , tome L(i 0) , com i =1,2,3, como as linhas de B0 e a11
como pivô e
calculam-se os multiplicadores m(i10) ( i =2,3).
• Etapa 2: Repete-se o processo para o próximo pivô, situado na diagonal da matriz B1 .
)
Em B1 , tome L(1i ) , com i =2,3 e a (1
22 como pivô.
• Método compacto para a TRIANGULAÇÃO U ⋅ x = c :
Linha
Multiplicador
m
Matriz Aumentada
B0 ⇒
(1)
2
3
-1
5
(0 )
(2) m21 =
4
4
-3
3
(3)
(0 )
m31
(2)
(3)
(1)
m32
=
2
-3
1
Transformação
-1
B1 ⇒
=
(3)
B2 ⇒
As linhas contendo os pivôs formam o sistema U ⋅ x = c .
Exercício 55
Resolver o sistema S4 com arredondamento em duas casas decimais, na matriz
aumentada.
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
 8,7 x1
24,5 x

1
S4 ⇒ A ⋅ x = b ⇒ 
52
,
3
x
1

 21,0 x1
+ 3,0 x 2
− 8,8 x 2
− 84,0 x 2
− 81,0 x 2
+
+
−
−
Resolução de sistemas de equações lineares
9,3x 3 + 11,0 x4 = 16,4
11,5 x 3 − 45,1x4 = − 49,7
.
23,5x 3 + 11, 4 x 4 = − 80,8
13,2 x3 + 21,5 x 4 = − 106,3
Resolução:
Linha
(1)
(2)
(3)
(4)
Multiplicador
B0 ⇒
(0 )
m21
(0 )
m31
(0 )
m41
m
Matriz Aumentada
8,70
3,00
9,30
11,00
16,40
=
24,50
-8,80
11,50
-45,10
-49,70
=
52,30
-84,00
-23,50
11,40
-80,80
=
21,00
-81,00
-13,20
21,50
-106,30
B1 ⇒
(2)
(3)
(1)
m32
=
(4)
(1)
=
m42
(3)
(4)
(2 )
m43
B2 ⇒
=
B3 ⇒
(4)
Então A ⋅ x = b ≈ U ⋅ x = c ⇒ [ A : b ] ≈ [ U : c ].
U ⋅ x=c ⇒
Logo: x =
3.2.1.1 Cálculo do Resíduo
Uma medida para avaliar a precisão dos cálculos é o resíduo, que é dado por:
(Eq.14)
r = b − Ax .
Exercício 56
Com base no exercício anterior, calcular o resíduo r do sistema A ⋅ x = b .
Resolução:
r = b − Ax .
r=
3.2.1.2 Algoritmo de Eliminação de Gauss
Seja o sistema A ⋅ x = b , com An×n , xn×1 e bn×1 .
Sempre supor que akk ≠0 na etapa k .
TRIANGULARIZAÇÃO: A ⋅ x = b ≈ U ⋅ x = c .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
3-35
Cálculo Numérico
Resolução de sistemas de equações lineares
3-36
Para k =1, 2, …, ( n −1)
Para i =( k +1), …, n
m=
aik
akk
aik =0
Para j =( k +1), …, n
aij = aij − m ∗ akj
bi = bi − m ∗ bk
FIM
FIM
FIM
RESOLUÇÃO DO SISTEMA U ⋅ x = c .
xn =
bn
ann
Para k =( n −1), …, 2, 1
s =0
Para j =( k +1), …, n
s = s + akj ∗ x j
FIM
xk =
bk − s
akk
FIM
3.2.2 Estratégia de Pivoteamento Completo
No momento de se calcular o multiplicador mik , se o pivô estiver próximo de zero, o
método pode ampliar os erros de arredondamento. Para se contornar estes problemas, escolhese como pivô MAX aij , com i , j =1, 2, …, n .
Dado A ⋅ x = b , tome B =[ A : b ].
 a11 a12

 a21 a22
 M
M
B=
a p1 a p 2
 M
M

 an1 an2
L
a1q
L
L
a2 q
L a2 n
M
a1n
M
L
a pq
M
L a pn
M
L
anq
L ann
b1 

b2 
M
.
bp 
M

bn 
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Resolução de sistemas de equações lineares
3-37
Seja a pq = MAX aij , (i , j =1, 2, …, n ) o pivô da linha p . Então, calcula-se o
multiplicador m(iq0) = −
aiq( 0)
a(pq0)
, em cada linha, ∀i ≠ p com i =1, 2, …, n . Assim, anulam-se os
elementos aij da coluna q através da operação:
L(1i ) ← m(iq0) ∗ L(p0) + L(i 0) .
Eliminando-se a linha pivotal p , repete-se o processo até que se obtenha L(i k ) com k
conjuntos de operações elementares aplicadas sobre B , onde k =1, 2, …, ( n −1).
Resolva S4 com arredondamento em duas casas decimais, utilizando
eliminação de Gauss com pivoteamento completo.
Exercício 57
 8,7 x1
24,5 x

1
S4 ⇒ A ⋅ x = b ⇒ 
52,3 x1
 21,0 x1
+ 3,0 x 2
− 8,8 x 2
− 84,0 x 2
− 81,0 x 2
+ 9,3x 3
+ 11,5 x 3
− 23,5x 3
− 13,2 x3
+ 11,0 x4
− 45,1x4
+ 11, 4 x 4
+ 21,5 x 4
= 16,4
= − 49,7
.
= − 80,8
= − 106,3
Resolução:
Linha
(0 )
(1) m12
=
(2)
(0 )
m22
=
B0 ⇒
(3)
(4)
(1)
(2)
Multiplicador
(0 )
m42
(1)
m14
=
m
8,70
Matriz Aumentada
3,00
9,30
11,00
16,40
24,50
-8,80
11,50
-45,10
-49,70
52,30
-84,00
-23,50
11,40
-80,80
21,00
-81,00
-13,20
21,50
-106,30
=
(4)
(1)
=
m44
(1)
(4)
(2 )
m11
=
B1 ⇒
B2 ⇒
B3 ⇒
(1)
Então A ⋅ x = b ≈ U ⋅ x = c ⇒ [ A : b ] ≈ [ U : c ].
U ⋅ x=c ⇒
3.2.3 Refinamento de Soluções
Seja x (0) a solução aproximada para A ⋅ x = b . Obtém-se a solução melhorada x (1)
aplicando-se a correção δ( 0) em x (0) .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Resolução de sistemas de equações lineares
3-38
x (1) = x (0) + δ( 0)
Se A ⋅ x (1) = b , então
A ⋅( x (0) + δ( 0) )= b
A ⋅ x (0) + A ⋅ δ( 0) = b
⇒ A ⋅ δ( 0) = b − A ⋅ x (0)
⇒ A ⋅ δ( 0) = r (0 ) . Assim, δ( 0) vem de [ A : r (0 ) ].
Obtido o δ( 0) , calcula-se x (1) = x (0) + δ( 0) .
Repete-se o processo para se obter x (2) , x (3) , …, x ( k ) , até que se tenha a precisão
desejada. Logo, obtém-se o refinamento de forma iterativa pela seguinte equação:
x (i ) = x ( i −1) + δ(i −1) , com i =1, 2, … k .
(Eq.15)
Considerando a resposta x do Exercício 55 , faça o refinamento de x até que
se obtenha o resíduo r (k ) =0, considerando precisão dupla ( 10−4 =0,0001), quatro casas
decimais.
Exercício 58
 8,7 x1
24,5 x

1
A⋅ x=b⇒
52,3 x1
 21,0 x1
+ 3,0 x 2
− 8,8 x 2
− 84,0 x 2
− 81,0 x 2
+ 9,3x 3
+ 11,5 x 3
− 23,5x 3
− 13,2 x3
+ 11,0 x4
− 45,1x4
+ 11, 4 x 4
+ 21,5 x 4
= 16,4
= − 49,7
= − 80,8
= − 106,3
x (0) = [1,01 2,01 − 1,01 1,00 ]T
r (0 ) = b − A ⋅ x (0) ⇒ r (0 ) = [− 0,024 − 0,042 0,082 0,468]
REFINAMENTO:
x ( k ) = x (k −1) + δ( k −1)
A ⋅ δ( k −1) = r (k −1) ⇒ [ A : r (k −1) ] ⇒ δ( k −1)
T
Resolução:
[ A : r (0 ) ] ⇒ δ( 0) ⇒ x (1) = x (0) + δ( 0)
• k =1
Linha
Multiplicador
(1)
B0 ⇒
m
Matriz Aumentada
8,7000
3,0000
9,3000
11,0000
-0,0240
(2)
(0 )
=
m21
24,5000
-8,8000
11,5000
-45,1000
-0,0420
(3)
(0 ) =
m31
52,3000
-84,0000
-23,5000
11,4000
0,0820
(4)
(0 ) =
m41
21,0000
-81,0000
-13,2000
21,5000
0,4680
B1 ⇒
(2)
(3)
(4)
)
m(1
32
)
m(1
42
=
=
B2 ⇒
(3)
(4)
(2 )
m43
=
B3 ⇒
Considerando 4 casas decimais:
(4)
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Resolução de sistemas de equações lineares
3-39
[ A : r (0 ) ] ≈
Então:
[ A : r (0 ) ] ⇒ δ( 0) ⇒
Como:
x (1) = x (0) + δ( 0) ⇒
r (1 ) = b − A ⋅ x (1) ⇒
• Logo,
3.3 Métodos iterativos
A solução x de um sistema de equações lineares A⋅ x = b pode ser obtido resolvendo,
de forma iterativa, o sistema equivalente da forma x = F ⋅ x + d , onde F é uma matriz n × n ,
x
e
d
vetores
n × 1 . Isto pode ser feito tomando
φ( x ) = F ⋅ x + d ,
x ( k +1) = φ( x ( k ) ) = F ⋅ x ( k ) + d , onde k =0, 1, …, M , e M é o número máximo de iterações e
x ( 0) é o vetor inicial.
3.3.1 Testes de parada
O processo iterativo x ( k +1) gera aproximações até que:
• máx x (i k +1) − xi( k ) ≤ ε , sendo ε a tolerância; ou
1≤i≤ n
• k > M , sendo M o número máximo de iterações.
3.3.2 Método de Gauss-Jacobi.
Adaptação de A⋅ x = b para x = F ⋅ x + d :
 a11x1
a x

A⋅ x = b ⇒  21 1
 M
 an1 x1

 x1

 x
x = F ⋅ x +d ⇒ 2
M

xn

+ a12 x2
+ a 22 x2
M
+ an 2 x2
=
=
=
+ L +
+ L +
+ L +
a1n xn
a2n xn
M
ann x n
= b1
= b2
M
= bn
b1 − ( a12 x 2 + a13x3 + a14 x4 + L + a1n xn )
a11
b2 − ( a21x1 + a 23x3 + a24 x 4 + L + a 2n x n )
a22
M
bn − ( an1x1 + an 2 x2 + a n3 x3 + L + an( n−1) x(n−1) )
ann
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Resolução de sistemas de equações lineares 3-40
OBS. 18:
Para o sistema x = F ⋅ x + d , é necessário que a ii ≠ 0, ∀i . Caso isto não ocorra,
o sistema A⋅ x = b deve ser reagrupado.
Assim, a fórmula recursiva x = F ⋅ x + d é dada na forma matricial por:

 0
 x1   a
 x  − 21
 2   a22
 x 3  =  a31
  −
 M   a33
 xn   M
 a
− n1
 ann
−
a12
a11
0
a32
a 33
M
an 2
−
a nn
−
a13
a11
a
− 23
a 22
−
0
M
a n3
−
a nn
a1n 
 b1 

a 
a11 
11
a2n   x1   b2 
 
L −
a22   x 2   a22 
a  ∗ x  +  b 
L − 3n   3   3 
a33   M   a33 
M  x   M 
  n   bn 
L
0 


 ann 

L −
ou ainda x ( k +1) = F ⋅ x ( k ) + d o que é equivalente a:
 ( k +1)
x1

 ( k +1)
x2

 M

 ( k +1)
xn

b1 − ( a12 x 2 + a13 x3 + a14 x4 + L + a1n x (nk ) )
a11
(k)
(k)
b2 − ( a21x1 + a 23x3 + a24 x (4k ) + L + a2n xn(k ) )
(k)
=
=
=
(k)
(k )
a22
M
(k )
(k )
(k )
(k)
bn − ( an1 x1 + an2 x2 + a n3 x3 + L + a n(n−1) x( n−1) )
ann
Resolva o sistema a seguir, utilizando o método de Gauss-Jacobi, com
= 0n×1 e ε = 10− 2 =0,01.
= 7
10 x1 + 2 x2 + x3

A⋅ x = b ⇒  x1 + 5 x2 + x3 = − 8 ⇒ x = F ⋅ x + d
 2 x + 3x + 10 x = 6
 1
2
3
Exercício 59
x
(0 )
Resolução:
F=
e d=
Neste caso a fórmula de recorrência fica:
x ( k +1) = F ⋅ x ( k ) + d ⇒
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Resolução de sistemas de equações lineares
k
x1( k )
x (2k )
x3( k )
max xi(k ) − xi(k −1)
0
1
2
3
4
5
6
0
0
0
-
3-41
1≤i ≤3
Com x ( 0) = [0 0 0]T e ε=0,01, o processo convergiu com ........... iterações para:
x=
3.3.2.1 Critério das linhas
Uma condição suficiente (mas não necessária) para garantir a convergência do método
de Gauss-Jacobi aplicado ao sistema A⋅ x = b , com a ii ≠ 0, ∀i , é
n
(Eq.16)
∑ aij
j =1
j ≠i
< aii , i = 1, 2, 3, … , n .
Neste caso, a matriz dos coeficientes das incógnitas A é dita estritamente diagonal
dominante.
Verificar se o
10 x1

que segue: A⋅ x = b ⇒  x1
 2x
 1
Exercício 60
critério das linhas é satisfeito no sistema de equações A⋅ x = b ,
+ 2 x2 + x3
= 7
+ 5 x2
+ 3x 2
+ x3
+ 10 x3
= −8
= 6
10 2 1 


Resolução:
A= 1 5 1 ⇒
 2 3 10 
Logo, a matriz dos coeficientes A é estritamente diagonal dominante, o que garante a
convergência do método de Gauss-Jacobi aplicado a este sistema com esta ordem de
equações e incógnitas.
Verificar se o critério
 x1 + 3x 2

que segue: A⋅ x = b ⇒ 5 x1 + 2 x2

6 x2

Exercício 61
das linhas é satisfeito no sistema de equações A⋅ x = b ,
+ x3 = − 2
+ 2 x3
+ 8 x3
= 3
= −6
1 3 1


Resolução:
A = 5 2 2 ⇒
0 6 8
Logo a matriz dos coeficientes A não é estritamente diagonal dominante. Isto significa
que não é garantida a convergência do método de Gauss-Jacobi aplicado a este sistema
com esta ordem de equações e incógnitas.
Mas permutando adequadamente as equações do sistema, obtém-se o sistema equivalente:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Resolução de sistemas de equações lineares
3-42
A=
Logo, esta nova matriz dos coeficientes A é estritamente diagonal dominante, o que
garante a convergência do método de Gauss-Jacobi aplicado a este sistema com esta nova
ordem de equações e incógnitas.
3.3.3 Método de Gauss-Seidel.
É semelhante ao método de Gauss-Jacobi, com a diferença de utilizar xi( k +1) , 1 ≤ i < p ,
para o cálculo de x (pk +1) . Desta forma, as equações recursivas ficam:
 ( k +1)
 x1

 ( k +1)
x2

 ( k +1)
 x3

 M

 x (nk +1)

Exercício 62
x
(0 )
= 0n×1
=
=
=
=
b1 − ( a12 x (2k ) + a13 x3(k ) + a14 x4(k ) + L + a1n x (nk ) )
a11
( k +1)
(k)
b2 − ( a21x1
+ a 23x3 + a 24x (4k ) + L + a2 n xn( k ) )
a22
( k +1)
( k +1)
b3 − ( a31x1
+ a32 x 2
+ a34 x (4k ) + L + a3n x (nk ) )
a33
M
(k +1)
(k +1)
bn − ( an1 x1
+ an 2 x2
+ an3 x3( k +1) + L + a n(n−1) x((kn−+11)) )
ann
Resolva o sistema a seguir, utilizando o método de Gauss-Seidel, com
e ε = 10− 2 =0,01.
10 x1

A⋅ x = b ⇒  x1
 2x
 1
+ 2 x2
+
+ 5 x2
+ 3x 2
+ x3
+ 10 x3
x3
=
7
= −8
= 6
Resolução:
Neste caso a fórmula de recorrência fica:
 ( k +1)
=
 x1

 ( k +1)
=
x2

 ( k +1)
=
 x3

Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Resolução de sistemas de equações lineares
k
x1( k )
x (2k )
x3( k )
max xi(k ) − xi(k −1)
0
1
2
3
4
0
0
0
-
3-43
1≤i ≤3
Com x ( 0) = [0 0 0]T e ε=0,01, o processo convergiu com ......... iterações para:
x=
3.3.4 Comparação entre os métodos
Exercício 63
x
(0 )
= 0n×1
Resolva o sistema A⋅ x = b , utilizando o método de Gauss-Jacobi, com
e ε=0,05.
5 x1

A⋅ x = b ⇒ 3 x1
3 x
 1
+
+
x2
+ 4 x2
+ 3x 2
x3
= 5
+ x3
+ 6 x3
= 6
= 0
Resolução:
F=
e d=
Neste caso a fórmula de recorrência fica:
 ( k +1)
 x1


( k +1)
(k)
x
= F ⋅ x + d ⇒  x (2k +1)

 ( k +1)
 x3

=
=
=
k
x1( k )
x (2k )
x3( k )
max xi(k ) − xi(k −1)
0
1
2
3
4
5
6
7
8
9
10
11
0
0
0
-
Universidade Tecnológica Federal do Paraná (UTFPR)
1≤i ≤3
LAURO / NUNES
Cálculo Numérico
Resolução de sistemas de equações lineares
3-44
Com x ( 0) = [0 0 0]T e ε=0,05, o processo convergiu com ......... iterações para:
x=
Resolva o sistema
= 0n×1 e ε=0,05.
5 x1 + x2 +

A⋅ x = b ⇒ 3 x1 + 4 x2 +
3 x + 3x +
 1
2
Exercício 64
x
(0 )
A⋅ x = b , utilizando o método de Gauss-Seidel, com
x3
= 5
x3
6 x3
= 6
= 0
Resolução:
Neste caso a fórmula de recorrência fica:
 ( k +1)
=
 x1

 ( k +1)
=
x2

 ( k +1)
=
 x3

k
x1( k )
x (2k )
x3( k )
max xi(k ) − xi(k −1)
0
1
2
3
0
0
0
-
1≤i ≤3
Com x ( 0) = [0 0 0]T e ε=0,05, o processo convergiu com ......... iterações para:
x=
3.3.5 Critério de Sassenfeld
Uma condição suficiente para garantir a convergência do método de Gauss-Seidel
aplicado ao sistema A⋅ x = b , com a ii ≠ 0, ∀i , é M <1, sendo M = max βi , onde:
1≤i≤ n
β1 =
n
1
⋅ ∑ a1 j
a11 j= 2
βi =
1
aii
OBS. 19:
 i −1
⋅ ∑ aij ⋅ β j +
 j =1
n

j =i +1

∑ aij  , i = 2, 3, … , n .
Se o critério das linhas é satisfeito, então o critério de Sassenfeld também será
satisfeito.
Verificar se o critério de Sassenfeld é
+ 0,5 x2 −
 x1
 0, 2 x
+
x2
−

1
A⋅ x = b , que segue: A⋅ x = b ⇒ 
− 0,1x1 − 0,7 x 2 +
 0,1x1 + 0,3 x2 +
Exercício 65
Universidade Tecnológica Federal do Paraná (UTFPR)
satisfeito
0,1x3 +
0,2 x 3 −
x3
+
0,2 x 3 +
no sistema
0,1x4 =
0,1x4 =
0,2 x4 =
x4
=
de equações
0, 2
− 2,6
1,0
− 2,5
LAURO / NUNES
Cálculo Numérico
Resolução de sistemas de equações lineares
 1
0,5 − 0,1 0,1 
 0, 2
1
− 0,2 − 0,1

A=
− 0,1 − 0,7
1
0,2 


0,3
0,2
1 
 0,1
Resolução:
β1 =
1
⋅ [ a12 + a13 + a14 ] =
a11
β2 =
1
⋅ [ a 21 ⋅ β1 + a23 + a 24 ] =
a 22
β3 =
1
⋅ [ a 31 ⋅ β1 + a 32 ⋅ β 2 + a 34 ] =
a33
β4 =
1
⋅ [ a 41 ⋅ β1 + a42 ⋅ β2 + a 43 ⋅ β3 ] =
a44
3-45
Então, M = max βi = max { ........ ; ........ ; ........ ; ........ } = .................... Logo o critério de
1≤i ≤4
Sassenfeld está satisfeito, o que garante a convergência do método de Gauss-Seidel
aplicado a este sistema.
Verificar se o critério de Sassenfeld é satisfeito no sistema de equações
2 x1 + x2 + 3 x3 = 9

A⋅ x = b , que segue: A⋅ x = b ⇒ 
− x2 + x3 = 1
x
+ 3 x3 = 3
 1
Exercício 66
Resolução:
Com esta disposição de linhas e colunas, tem-se que:
β1 =
1
⋅ [ a12 + a13 ] =
a11
β1 =
1
⋅ [ a12 + a13 ] =
a11
β2 =
1
⋅ [ a 21 ⋅ β1 + a 23 ] =
a 22
β3 =
1
⋅ [ a31 ⋅ β1 + a32 ⋅ β 2 ] =
a33
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Resolução de sistemas de equações lineares
Então, M = max βi =
1≤i ≤3
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
3-46
Cálculo Numérico
Interpolação
4-47
4 Interpolação
4.1 Interpolação polinomial
Uma função f ( x ) pode ser conhecida por um conjunto finito e discreto de n +1
pontos.
y
( x4 , y4 )
( x0 , y0 )
( x1 , y1 )
( x3 , y3 )
( x2 , y2 )
x0
x1
[Fig. 24]: Interpolação de
xi
yi
x0
y0
x1
x2
x3
x4
x5
y1
y2
y3
y4
y5
x2
x3
f (x)
( x5 , y5 )
P(x )
x4
x5
x
f ( x ) pelo polinômio P ( x ).
Para se INTERPOLAR os n +1 pontos obtidos da tabela, é utilizado um polinômio
Pn ( x ) de tal forma que:
Pn ( xi )= f ( xi ) para i =0, 1, …, n .
(Eq.17)
4.1.1 Existência e Unicidade do Polinômio Interpolador
Pn(x)
Existe um único polinômio Pn ( x ), de grau ≤ n , tal que: Pn ( xi )= f ( xi ), com
Teorema 4
i =0,1,…, n , desde que xi ≠ x j , i ≠ j .
n
Tome
Pn ( xi )= ∑ ak xik = f ( xi ) para
i =0,1,…, n .
Desenvolvendo
o
sistema
k =0
n
f ( xi )= ∑ ak xik ( i =0,1,…, n ), obtém-se:
k =0
 a0

 a0

M
a
 0
+
+
a1x0
a1x1
+ a2 x0
+ a2 x12
2
M
+
a1xn
+ L +
+ L +
M
+
a2 xn2
n
an x0
an x1n
=
=
M
+ L +
an xnn
f ( x0 )
f ( x1 )
(i = 0 )
( i = 1)
M
=
f ( xn )
Universidade Tecnológica Federal do Paraná (UTFPR)
(i = n )
LAURO / NUNES
Cálculo Numérico
Interpolação 4-48
Daí, retira-se a matriz dos coeficientes A para se calcular as incógnitas a0 , a1 ,…, a n .
1 x0

1 x1
A= 
M M

1 xn
2
n
x0 L x0 

x12 L x1n 
.
M
M

xn2 L xnn 
A é uma matriz de VANDERMONDE e, sendo xi com i =0,1,…, n , pontos distintos,
o det A ≠0. Assim o sistema admite solução única.
OBS. 20:
det A =( xn − xn−1 )∗( xn − xn −2 )∗…∗( xn − x0 )∗( xn−1 − xn −2 )∗( xn−1 − xn −3 )∗…∗( xn−1 − x0 )
∗…∗
∗( x3 − x2 )∗( x3 − x1 )∗( x3 − x0 )∗( x2 − x1 )∗( x2 − x0 )∗( x1 − x0 ) ⇒ det A = ∏ (xi − x j ) .
i> j
ENTÃO: O polinômio Pn ( x ) existe e é único.
4.1.2 Forma de Lagrange
Seja f uma função tabelada em ( n +1) pontos distintos x0 , x1 ,…, xn e seja Li ( x )
polinômios de Lagrange de grau n , onde Li é dado por:
(x − xj)
1 , se i = k
de tal forma que Li ( xk )= 
0 , se i ≠ k
j = 0 ( xi − x j )
n
Li ( x )= ∏
j≠i
Exercício 67
Determine Li ( xk ) para i =0,1,2, k =0,1,2 e n =2.
Resolução:
• i =0 ⇒ L0 ( x )=
( x − x1 )( x − x2 )
( x0 − x1 )( x0 − x2 )
k =0 ⇒ L0 ( x0 )= .......... .
k =1 ⇒ L0 ( x1 )= .......... .
k =2 ⇒ L0 ( x2 )= .......... .
• i =1 ⇒ L1 ( x )=
( x − x0 )( x − x2 )
( x1 − x0 )( x1 − x2 )
k =0 ⇒ L1 ( x0 )= .......... .
k =1 ⇒ L1 ( x1 )= .......... .
k =2 ⇒ L1 ( x2 )= .......... .
• i =2 ⇒ L2 ( x )=
( x − x0 )( x − x1)
( x2 − x0 )( x2 − x1 )
k =0 ⇒ L2 ( x0 )= .......... .
k =1 ⇒ L2 ( x1 )= .......... .
k =2 ⇒ L2 ( x2 )= .......... .
Para x = xk , com k =0,1,2,…, n , temos:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Interpolação
4-49
n
Pn ( xk )= ∑ yi Li ( xk )
i =0
⇒ i ≠ k ⇒ yi Li ( xk ) =0
123
=0
⇒ i = k ⇒ yi Li ( xi ) = yi
123
=1
A forma de Lagrange para o polinômio interpolador é:
n
(Eq.18)
Pn ( x )= ∑ yi Li ( x ) ou
i =0
n
Pn ( x )= ∑ yi
i =0
n
(x − xj)
∏ (x
j=0
j≠i
i
− xj)
Interpolar o ponto x =1,5 na tabela abaixo, empregando o polinômio
interpolador de Lagrange.
i
0
1
2
3
xi
−1
0
1
2
Exercício 68
yi
1
3
1
1
n =3 é o grau máximo de P3 ( x ).
Resolução:
3
P3 ( x )= ∑ yi Li ( x ) ⇒ P3 ( x )= .......... ⋅ L0 ( x )+ .......... ⋅ L1 ( x )+ .......... ⋅ L2 ( x )+ .......... ⋅ L3 ( x )
i =0
3
Li ( x )= ∏
(x − x j )
j = 0 ( xi
j≠i
L0 ( x )=
− xj)
( x − x1)( x − x2 )( x − x3 )
=
( x0 − x1 )( x0 − x2 )( x0 − x3 )
( x − x0 )( x − x2 )( x − x3 )
=
( x1 − x0 )( x1 − x2 )( x1 − x3 )
( x − x0 )( x − x1 )( x − x3 )
L2 ( x )=
=
( x2 − x0 )( x2 − x1 )( x2 − x3 )
( x − x0 )( x − x1 )( x − x2 )
L3 ( x )=
=
( x3 − x0 )( x3 − x1 )( x3 − x2 )
L1 ( x )=
Logo:
P3 ( x )=
⇒ P3 ( x )=
P3 (1,5)= P3 ( 32 )=
P3 (1,5)=
P3 (1,5)=
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
4-50
Interpolação
y
3
P3 (x )
2
1
3
8
-1
0
1
3
2
2
x
[Fig. 25]: Interpolação por Lagrange.
4.1.3 Forma de Newton
A forma de Newton para o polinômio Pn ( x ) que interpola f ( x ) em x0 , x1 ,…, xn ,
( n +1) pontos distintos é a seguinte:
Pn ( x )= f [ x0 ]+( x − x0 )⋅ f [ x0 , x1 ]+( x − x0 )⋅( x − x1 )⋅ f [ x0 , x1 , x2 ]+…
…+( x − x0 )⋅( x − x1 )⋅…⋅( x − xn−1 )⋅ f [ x0 , x1 ,…, xn ].
Onde
(Eq.19)
ORDEM
f [ x0 ]= f ( x0 )= y0
f [ x0 , x1 ]=
0
f [ x1 ] − f [ x0 ] f ( x1 ) − f ( x0 ) y1 − y0
=
=
x1 − x0
x1 − x0
x1 − x0
f [ x0 , x1 , x2 ]=
f [ x1, x2 ] − f [ x0 , x1]
x2 − x0
f [ x0 , x1 , x2 , x3 ]=
2
f [ x1 , x2 , x3 ] − f [ x0 , x1, x2 ]
x3 − x0
M
f [ x0 , x1 ,…, xn ]=
1
f [ x1, x2 ,L , xn ] − f [ x0 , x1,L, xn−1]
xn − x0
3
M
n
f [ x0 , x1 ,…, xn ] é a DIFERENÇA DIVIDIDA de ordem n da função f ( x ) sobre os n +1
pontos x0 , x1 ,…, xn .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Interpolação
4-51
4.1.3.1 Tabela Prática (DIFERENÇAS DIVIDIDAS)
x
ordem 0
x0
f [ x0 ]
ordem 1
ordem 2
… ordem n
ordem 3
f [ x0 , x1 ]
x1
f [ x1 ]
f [ x0 , x1 , x2 ]
f [ x1 , x2 ]
x2
f [ x0 , x1 , x2 , x3 ]
f [ x2 ]
f [ x1 , x2 , x3 ]
f [ x2 , x3 ]
x3
f [ x1 , x2 , x3 , x4 ]
f [ x2 , x3 , x4 ]
f [ x3 ]
f [ x3 , x4 ]
x4
O
f [ x0 ,…, xn ]
M
f [ x4 ]
M
N
f [ xn −3 , xn −2 , xn−1 , xn ]
M
M
O
f [ xn −2 , xn−1 , xn ]
M
f [ xn−1 , xn ]
xn
f [ xn ]
i
0
Interpolar o ponto x =1,5 na tabela abaixo, empregando a forma de Newton.
1
2
3
xi
−1
1
0
3
Exercício 69
yi
Resolução:
x
1
1
2
1
n =3 é o grau máximo de P3 ( x ). Tabela de diferenças divididas:
ordem 0
ordem 1
ordem 2
ordem 3
−1
0
1
2
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
4-52
Interpolação
P3 ( x )= f [ x0 ]+( x − x0 )⋅ f [ x0 , x1 ]+( x − x0 )⋅( x − x1 )⋅ f [ x0 , x1 , x2 ]+
+( x − x0 )⋅( x − x1 )⋅( x − x 2 )⋅ f [ x0 , x1 , x2 , x3 ]
P3 ( x )=
P3 ( x )=
P3 ( x )=
4.2 Estudo de erro na interpolação
Sejam x 0 < x1 < x 2 <…< x n , ( n +1) pontos. Seja f ( x ) com derivadas até ordem ( n +1)
para x pertencente ao intervalo [ x 0 , x n ].
Seja Pn ( x ) o polinômio interpolador de f ( x ) nos pontos x 0 , x1 , x 2 ,…, x n .
Então, em qualquer ponto x pertencente ao intervalo [ x 0 , x n ], o erro é dado por:
En ( x )= f ( x )− Pn ( x )
En ( x )=( x − x0 )⋅( x − x1 )⋅…⋅( x − x n )⋅
(Eq.20)
f ( n+1) ( ξ x )
( n + 1)!
onde ξ x ∈( x 0 , x n ).
Esta fórmula tem uso limitado, pois são raras as situações em que
conhecida e o ponto ξ x nunca é conhecido.
f ( n+1) ( x ) é
4.2.1 Estimativa para o Erro
Utilizando a (Eq.20), sendo f ( n+1) ( x ) contínua em I =[ x 0 , x n ], pode-se escrever:
| En ( x )|=| f ( x )− Pn ( x )|
n
| En ( x )|≤ ∏ ( x − xi ) ⋅
i =0
M n+1
, onde M n+1 = max f (n+1) ( x ) .
x∈I
( n + 1)!
Ao se construir a tabela de diferenças divididas até ordem n +1, pode-se usar o maior
M n+1
valor em módulo desta ordem como aproximação para
no intervalo I =[ x 0 , x n ].
( n + 1)!
Então:
| En ( x )|≈ ∏ ( x − xi ) ⋅ max ( Dd )
n
(Eq.21)
i =0
sendo Dd os valores da tabela de diferenças divididas de ordem ( n +1).
Exercício 70
Seja f ( x ) dada em forma de tabela de valores, como segue:
x
0,2
0,34
0,4
0,52
0,6
0,72
f (x)
0,16
0,22
0,27
0,29
0,32
0,37
a) Obter f (0,47) usando um polinômio de grau 2;
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
b) Dar uma estimativa para o erro.
4-53
Tabela de diferenças divididas:
Resolução:
x
Interpolação
ordem 0
ordem 1
ordem 2
ordem 3
0,2
0,34
0,4
0,52
0,6
0,72
Deve-se escolher 3 pontos próximos de 0,47 para a obtenção de P2 ( x ).
P2 ( x )= f [ x0 ]+( x − x0 )⋅ f [ x0 , x1 ]+( x − x0 )⋅( x − x1 )⋅ f [ x0 , x1 , x2 ]
P2 ( x )=
P2 ( x )=
a) P2 (0,47)=
.......... .......... ≈
f (0,47)
b) | En (0,47)|≈
| En (0,47)|≈ ..........
Prove a igualdade seguinte.
x − x1
x − x0
P1 ( x )= f ( x 0 )⋅
+ f ( x1 )⋅
= f [ x0 ]+( x − x0 )⋅ f [ x0 , x1 ]
x0 − x1
x1 − x0
Exercício 71
Resolução:
x
x0
ordem 0
ordem 1
f [ x0 ]= ..........
f [ x0 , x1 ]=
x1
f [ x1 ]= ..........
..........
⇒
P1 ( x )= f [ x0 ]+( x − x0 )⋅ f [ x0 , x1 ]
P1 ( x )=
⇔ P1 ( x )=
⇔ P1 ( x )=
⇔ P1 ( x )=
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
4-54
Interpolação
⇔ P1 ( x )=
⇔ P1 ( x )=
⇔ P1 ( x )=
⇔
P1 ( x )= f ( x 0 )⋅
x − x1
x − x0
+ f ( x1 )⋅
x0 − x1
x1 − x0
4.3 Interpolação inversa: casos existentes
O problema da interpolação inversa consiste em: dado y ∈( f ( x0 ), f ( xn )), obter x ,
tal que f ( x )= y .
São duas, as formas de se obter x . A primeira é encontrar x tal que Pn ( x )= y ; A
segunda é fazer a própria interpolação inversa, utilizando para isso, os valores de y .
4.3.1 Encontrar
x
tal que
Pn ( x )
Obter Pn ( x ) que interpola f ( x ) em x 0 , x1 , x 2 ,…, x n e em seguida encontrar x , tal
que f ( x )= y .
OBS. 21:
x obtido desta forma não permite se estimar o erro.
Exercício 72
Encontre x tal que f ( x )=2 pela tabela abaixo:
x
0,5
0,6
0,7
0,8
0,9
1,0
f (x)
1,65
1,82
2,01
2,23
2,46
2,72
Resolução:
Fazendo interpolação linear por x 0 =0,6 e x1 =0,7:
x − x1
x − x0
P1 ( x )= f ( x 0 )⋅
+ f ( x1 )⋅
x0 − x1
x1 − x0
P1 ( x )=
P1 ( x )=
P1 ( x )=
P1 ( x )=2
x = .......... .......... .......... .
4.3.2 Interpolação inversa
Se f ( x ) for inversível num intervalo contendo y , então x = f −1 ( y )= g ( y ).
Condição para a inversão de f ( x ): f é contínua e monótona crescente (decrescente)
num intervalo [ a , b ].
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Interpolação 4-55
Dado f ( x ) contínua em ( x 0 , x n ), então f ( x ) será admitida monótona crescente se
f ( x 0 ) < f ( x1 ) < … < f ( x n ) e monótona decrescente se f ( x 0 ) > f ( x1 ) > … > f ( x n ).
Respeitadas as condições dadas acima, será obtido o polinômio Pn ( y ) que interpola
g ( y )= f −1 ( y ) sobre [ y0 , yn ].
Exercício 73
Considere a tabela a seguir:
x
0
0,1
0,2
0,3
0,4
0,5
y =ex
1
1,1052
1,2214
1,3499
1,4918
1,6487
Obter x , tal que e x =1,3165, usando um processo de interpolação quadrática. Usar a forma de
Newton para obter P2 ( y ). Construir a tabela de diferenças divididas.
Resolução:
y
ordem 0
ordem 1
ordem 2
ordem 3
1
1,1052
1,2214
1,3499
1,4918
1,6487
P2 ( y )= g [ y0 ]+( y − y0 )⋅ g [ y0 , y1 ]+( y − y0 )⋅( y − y1 )⋅ g [ y0 , y1 , y2 ]
P2 ( y )=
P2 (1,3165)=
Assim, e .............................≈1,3165
Erro cometido:
| E2 ( y )| ≤ |( y − y0 )⋅( y − y1 )⋅( y − y2 )|⋅
Na calculadora = 1,316359.
M3
3!
| E2 (1,3165)| ≤
| E2 (1,3165)| ≤
⇒ M 3 = max g' ' ' ( y ) , y ∈[ y0 , y2 ].
M3
pode ser aproximado por .......... (tabela de diferenças divididas de ordem 3).
3!
| E2 (1,3165)| ≈ .......... .......... .......... .......... ⇒ | E2 ( y )| ≈ .......... .......... .......... .
1o Caso:
2o Caso: f ( x )= e x ⇒ g ( y )= f −1 ( y )= ln y
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
⇒
Interpolação
4-56
Logo: M 3 =
| E2 (1,3165)| ≤
4.4 Funções spline em interpolação
Considere f ( x )=
i =0,1,…, n .
1
2i
tabelada no intervalo [−1,1] nos pontos xi =−1+ , com
2
n
1 + 25 x
No gráfico abaixo, pode ser observada a função f ( x ) e o polinômio P10 ( x ) que
interpola o conjunto discreto de pontos para n =10.
−1,0
−0,8
−0,6
−0,4
−0,2
0
0,2
0,4
0,6
0,8
1,0
f ( x ) 0,038
0,059
0,1
0,2
0,5
1,0
0,5
0,2
0,1
0,059
0,038
x
3
2
y
1
1
2
P10(x )
f (x )
-1
- 12
[Fig. 26]: Gráfico do polinômio
0
1
2
1
x
P10( x ) interpolando f ( x ) .
Em certos casos, a aproximação por Pn ( x ) pode ser desastrosa. Uma alternativa é
interpolar f ( x ) em grupos de poucos pontos, obtendo-se polinômios de graus menores, e
impor condições para que a função de aproximação seja contínua e tenha derivadas contínuas
até uma certa ordem.
4.4.1 Função Spline
Considere a função f ( x ) tabelada nos pontos x 0 < x1 < x 2 <…< x n .
Uma função S p ( x ) é denominada SPLINE DE GRAU p com nós nos pontos xi ,
com i =0,1,…, n , se satisfaz as 3 seguintes condições:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Interpolação 4-57
1) Em cada subintervalo [ xi , xi +1 ], com i =0,1,…,( n −1), S p ( x ) é um polinômio de grau p
representado por si ( x ).
2) S p ( x ) é contínua e tem derivada contínua até ordem ( p −1) em [ a , b ].
3) S p ( xi )= f ( xi ), com i =0,1,…, n .
Nestes termos, S p ( x ) é denominada SPLINE INTERPOLANTE.
4.4.2 Spline linear interpolante
É representada por S1 ( x ) .
S1 ( x ) pode ser escrita em cada subintervalo [ xi −1 , xi ], com i =1,2,…, n como:
si ( x )= f ( xi −1 )
(Eq.22)
xi − x
x − xi −1
+ f ( xi )
, ∀x ∈[ xi −1 , xi ].
xi − xi −1
xi − xi −1
S1 ( x ) definida dessa forma satisfaz as condições 1) , 2) e 3) .
Achar a função spline linear que interpola a função f ( x ) tabelada a seguir.
Exercício 74
x0
x1
x2
x3
x
1
2
5
7
y = f (x)
1
2
3
2,5
3
2,5
2
y
s3 (x)
f (x )
s2 (x)
s1(x)
1
0
1
2
3
4
5
6
7
x
[Fig. 27]: Spline linear interpolando 4 pontos.
Pela definição, pode-se definir 3 splines lineares para os 4 pontos: s1 ( x ),
s2 ( x ) e s3 ( x ).
Resolução:
s1 ( x )= y0
x1 − x
x − x0
+ y1
x1 − x0
x1 − x0
⇒ s1 ( x )= .......... , x ∈[ .......... , .......... ].
s1 ( x )=
s2 ( x )= y1
x2 − x
x − x1
+ y2
x2 − x1
x2 − x1
s2 ( x )=
s3 ( x )= y2
⇒ s2 ( x )= .......... .......... .......... , x ∈[ .......... , .......... ].
x3 − x
x − x2
+ y3
x3 − x2
x3 − x2
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
s3 ( x )=
Interpolação
4-58
⇒ s3 ( x )= .......... .......... .......... , x ∈[ .......... , .......... ].
Então, no intervalo [ a , b ]=[1,7], a spline linear S1 ( x ) é dada por:
S1 ( x )=
4.4.3 Spline cúbica interpolante
É representada por S3 ( x ) .
A spline linear tem derivada primeira descontínua nos nós. A spline quadrática S2 ( x )
tem derivadas contínuas até ordem 1, portanto, pode ter picos ou troca abrupta de curvatura
nos nós.
A spline cúbica S3 ( x ) é mais utilizada por ter derivadas primeira e segunda contínuas,
que faz S3 ( x ) ser mais suave nos nós.
Suponha f ( x ) dada por xi , com i =0,1,…, n .
Tome S3 ( x ) como spline cúbica de f ( x ) nos nós xi , caso existam n polinômios de
grau 3 definidos em cada subintervalo k por sk ( x ), com k =1,2,…, n . Então a spline cúbica
S3 ( x ) deve satisfazer as 5 igualdades seguintes:
1) S3 ( x )= sk ( x ) para x ∈[ xk −1 , xk ], k =1,2,…, n .
2) S3 ( xi )= f ( xi ), com i =0,1,…, n .
3) sk ( xk )= sk +1 ( xk ), k =1,2,…,( n −1).
4) sk, ( xk )= sk, +1 ( xk ), k =1,2,…,( n −1).
5) sk,, ( xk )= sk,,+1 ( xk ), k =1,2,…,( n −1).
Em cada intervalo [ xk −1 , xk ], sk ( x ) será dada por:
(Eq.23)
sk ( x )= a k ( x − xk )3 + bk ( x − xk )2 + ck ( x − xk )+ d k , com k =1,2,…, n .
São 4 coeficientes para cada k à serem determinados.
Tome a notação h k = xk − xk −1 , para x = xk −1 .
Condição 1) : é satisfeita pela definição de sk ( x ).
Para a condição 2) , tem-se as equações:
(Eq.24)
d k = f ( xk )= sk ( xk ), k =1,2,…, n .
(Eq.25)
s1 ( x0 )= f ( x0 ) ⇒ − a1 h13 + b1 h12 − c1 h1 + d 1 = f ( x0 ), k =1.
Condição 3) para k =1,2,…,( n −1).
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
sk +1 ( xk )= f ( xk )
Interpolação
4-59
− a k +1 h 3k +1 + bk +1 h 2k +1 − ck +1 h k +1 + d k +1 = f ( xk ).
(Eq.26)
Para as condições 4) e 5) , tome as derivadas:
(Eq.27)
sk, ( x )=3 a k ( x − xk )2 +2 b k ( x − xk )+ c k .
(Eq.28)
sk,, ( x )=6 a k ( x − xk )+2 b k .
Para x = xk ⇒ sk,, ( xk )=2 b k . Assim, o coeficiente b k é dado por:
bk =
(Eq.29)
sk,, ( xk )
.
2
Para x = xk −1 ⇒ sk,, ( xk −1 )=−6 a k h k +2 b k .
ak =
2bk − sk,, ( xk −1 )
6h k
=
sk,, ( xk ) − s ,k, ( xk −1 )
.
6 hk
Impondo a condição 5) , sk,, ( xk −1 )= sk,,−1 ( xk −1 ), obtém-se:
sk,, ( xk ) − s ,k,−1( xk −1 )
ak =
, com s0,, ( x0 ) arbitrária.
6h k
(Eq.30)
Na obtenção de c k , utilizam- se as equações (Eq.25) e (Eq.26):
ck =
ck =
− f ( xk −1) − ak h3k + bk h2k + d k
hk
, d k = f ( xk )
f ( xk ) − f ( xk −1 )
−( a k h 2k − b k h k ), substituindo a k e b k obtém-se:
hk
f ( xk ) − f ( xk −1 )  sk,, ( xk ) − sk,,−1( xk −1)
sk,, ( xk )
ck =
− 
hk −
hk
2
hk
6




Daí, c k pode ser dado por:
(Eq.31)
ck =
,,
,,
f ( xk ) − f ( xk −1 ) 2s k ( xk ) ⋅ h k + sk −1( xk −1) ⋅ h k
+
.
6
hk
Na obtenção dos coeficientes, tome yk = f ( xk ) e g k = sk,, ( xk ).
(Eq.32)
ak =
gk − g k −1
6h k
(Eq.33)
bk =
gk
2
(Eq.34)
ck =
y k − yk −1 2hk gk + g k −1hk
+
hk
6
(Eq.35)
d k = yk .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Interpolação
4-60
Impondo a última condição 4) , sk, ( xk )= sk, +1 ( xk ), com k =1,2,…,( n −1), conclui-se
que:
Para x = xk ⇒ sk, ( xk )= c k , então:
c k =3 a k +1 h 2k +1 −2 b k +1 h k +1 + c k +1
⇒ c k +1 = c k −3 a k +1 h 2k +1 +2 b k +1 h k +1 .
Fazendo-se algumas substituições, através das equações (Eq.32), (Eq.33) e (Eq.34):
gk +1h k +1
y k +1 − yk 2h k +1g k +1 + g k hk +1 y k − yk −1 2hk gk + g k −1hk
g − gk
+
=
+
−3 k +1
h k +1 +2
6
6
2
hk +1
hk
6
Daí, chega-se a (Eq.36):
(Eq.36)
 y − yk yk − yk −1 
 , com k =1,2,…,( n −1).
h k g k −1 +2( h k + h k +1 ) g k + h k +1 gk +1 =6  k +1
−
 hk +1

h
k


A equação (Eq.36) é um sistema de equações lineares A g = b , onde k =1,2,…,( n −1).
A ordem do sistema é: A(n −1 )× (n +1 ) , g (n +1)×1 e b( n −1)×1 .
Pela variação de k , o sistema A g = b é indeterminado. Para se resolver o sistema, de
forma única, é necessário impor mais duas condições, apresentadas nas três alternativas a
seguir.
Spline Natural ⇒ nos extremos, S3 ( x0 ) é aproximadamente linear.
(1a)
S3" ( x0
S3" ( xn
(2a)
)= g 0 =0
)= g n =0
Nos extremos, S3 ( x ) é aproximadamente parábola.
g 0 = g1
g n = g n−1
(3a)
Nos extremos, é dada uma inclinação I 0 e I n para S3 ( x ).
S3' ( x0
)= I 0 ⇒ s1, ( x0 )= I 0 ⇒
3 a1 h12 −2 b1 h1 + c1 = I 0
S3' ( xn )= I n ⇒ s ,n ( xn )= I n ⇒
c n = In .
Nas alternativas (1a) e (2a), são eliminadas duas variáveis, g 0 e g n . Assim A g = b é
SPD, sendo que, o sistema é dado na ordem: A(n −1 )× (n −1) , g (n −1)×1 e b( n −1)×1 .
Na alternativa (3a), são acrescentadas duas equações. Assim A g = b é SPD, sendo
que, o sistema é dado na ordem: A(n +1)×(n +1) , g (n +1)×1 e b(n +1)×1 .
Encontrar uma aproximação para
interpolando a tabela:
Exercício 75
f (0,25) por spline cúbica natural,
x0
x1
x2
x3
x4
x
0
0,5
1,0
1,5
2,0
y = f (x)
3
1,8616
−0,5571
−4,1987
−9,0536
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Resolução:
n =4, logo, procura-se s1 ( x ), s 2 ( x ), s 3 ( x ) e s 4 ( x ).
Interpolação
4-61
Spline Natural ⇒ k =1,2,…,( n −1) ⇒ k =1,2,3 ⇒ Utilizando a (Eq.36), segue que:
 y − yk yk − yk −1 

(Eq.36) ⇒ h k g k −1 +2( h k + h k +1 ) g k + h k +1 gk +1 =6  k +1
−
 hk +1

h
k


Desenvolvendo o sistema A g = b :

 .......... g 0 + .......... g1 + ..........g 2 =
..........


 .......... g1 + .......... g 2 + .......... g3 =
..........


 .......... g 2 + ..........g 3 + ..........g 4 =
..........

g 0 = g 4 = .......... (Spline Natural).
Então,
 .......... .......... ..........  g1 
 ..........






.
A g = b ⇒  .......... .......... .......... ⋅  g 2  = .......... ∗  ..........

 .......... .......... ..........  g 3 
 ..........

Substituindo os valores:
 .......... .......... ..........  g1   ..........

 ..........

 ⋅ g  = 
 ⇒ g =
 .......... .......... ..........  2   ..........

 ..........
 .......... .......... ..........  g 3   ..........

 ..........
Forma geral de s i ( x ) ⇒ s i ( x )= ai ( x − xi )3 + bi ( x − xi )2 + ci ( x − xi )+ d i , com

.


i =1,2,3,4.
f (0,25) ≈ s1 (0,25)
g1 − g0
= ..........
6h
g
b1 = 1 = ..........
2
y − y0 2hg1 + g0h
c1 = 1
+
= ..........
h
6
d 1 = y1 = ..........
a1 =
Logo, s1 (0,25)=
⇒ s1 (0,25)=
⇒ a1 = ..........
⇒ b1 = ..........
⇒ c1 = ..........
⇒ d 1 = ..........
..........
.......... ..........
≈ f (0,25) .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Interpolação
4-62
Considerando os próximos 5 exercícios, encontrar uma aproximação para f ( x ) por
spline cúbica natural, interpolando a tabela:
x0
x1
x2
x3
x4
0
0,5
1,0
1,5
2,0
x
y = f (x)
3
1,8616
−0,5571 −4,1987 −9,0536
n =4, logo, procura-se s1 ( x ), s 2 ( x ), s 3 ( x ) e s 4 ( x ).
Do exercício anterior, a forma geral de s i ( x ) é dada por:
s i ( x )= ai ( x − xi )3 + bi ( x − xi )2 + ci ( x − xi )+ d i , com i =1,2,3,4.
Exercício 76
f (0,8).
Resolução:
f (0,8) ≈ s 2 (0,8)
g2 − g1
= ..........
6h
g
b 2 = 2 = ..........
2
y − y1 2hg 2 + g1h
c2= 2
+
= ..........
h
6
d 2 = y2 = ..........
a2=
Logo, s 2 (0,8)=
⇒ s 2 (0,8)=
Exercício 77
⇒ a 2 = ..........
⇒ b 2 = ..........
⇒ c 2 = ..........
⇒ d 2 = ..........
..........
.......... ..........
≈ f (0,8) .
f (1,1).
Resolução:
f (1,1) ≈ s 3 (1,1)
g3 − g 2
= ..........
6h
g
b3 = 3 = ..........
2
y − y2 2hg 3 + g 2h
c3 = 3
+
= ..........
h
6
d 3 = y3 = ..........
a3 =
⇒ a 3 = ..........
⇒ b3 = ..........
⇒ c 3 = ..........
⇒ d 3 = ..........
Logo, s 3 (1,1)=−0,7137(−0,4) −3,1260(−0,4)2 −8,6678(−0,4)−4,1987
3
⇒ s 3 (1,1)=
.......... ..........
≈ f (1,1) .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Exercício 78
Interpolação
4-63
f (1,2).
Resolução:
f (1,2) ≈ s 3 (1,2)
g3 − g 2
= ..........
6h
g
b3 = 3 = ..........
2
y − y2 2hg 3 + g 2h
c3 = 3
+
= ..........
h
6
d 3 = y3 = ..........
a3 =
Logo, s 3 (1,2)=
⇒ s 3 (1,2)=
Exercício 79
⇒ a 3 = ..........
⇒ b3 = ..........
⇒ c 3 = ..........
⇒ d 3 = ..........
..........
.......... ..........
≈ f (1,2) .
f (1,3).
Resolução:
f (1,3) ≈ s 3 (1,3)
g3 − g 2
= ..........
6h
g
b3 = 3 = ..........
2
y3 − y2 2hg 3 + g 2h
c3 =
+
= ..........
h
6
d 3 = y3 = ..........
a3 =
Logo, s 3 (1,3)=
⇒ s 3 (1,3)=
Exercício 80
⇒ a 3 = ..........
⇒ b3 = ..........
⇒ c 3 = ..........
⇒ d 3 = ..........
..........
.......... ..........
≈ f (1,3) .
f (1,7).
Resolução:
f (1,7) ≈ s 4 (1,7)
g4 − g3
= ..........
6h
g
b 4 = 4 = ..........
2
y − y3 2hg4 + g3h
c4= 4
+
= ..........
h
6
d 4 = y4 = ..........
a4=
Logo, s 4 (1,7)=
⇒ s 4 (1,7)=
⇒ a 4 = ..........
⇒ b 4 = ..........
⇒ c 4 = ..........
⇒ d 4 = ..........
..........
.......... ..........
≈ f (1,7) .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Ajuste de curvas pelo método dos mínimos quadrados
5-64
5 Ajuste de curvas pelo método dos
mínimos quadrados
5.1 Introdução
Uma forma de se trabalhar com uma função definida por uma tabela de valores é a
interpolacão. Contudo, a interpolação pode não ser aconselhável quando:
É preciso obter um valor aproximado da função em algum ponto fora do intervalo de
tabelamento (extrapolação).
Os valores tabelados são resultado de experimentos físicos, pois estes valores poderão conter
erros inerentes que, em geral, não são previsíveis.
Surge então a necessidade de se ajustar a estas funções tabeladas uma função que seja
uma “boa aproximação” para as mesmas e que nos permita “extrapolar” com certa margem de
segurança.
Assim, o objetivo deste processo é aproximar uma função f por outra função g ,
escolhida de uma família de funções em duas situações distintas:
Domínio discreto: quando a função f é dada por uma tabela de valores.
y
x
[Fig. 28]: Domínio discreto
Domínio contínuo: quando a função f é dada por sua forma analítica.
y
y = f (x)
x
a
b
[Fig. 29]: Domínio contínuo
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Ajuste de curvas pelo método dos mínimos quadrados
5-65
5.2 Caso Discreto
O problema do ajuste de curvas no caso em que se tem uma tabela de pontos:
x1
x2
x3
…
xm
f ( x1 )
f ( x2 )
f ( x3 )
…
f ( xm )
com x1 , x 2 , x 3 , … , x m ∈[ a , b ], consiste em: “escolhidas” n funções contínuas g 1 ( x ),
g 2 ( x ), g 3 ( x ), … , g n ( x ), contínuas em [a , b ], obter n constantes α 1 , α 2 , α 3 , … , α n
tais que a função g ( x )= α 1 g 1 ( x )+ α 2 g 2 ( x )+ α 3 g 3 ( x )+ … + α n g n ( x ) se aproxime ao
máximo de f ( x ).
Este modelo matemático é linear pois os coeficientes que devem ser determinados α 1 ,
α 2 , α 3 , … , α n aparecem linearmente, embora as funções g 1 ( x ), g 2 ( x ), g 3 ( x ), …
, g n ( x ) possam ser não lineares.
Surge então a primeira pergunta: Como escolher as funções contínuas g 1 ( x ), g 2 ( x ),
g 3 ( x ), … , g n ( x ) ?
Esta escolha pode ser feita observando o gráfico dos pontos tabelados (diagrama de
dispersão) ou baseando-se em fundamentos teóricos do experimento que forneceu a tabela.
Seja d k = f ( x k )− g ( x k ) o desvio em x k .
O método dos mínimos quadrados consiste em escolher os coeficientes α 1 , α 2 , α 3 ,
… , α n de tal forma que a soma dos quadrados dos desvio s seja mínima, isto é:
m
m
k =1
k =1
∑ d k2 = ∑ [ f ( xk ) − g (x k )]2
deve ser mínimo.
Assim, os coeficientes α 1 , α 2 , α 3 , … , α n que fazem com que g ( x ) se aproxime ao
máximo de f ( x ), são os que minimizam a função:
m
F ( α 1 , α 2 , α 3 ,…, α n )= ∑ [ f ( xk ) − g ( x k )] 2 =
k =1
m
∑ [ f ( xk ) − α1 g1 ( xk ) − α 2 g 2 (x k ) − α 3 g 3 ( x k ) − L − α n g n ( xk )]2 .
k =1
g
y
f ( x k)
dk
g( xk)
x
xk
[Fig. 30]: O método do mínimos quadrados
Para isto é necessário que:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Ajuste de curvas pelo método dos mínimos quadrados
∂F
(α 1 , α 2 , α 3 , L, α n ) =0, j =1, 2, 3, …, n , isto é:
∂α j
5-66
∂F
(α 1 , α 2 , α 3 ,L , α n ) =
∂α j
m
2· ∑ [ f ( x k ) − α1 g1 ( xk ) − α 2 g 2 ( xk ) − L − α n g n ( x k )] ⋅ [− g j ( x k )] =0, j =1, 2, 3, …, n
k =1
ou
m
∑ [ f (x k ) − α1g1 ( xk ) − α 2 g 2 ( xk ) − L − α n g n ( xk )] ⋅ [g j ( xk )] =0,
k =1
j =1, 2, 3, …, n
Assim, tem-se o seguinte sistema de n equações lineares com n incógnitas α 1 , α 2 ,
α3 , … , αn :
(Eq.37)
m
∑ [ f ( x k ) − α1g1 ( x k ) − α 2 g 2 ( xk ) − L − α n g n ( xk )] ⋅[ g1( xk )] = 0
 k =1
m
∑ [ f ( x k ) − α1g1 ( x k ) − α 2 g 2 ( xk ) − L − α n g n ( xk )] ⋅[ g 2 ( x k )] = 0
 k =1
M
M

m
∑ [ f ( x k ) − α1g1 ( x k ) − α 2 g 2 ( xk ) − L − α n g n ( xk )] ⋅[ g n ( x k )] = 0
 k =1
Que é equivalente a:
(Eq.38)
m
 m

m

g
(
x
)
⋅
g
(
x
)
⋅
α
+
L
+
g
(
x
)
⋅
g
(
x
)
⋅
α
=
 ∑ 1 k
∑ 1 k
1 k 
1
n k 
n
∑ g1( xk ) ⋅ f ( xk )
  k =1

 k =1

k =1
 m
m

m


g
(
x
)
⋅
g
(
x
)
⋅
α
+
L
+
g
(
x
)
⋅
g
(
x
)
⋅
α
=
∑ 2 k
∑ 2 k
∑ g 2 ( xk ) ⋅ f ( x k )
1 k 
1
n k 
n
  k =1

 k =1

k =1

M
M
 m
m

m

 ∑ g n ( xk ) ⋅ g1 ( x k ) ⋅ α1 + L + ∑ g n ( xk ) ⋅ g n ( xk )  ⋅ α n = ∑ g n ( x k ) ⋅ f ( xk )
  k =1
k =1

 k =1

As equações deste sistema linear são chamadas de equações normais.
Este sistema pode ser escrito na forma matricial A ⋅ α = b :
 a11α1
a
 21α1

M

 an1α

1
+
+
a12α 2
a 22α 2
+ L + a1n α n
+ L + a 2n α n
M
+
an 2 α 2
M
+ L +
ann α n
= b1
= b2
M
= bn
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Ajuste de curvas pelo método dos mínimos quadrados
m
m
k =1
k =1
5-67
onde A = ( aij ) tal que aij = ∑ g i ( xk ) ⋅ g j ( xk ) = ∑ g j ( xk ) ⋅ g i ( xk ) = a ji , ou seja, A é
uma matriz simétrica;
m
α = [α1 , α 2 ,L , α n ]T e b = [b1 , b2 ,L , bn ]T é tal que bi = ∑ gi ( x k ) ⋅ f ( x k ) .
k =1
m
Lembrando que, dados os vetores x e y ∈ ℜm o número real 〈 x, y〉 = ∑ xk ⋅ y k é
k =1
chamado de produto escalar de x por y , e usando esta notação no sistema normal A ⋅ α = b ,
tem-se: a ij = 〈 gi , g j 〉 e bi = 〈 g i , f 〉 onde:
g l é o vetor [ g l ( x1 ) g l ( x2 ) g l ( x3 ) L g l ( x m )]T e
f é o vetor [ f ( x1 ) f ( x 2 ) f ( x 3 ) L f ( x m )]T .
Desta forma o sistema na forma matricial fica:
(Eq.39)
 〈 g1 , g1〉 〈 g1 , g 2 〉 L 〈 g1 , g n 〉   α1   〈 g1 , f 〉 

〈 g , g 〉 〈 g , g 〉 L 〈 g , g 〉  α  
〈g 2 , f 〉
2
2
2
n   2

 2 1
⋅
=
 M
M
M   M   M 


   
〈 g n , g1 〉 〈 g n , g 2 〉 L 〈 g n , g n 〉  α n  〈 g n , f 〉 
Demonstra-se que, se as funções g 1 ( x ), g 2 ( x ), g 3 ( x ), … , g n ( x ) forem tais que os
vetores g1 , g 2 , g 3 ,L , g n , sejam linearmente independentes (LI), então det A ≠ 0 e o sistema
de equações é possível e determinado (SPD). Demonstra-se ainda que a solução única deste
sistema, α 1 , α 2 , α 3 , … , α n é o ponto em que a função F ( α 1 , α 2 , α 3 ,…, α n ) atinge seu
valor mínimo.
Se os vetores g1 , g 2 , g 3 ,L , g n , forem ortogonais entre si, isto é, se
〈 g i , g j 〉 = 0 se i ≠ j e 〈 g i , g j 〉 ≠ 0 se i = j , a matriz dos coeficientes A será uma matriz
diagonal, o que facilita a resolução do sistema A ⋅ α = b .
OBS. 22:
Exercício 81
(Regressão Linear) Ajustar os dados da tabela abaixo através de uma reta.
i
1
2
3
4
5
xi
1,3
3,4
5,1
6,8
8,0
f ( xi )
2,0
5,2
3,8
6,1
5,8
Fazendo g ( x ) = α1 ⋅ g1( x ) + α 2 ⋅ g 2 ( x) e considerando g1 (x ) = .......... e
g 2 (x ) = .......... , tem-se: g ( x ) = ..................................................
Assim, a reta que melhor se ajusta aos valores da tabela terá coeficientes α1 e α 2 , que
são solução do seguinte sistema na forma matricial:
 〈 g1, g1〉 〈 g1 , g 2 〉   α1   〈 g1, f 〉 

〈 g , g 〉 〈 g , g 〉  ⋅ α  = 
 2 1
2
2   2
〈 g 2 , f 〉 
g1 =[ .......... .......... .......... .......... .......... ]T
g 2 =[ .......... .......... .......... .......... .......... ]T
Resolução:
f =[ ..........
T
.......... .......... .......... .......... ]
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Ajuste de curvas pelo método dos mínimos quadrados 5-68
〈 g1 , g1 〉 = ..............................................................................................................................................................................
〈 g1 , g 2 〉 = ..............................................................................................................................................................................
〈 g 2 , g1 〉 = ..............................................................................................................................................................................
〈 g 2 , g 2 〉 = ..............................................................................................................................................................................
〈 g1 , f 〉 = ..............................................................................................................................................................................
〈 g2 , f 〉 =
Assim,
..............................................................................................................................................................................
Logo a equação da reta procurada é:
g ( x ) = .................................................
Exercício 82
Ajustar os dados da tabela através da parábola g ( x ) = x 2 :
i
1
2
3
4
5
6
7
8
9
10
11
xi
−1
−0,75
−0,6
−0,5
−0,3
0
0,2
0,4
0,5
0,7
1
f ( xi )
2,05
1,153
0,45
0,4
0,5
0
0,2
0,6
0,512
1,2
2,05
y
2
1
-1
1
x
[Fig. 31]: Diagrama de dispersão.
Fazendo
g ( x ) = α1 ⋅ g1 ( x)
e
considerando
g1 (x ) = x 2 ,
obtém-se
g ( x ) = .................... . Assim, para se obter a parábola que melhor se ajusta aos pontos da
tabela, será necessário encontrar α1 do sistema:
[〈 g1 , g1 〉 ] ⋅ [α 1 ] = 〈 f , g1 〉
g1 =[ .................... .................... .................... … .................... .................... ]T
Resolução:
[
f =[ ....................
]
.................... ....................
…
T
.................... .................... ]
〈 g1 , g1 〉 = ..............................................................................................................................................................................
..............................................................................................................................................................................
〈 g1 , f 〉 = ..............................................................................................................................................................................
..............................................................................................................................................................................
Assim, α1 = .................... .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Ajuste de curvas pelo método dos mínimos quadrados
Logo a equação da parábola procurada é: g ( x ) = .................................................
5-69
Ajustar os dados da tabela abaixo por um polinômio do segundo grau
g ( x ) = α1 + α2 ⋅ x + α 3 ⋅ x 2 .
Exercício 83
i
1
2
3
4
xi
−2
−1
1
2
f ( xi )
1
−3
1
9
Neste caso tem-se que: g1 (x ) = .......... , g 2 (x ) = .......... e g 3 ( x) = ..........
Resolução:
 〈 g1 , g1〉 〈 g1 , g 2 〉 〈 g1 , g 3 〉   α1   〈 g1 , f 〉 


   
〈 g 2 , g1〉 〈 g 2 , g 2 〉 〈 g 2 , g 3 〉  ⋅ α 2  = 〈 g 2 , f 〉 
〈 g 3 , g1 〉 〈 g 3 , g 2 〉 〈 g 3 , g 3 〉   α 3   〈 g 3 , f 〉 


T
g1 =[ .......... .......... .......... .......... ]
g 2 =[ ..........
g 3 =[ ..........
T
.......... .......... .......... ]
T
.......... .......... .......... ]
=[ .......... .......... .......... .......... ]T
f
〈 g1 , g1 〉 = ..............................................................................................................................................................................
〈 g1 , g 2 〉 = ..............................................................................................................................................................................
〈 g 2 , g1 〉 = ..............................................................................................................................................................................
〈 g1 , g 3 〉 = ..............................................................................................................................................................................
〈 g 3 , g1〉 = ..............................................................................................................................................................................
〈 g 2 , g 2 〉 = ..............................................................................................................................................................................
〈 g 2 , g 3 〉 = ..............................................................................................................................................................................
〈 g 3 , g 2 〉 = ..............................................................................................................................................................................
〈 g 3 , g 3 〉 = ..............................................................................................................................................................................
〈 g1 , f 〉 = ..............................................................................................................................................................................
〈 g2 , f 〉 =
..............................................................................................................................................................................
〈 g 3 , f 〉 = ..............................................................................................................................................................................
Assim,
Logo a equação da parábola procurada é:
g ( x ) = .................................................
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Ajuste de curvas pelo método dos mínimos quadrados
5-70
5.3 Caso Contínuo
No caso contínuo, o problema de ajuste de curvas consiste em: dada uma função
f ( x ) , contínua em [a , b ] e escolhidas as funções g1 ( x ), g 2 ( x ), g 3 ( x ),…, g n ( x ), todas
contínuas em [ a , b ], determinar constantes α 1 , α 2 , α 3 ,…, α n de modo que a função
g ( x )= α 1 g 1 ( x )+ α 2 g 2 ( x )+ α 3 g 3 ( x )+…+ α n g n ( x ) se aproxime ao máximo de f ( x ) no
intervalo [ a , b ].
Seguindo o critério dos mínimos quadrados para o conceito de proximidade entre
f ( x ) e g ( x ), os coeficientes α 1 , α 2 , α 3 ,…, α n a serem obtidos são tais que
b
∫a [ f (x ) − g (x )] dx seja o menor possível.
2
Para achar α tal que g ( x )≈ f ( x ), tome:
b
∫a [ f (x ) − g (x )] dx = F (α)= F ( α 1 , α 2 , α 3 ,…, α n ).
2
Encontram-se os pontos críticos de F (α):
∂F
(α)=0, j =1,2,…, n .
∂α j
b
b
a
a
Mas, F (α)= ∫ [ f ( x ) − g ( x )]2dx = ∫ [ f ( x )2 − 2 f ( x ) g ( x ) + g ( x )2 ]dx
b
b
b
a
a
a
⇒ F (α)= ∫ f ( x) 2dx −2 ∫ f ( x ) g ( x)dx + ∫ g ( x) 2dx .
Ao desenvolver
∂F
(α)=0, j =1,2,…, n , obtém-se:
∂α j
b
  bg 2 ( x) dx α
+ L +  ∫ g1( x ) g n ( x )dx  α n
1
  ∫a 1

 a

 b
b
 ∫ g2 ( x ) g1 ( x )dx  α1 + L +  ∫ g 2 ( x ) g n ( x )dx  α n

 a

  a

M
O
M
 b
 b g 2( x )dx  α
 ∫a gn ( x ) g1( x ) dx α1 + L +
 ∫a n
 n


b
=
∫a f (x )g1( x)dx
=
∫a f (x )g 2 ( x)dx .
b
M
=
b
∫a f (x )g n ( x )dx
Este é um sistema linear A α= b de ordem n .
b
A =( aij ) tal que aij = ∫ gi ( x) g j ( x )dx = a ji ⇒ aij = a ji .
a
A
é
SIMÉTRICA. α=( α 1 , α 2 , α 3 ,…, α n ) e
b =( b1 , b2 , b3 ,…, bn ),
tal
que
b
bi = ∫ f ( x) gi ( x) dx .
a
Usando a definição de produto escalar de duas funções p ( x ) e q ( x ) no intervalo
b
[ a , b ] por 〈 p , q 〉 = ∫ p( x) ⋅ q ( x) dx , o sistema A α= b fica:
a
(Eq.40)
A =( aij )= gi , g j e b =( bi )= f , gi .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Ajuste de curvas pelo método dos mínimos quadrados
5-71
Aproximar a função f ( x )=4 x 3 por um polinômio do primeiro grau, uma reta,
no intervalo [0,1].
Exercício 84
Resolução:
g ( x )= α 1 g 1 ( x )+ α 2 g 2 ( x )= ................................................., isto é, g1 ( x )= .......... e g 2 ( x )= .......... .
 a11 a12   α1   b1 
〈 g , g 〉
A α= b ⇒ 
⋅  =  ⇒  1 1

a21 a22  α 2  b2 
〈 g2 , g1〉
〈 g1 , g 2 〉   α1   〈 f , g1〉 
⋅
=
〈 g 2 , g 2 〉  α 2  〈 f , g2 〉 
a11 = 〈 g1, g1〉 = ..........
a12 = 〈 g1 , g 2 〉 = 〈 g 2 , g1 〉 = a21 = ..........
a22 = 〈 g 2 , g 2 〉 = ..........
b1 = 〈 f , g1〉 = ..........
b2 = 〈 f , g 2 〉 = ..........
A α= b ⇒
Logo:
g ( x )=
Exercício 85
.................................................
≈ f ( x )=4 x 3 em [0,1].
Aproximar a função f ( x )= e x no intervalo [0,1] por uma reta.
Resolução:
g ( x )= α 1 g 1 ( x )+ α 2 g 2 ( x )= ................................................., isto é, g1 ( x )= .......... e g 2 ( x )= .......... .
 a11 a12   α1   b1 
〈 g , g 〉
A α= b ⇒ 
⋅  =  ⇒  1 1

a21 a22  α 2  b2 
〈 g2 , g1〉
〈 g1 , g 2 〉   α1   〈 f , g1〉 
⋅
=
〈 g 2 , g 2 〉  α 2  〈 f , g2 〉 
a11 = 〈 g1, g1〉 = ..........
a12 = 〈 g1 , g 2 〉 = 〈 g 2 , g1 〉 = a21 = ..........
a22 = 〈 g 2 , g 2 〉 = ..........
b1 = 〈 f , g1〉 = ..........
b2 = 〈 f , g 2 〉 = ..........
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Ajuste de curvas pelo método dos mínimos quadrados
5-72
Usando o método de integração por partes em b2 : ∫ u ⋅ dv = u ⋅ v − ∫ v ⋅ du
g ( x )=
.................................................
5.4 Família de
Parâmetros
≈ f ( x )= e x em [0,1].
Funções
Não
Lineares
nos
Em alguns casos, a família de funções escolhidas pode ser não linear nos parâmetros,
m
isto é, g ( x ) não é da forma
∑ α k ⋅ gk (x ) . Nestes casos é preciso efetuar uma “linearização”,
k =1
através de transformações convenientes.
Exemplos:
1o ) f ( x )≈ α1 ⋅ e α2 x = g ( x )
ln f ( x )≈ ln α1 ⋅ e α2 x = ln α1 + α 2 ⋅ x = G ( x ).
Fazendo ln α1 = a1 e α 2 = a2 , tem-se: G ( x )= a1 + a 2 ⋅ x ,
Desta forma G ( x )≈ ln f ( x ), sendo que G ( x ) é linear nos parâmetros a1 e a 2 .
2o ) f ( x )≈
1
=g (x)
α1 + α 2 ⋅ x
1
≈ α1 + α 2 ⋅ x = G ( x ).
f ( x)
Fazendo α1 = a1 e α 2 = a2 , tem-se: G ( x )= a1 + a 2 ⋅ x ,
Desta forma G ( x )≈
1
, sendo que G ( x ) é linear nos parâmetros a1 e a 2 .
f ( x)
3o ) f ( x )≈ α1 + α 2 ⋅ x = g ( x )
f 2 ( x )≈ α1 + α 2 ⋅ x = G ( x ).
Fazendo α1 = a1 e α 2 = a2 , tem-se: G ( x )= a1 + a 2 ⋅ x ,
Desta forma G ( x )≈ f 2 ( x ), sendo que G ( x ) é linear nos parâmetros a1 e a 2 .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Ajuste de curvas pelo método dos mínimos quadrados 5-73
Exercício 86
Ajustar os dados da tabela que segue por uma função da forma
α2x
g ( x )= α1 ⋅ e .
x
0
1
2
f (x)
1
0,5
0,7
Desta forma, “linearizando” a função g ( x )= α1 ⋅ e α 2 x , como no primeiro
exemplo anterior, tem-se:
Resolução:
 〈 g1, g1〉 〈 g1, g 2 〉   a1   〈
〈 g , g 〉 〈 g , g 〉  ⋅ a  = 〈
 2 1
2
2   2 
, g1 〉 
, g 2 〉 
g1 =[ .......... .......... .......... ]T
g 2 =[ .......... .......... .......... ]T
T
.................... =[ .................... .................... .................... ]
〈 g1 , g1〉 = ..............................................................................................................................................................................
〈 g1 , g 2 〉 = ..............................................................................................................................................................................
〈 g 2 , g1〉 = 〈 g1 , g 2〉 = ..........
〈 g 2 , g 2 〉 = ..............................................................................................................................................................................
〈 .............. , g1 〉 = ..............................................................................................................................................................................
〈 ............. , g 2 〉 = ..............................................................................................................................................................................
⇒ g ( x )=
.................................................
≈ f ( x ).
Os parâmetros assim obtidos não são ótimos dentro do critério dos mínimos
quadrados, isto porque estamos ajustando o problema linearizado por mínimos quadrados e
não o problema original. Portanto, os parâmetros a1 e a 2 do exemplo, são os que ajustam a
função G ( x ) à função ln f ( x ), no sentido dos mínimos quadrados. Não se pode afirmar
que os parâmetros α1 e α 2 (obtidos de a1 e a 2 ) são os que ajustam g ( x )= α1 ⋅ e α 2 x à f ( x ),
dentro do critério dos mínimos quadrados.
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Integração Numérica
6-74
6 Integração Numérica
Se uma função f ( x ) é contínua em um intervalo [a , b ] e sua primitiva F ( x ) é
conhecida, então
(Eq.41)
b
∫a f (x )dx = F ( b )− F ( a )
onde F ' ( x )= f ( x ).
Por outro lado, nem sempre se tem F ( x ) e em alguns casos, a função a ser integrada
é dada por meio de tabela de pontos. Neste caso, torna-se necessária a utilização de métodos
numéricos.
A idéia básica da integração numérica é a substituição da função f ( x ) por um
polinômio que a aproxime no intervalo [ a , b ]. Assim o problema fica resolvido pela
integração de polinômios, o que é trivial de se fazer.
6.1 Fórmulas de Newton-Cotes
Neste caso, o polinômio que interpola f ( x ) o faz em pontos igualmente espaçados de
[ a , b ].
Fórmulas fechadas: x 0 = a , x n = b e
b
∫a
n
f ( x ) dx = ∑ Ai f ( x i ) , sendo Ai coeficientes
i =0
determinados de acordo com o grau do polinômio aproximador.
6.1.1 Regra dos Trapézios
f ( x1)
f (x )
p 1(x )
y
f ( x0)
0
a= x 0
b = x1
x
h= b - a , h= x1- x 0
[Fig. 32]: Regra dos trapézio.
A integral de f ( x ) no intervalo [ a , b ] é aproximada pela área de um trapézio.
(Eq.42)
b
∫a f (x )dx ≈
h
[ f ( x 0 )+ f ( x1 )] = I T
2
A aproximação de f ( x ) pela fórmula de Lagrange é p 1 ( x )= y 0 L0 ( x )+ y1 L1 ( x )
com L0 ( x )=
(Eq.43)
x − x1
x 0 − x1
p 1 ( x )=
e L1 ( x )=
x − x0
x1 − x 0
, logo:
x − x1
x − x0
f ( x 0 )+
f ( x1 )
−h
h
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Integração Numérica
6-75
6.1.1.1 Estimativa para o Erro
E n ( x )=( x − x 0 )⋅…⋅( x − x n )⋅
f ( x )= p 1 ( x )+ E ( x )
E ( x )=( x − x 0 )⋅( x − x1 )⋅
f ( n+1) ( ξx )
( n + 1)!
f " (ξ x )
, ξ x ∈( x 0 , x1 )
2
f " (ξ x )
f ( x )= p 1 ( x )+( x − x 0 )⋅( x − x1 )⋅
.
2
Integrando f ( x ):
x1
∫x
x1
x1
x0
x0
f ( x ) dx = ∫ p1 ( x ) dx + ∫
0
( x − x 0 )( x − x1 )
x = a
f " (ξ x )
dx , com  0
2
 x1 = b
b
I T = ∫ p1 ( x ) dx
a
b
ET = ∫ ( x − a )( x − b )
a
f " (ξ x )
dx
2
b
1 "
f ( c ) ∫ ( x − a )( x − b ) dx
a
2
1
( b − a )3
ET = f " ( c )
2
6
ET =
ET =
(Eq.44)
h3 "
f ( c ) com c ∈( a , b )
12
ou
| ET |≤
(Eq.45)
h3
max | f " ( x )|
12 x∈[ a,b]
OBS. 23:
b
b
∫a
 x 3 ax 2 bx 2  b 3 − 3ab 2 + 3a 2 b − a 3 ( b − a )3
( x − a x − b x + a b ) dx =  −
=
.
−
 =
6
2
2 
6
3
a
2
Exercício 87
Calcular
9
∫1
6x − 5 dx , usando a regra dos trapézios.
Resolução:
9
∫1
6x − 5 dx ≈ I T = .................... .
O erro cometido será, no máximo:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Logo, | ET |≤
Integração Numérica
6-76
.................... .
6.1.2 Regra dos Trapézios repetida
y
f (x )
0
a = x0
x1
x2
x3
x n-1
b = xn x
[Fig. 33]: Regra dos trapézios repetida
h = x1 − x 0 = x 2 − x1 = x 3 − x 2 = … = x n − x n −1
h=
b−a
, com n sendo o número de subdivisões do intervalo [ a , b ].
n
b
∫a f (x )dx ≈
Ai =
(Eq.46)
A1 + A2 + A3 +…+ An tal que Ai =área do trapézio i , com i =1,2,…, n .
h
[ f ( x i −1 )+ f ( x i )]
2
b
∫a
f ( x ) dx ≈
n −1
h
[ f ( x 0 )+ f ( x n )+2⋅ ∑ f ( x i ) ]
2
i =1
6.1.2.1 Estimativa para o Erro
(Eq.47)
| ETR |≤
Exercício 88
(b − a ) 3
max | f " ( x )|
2
x∈
[ a ,b ]
12n
Calcular
9
∫1
6x − 5 dx empregando o método dos trapézios com 8 repetições.
Determine uma aproximação para o erro cometido.
Resolução:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Integração Numérica
6-77
x 0 = ....... x1 = ....... x 2 = ....... x 3 = ....... x 4 = ....... x 5 = ....... x 6 = ....... x 7 = ....... x 8 = .......
x
f (x)
9
∫1
6x − 5 dx ≈
.................... .
Erro cometido será, no máximo:
| ETR | ≤
.................... .
Neste caso em particular, f ( x ) pode ser integrada de forma exata:
9
∫1
6x − 5 dx =
Exercício 89
.................... .
1
Seja I = ∫ e x dx . Calcule uma aproximação para I usando 10 subintervalos e a
0
regra dos trapézios repetida. Estimar o erro cometido.
Resolução:
1 x
∫0 e dx ≈ .................... .
Erro cometido será, no máximo:
| ETR | ≤
Exercício 90
.................... .
1
Seja I = ∫ e x dx . Qual o número mínimo de subdivisões, para a regra dos
0
trapézios repetida aplicada em I , de modo que o erro seja inferior a 10−3 ?
Resolução:
n = .................... .
6.1.3 Regra 1/3 de Simpson
É obtida aproximando-se a função f ( x ) da (Eq.41) por um por um polinômio
interpolador de 2o grau, p2 ( x ), que é dado pela fórmula de Lagrange:
p2 ( x )= L0 ( x ) f ( x0 )+ L1 ( x ) f ( x1 )+ L2 ( x ) f ( x2 )
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
2
Integração Numérica
tal que Li ( x )= ∏
(x − x j )
j = 0 ( xi
j≠i
− xj)
6-78
, com i =0,1,2.
y
f (x)
p2 (x )
f ( x0)
f ( x2)
f ( x1)
0
a= x 0
h
m= x 1
h
b = x2
x
[Fig. 34]: Regra 1/3 de Simpson.
x0 = a , x1 = m e x2 = b
a+b
2
m = x1 =
h=
b−a
2
x0 − x1 =− h ,
x0 − x2 =−2 h ,
x1 − x0 = h ,
x1 − x2 =− h ,
x2 − x0 =2 h ,
x2 − x1 = h .
p 2 ( x )=
b
∫a
( x − x1 )( x − x 2 )
( − h)( −2h )
f ( x 0 )+
x2
x2
x0
x0
( x − x 0 )( x − x 2 )
( h )( −h )
f ( x1 )+
( x − x 0 )( x − x1 )
( 2 h)( h )
f ( x2 )
f ( x ) dx = ∫ f ( x )dx ≈ ∫ p2 ( x) dx
f ( x0 ) x2
f ( x1) x2
f ( x2 )
(
x
−
x
)(
x
−
x
)
dx
−
( x − x0 )( x − x2 )dx +
1
2
2 ∫x
2
∫
x
0
0
2h
h
2h2
h
= [ f ( x 0 )+4 f ( x1 )+ f ( x2 )]. Logo:
3
=
(Eq.48)
b
∫a
x2
f ( x ) dx = ∫ f ( x )dx ≈
x0
x2
∫x (x − x0 )( x − x1)dx
0
h
[ f ( x 0 )+4 f ( x1 )+ f ( x2 )].
3
6.1.3.1 Estimativa para o Erro
x2
∫x
0
(Eq.49)
x2
x2
x0
x0
f ( x )dx = ∫ p2 ( x) dx + ∫ R2 ( x )dx
b
b
a
a
ES = ∫ R2 ( x) dx = ∫ ( x − x0 )( x − x1 )( x − x2 )
f ''' ( ξx )
dx
3!
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Integração Numérica
6-79
6.1.3.2 Mudança de Variável
z=
x − x0
⇒ x = hz + x 0
h
x0 − x 0
⇒ z =0
h
x − x 0 2h
x = x2 = b ⇒ z = 2
=
⇒ z =2
h
h
x = x0 = a ⇒ z =
dx
⇒ dx = h dz
h
h ⋅ f ''' (ξ z ) 2
ES =
∫0 hz ( hz − h )( hz −2 h ) dz
6
2
h4 ⋅ f ''' ( ξ z ) 2
h4 '''
ES =
z
(
z
−1)(
z
−2)
dz
=
f ( ξ z ) ∫ ( z 3 − 3z 2 + 2 z ) dz
∫
0
0
6
6
dz =
2
z4

h4 '''
h4 '''
ES =
f (ξ z ) ⋅  − z3 + z2  =
f ( ξ z ) ⋅0 = 0.
6
4
0 6
1442443
=0
Logo, ES =0. Isso quer dizer que ES não depende de R2 (resíduo de 2o grau).
Então:
f 4 (ξ x )
ES = ∫ R3 ( x ) dx = ∫ ( x − x0 )( x − x1)( x − x2 )( x − x3 )
dx
a
a
4!
2
h5 ⋅ f 4 ( ξ z ) 2
h5 4
ES =
z
(
z
−1)(
z
−2)(
z
−3)
dz
=
f ( ξ z ) ∫ ( z 4 − 6 z 3 + 11z 2 − 6 z ) dz
∫
0
0
24
24
1
4444244443
b
b
=−
ES = −
(Eq.50)
h5 4
f (ξ) com ( a ≤ξ≤ b ).
90
| ES | =
h5
⋅ max | f 4 ( x ) |
90 x∈[ a,b]
Considerando h =
(Eq.51)
4
15
| ES | ≤
b−a
( b − a) 5
⇒ h5 =
, tem-se:
2
32
( b − a) 5
⋅ max | f 4 ( x ) |.
2880 x∈[ a,b]
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
6-80
Integração Numérica
6.1.4 Regra 1/3 de Simpson repetida
y
f (x)
0
a = x0
x1
h
x2
x3
x4
x5
x6
xm -2
xm -1 b = x m x
[Fig. 35]: Regra 1/3 de Simpson repetida
Na figura, tome h =
b−a
⇒ h = xi − xi −1 ( i =1,2,… m ), para m =2 n ⇒ m é par.
2n
Aplica-se a regra de Simpson repetidas vezes no intervalo [ a , b ]=[ x0 , xm ].
x0 , x1 ,…, xm são pontos igualmente espaçados.
Então:
b
∫a
f ( x ) dx = ∫
xm
x0
f ( x ) dx
h
h
h
[ y 0 +4 y1 + y2 ]+ [ y2 +4 y 3 + y4 ]+…+ [ ym − 2 +4 y m −1 + ym ]
3
3
3
b
h
∫a f (x )dx ≈ 3 [ y 0 + ym +2( y2 + y4 +…+ ym − 2 )+4( y1 + y 3 +…+ y m −1 )]
≈
−1
2
2

h 
f ( x ) dx ≈  y 0 + ym +2 ∑ y2i +4 ∑ y2i −1  .
3 
i =1
i =1

m
(Eq.52)
b
∫a
m
6.1.4.1 Estimativa para o erro: ESR
ESR
h5
≤ n ⋅ ⋅ max | f 4 ( x ) |
90 x∈[ a,b]
(Eq.53)
ESR ≤ n ⋅
h5
⋅ max | f 4 ( x ) |
90 x∈[ a,b]
Considerando h =
(Eq.54)
ESR ≤
b−a
( b − a)
⇒ h5 =
, tem-se:
2n
32n 5
5
( b − a) 5
⋅ max | f 4 ( x ) |
4 x∈[ a ,b ]
2880n
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Exercício 91
Integração Numérica
6-81
1
Seja I = ∫ e xdx . Calcule uma aproximação para I usando a regra 1/3 de
0
Simpson com m =10. Estime o erro cometido.
Resolução:
1 x
∫0 e dx ≈ ................................................. .
Estimativa do erro:
ESR ≤ ................................................. .
Observe que ESR ≤ ................................................. e ETR ≤ ..................................................
Exercício 92
1
Seja I = ∫ e xdx . Para que valor de m teríamos erro inferior a 10−3 ?
0
Resolução:
m = .................... ⇒ Para um erro inferior a 10−3 seriam necessários
Obs: na regra dos trapézios com repetição são necessários
Exercício 93
..........
..........
subintervalos.
intervalos.
10
Seja I = ∫ log xdx . Aproxime I com a regra dos trapézios com 8 repetições.
6
Estime o erro cometido.
Resolução:
h = ................................................. ⇒ h = .................... .
i
0
1
2
3
xi
4
5
6
7
8
f ( xi )
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
10
∫6
Integração Numérica
6-82
log xdx ≈ ................................................. .
Estimativa do erro:
⇒ ETR ≤ ................................................. .
Exercício 94
Seja
10
I = ∫ log xdx . Aproxime
6
I
com a regra de Simpson com 8
subintervalos. Estime o erro cometido.
Resolução:
h = ................................................. ⇒ h = .................... . m = .......... e n = .......... .
i
0
1
2
3
4
5
6
xi
7
8
f ( xi )
10
∫6
log xdx ≈ ................................................. .
Estimativa do erro:
⇒ E SR ≤ ................................................. .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Solução numérica de equações diferenciais ordinárias
7-83
7 Solução numérica de equações
diferenciais ordinárias
7.1 Introdução
Se uma equação diferencial tem apenas uma variável independente, então ela é uma
equação diferencial ordinária.
EXEMPLOS:
dy
= x + y ; y , = x 2 + y 2 ; y ,, +(1− y 2 ) y , + y =0.
dx
Se uma equação diferencial envolve mais que uma variável independente, então ela é
equação diferencial parcial.
EXEMPLO:
∂ 2u ∂ 2u
+
=0, com u ≡ u ( x , y ).
∂x 2 ∂y 2
A ordem de uma equação diferencial é a mais alta ordem de derivação que aparece na
equação.
Se, dada uma equação diferencial de ordem n , a função, assim como suas derivadas
até ordem n −1, são especificadas em um mesmo ponto, então temos um problema de valor
inicial (PVI).
Se, em problemas envolvendo equações diferenciais ordinárias de ordem n , n ≥2, as
n condições fornecidas não são dadas todas num mesmo ponto, então temos um problema de
valor de contorno (PVC).
Exercício 95
Resolver a seguinte EDO:
dy
=− xy .
dx
Resolução:
⇒ y = .................... , para k ∈ℜ. Que representa uma família de curvas em ℜ2 .
Para a mesma EDO anterior, y , =− xy , resolva considerando uma condição
inicial y ( x0 )= y0 , com x0 =0 e y0 =1.
Exercício 96
Resolução:
y = .................... .
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Solução numérica de equações diferenciais ordinárias
7-84
7.2 Problema de valor inicial (PVI)
Uma equação diferencial de ordem n se apresenta da seguinte forma:
(Eq.55)
y (n ) = f ( x , y , y , , y ,, , y (3 ) , y (4 ) ,…, y (n −1 ) )
onde
y (l ) =
dl y
, l =1,2,…, n , x ∈[ a , b ] e y :[ a , b ]→ℜ.
dx l
Associadas a (Eq.55), podem existir condições cujo número coincide com a ordem da
EDO. Se tais condições se referem a um único valor x , tem-se um PROBLEMA DE VALOR
INICIAL – PVI. Caso contrário, tem-se um problema de valores de conterno.
7.2.1 Solução numérica de um PVI de primeira ordem
Toma-se m subintervalos de [a , b ], (m ≥1), e faz-se x j = x0 + j h onde h =
b−a
,
m
j =0,1,2,…, m , x j ∈[ a , b ].
I h ={ x0 , x1 ,…, xm } é denominado REDE ou MALHA de [ a , b ]. A solução numérica
ym ( x ) é a função linear por partes.
y (xm )
Solução
Exata
ym
y (x3 )
y (x0 ) = y0
y (x0)
y0
y (x2)
y (x1 )
x0
y1
x1
y2
x2
Solução
Numérica
y3
x3
xm -1
xm
[Fig. 36]: Gráfico da solução numérica de um PVI.
NOTAÇÃO: y ( x j )≈ y j significa que y j é aproximação para y ( x j ), x j ∈ I h .
NO GRÁFICO:
y ( x j ) ⇒ valor exato;
y j ⇒ valor aproximado;
j =1,2,…, m .
7.2.2 Método de Euler
Seja o PVI de primeira ordem definido por:
(Eq.56)
 y , = f ( x , y )

 y( x0 ) = y0 = η, sendo η um número dado.
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Solução numérica de equações diferenciais ordinárias 7-85
Para se aproximar y j para as soluções exatas y ( x j ), com j =1,2,…, m , procura-se
inicialmente y1 .
T
y1
e1
y (x1 )
y (x )
y (x0 ) = y0
x0
x1
[Fig. 37]: Gráfico do método de Euler.
Traça-se a tangente T à curva y ( x ) no ponto ( x0 , y ( x0 )), cuja equação é:
(Eq.57)
y ( x )= y ( x0 )+( x − x0 ) y , ( x0 ).
Fazendo x = x1 e lembrando que y ( x0 )= y0 , x1 − x0 = h , y , ( x0 )= f ( x0 , y ( x0 )) e
y1 ≈ y ( x1 ), tem-se:
(Eq.58)
y1 = y0 + h ⋅ f ( x0 , y ( x0 )).
7.2.2.1 Erro cometido
e1 = y1 − y ( x1 ).
7.2.2.2 Aproximação e erro de y j de forma geral
(Eq.59)
 y j +1 = y j + h ⋅ f ( x j , y j )
, com j =0,1,2,…, m −1.

e j +1 = y j +1 − y ( x j +1)
O método de Euler consiste em calcular RECURSIVAMENTE a seqüência { y j }
através das fórmulas:
(Eq.60)
( A )

( B )
Exercício 97
y0 = y( a ) = η
y j +1 = y j + h ⋅ f ( x j , y j )
, com j =0,1,2,…, m −1.
 y, = x − y + 2
Achar aproximações para a solução do PVI 
na malha de [0,1]
y
(
0
)
=
2

com h =0,1.
Resolução:
1− 0
→ m =10.
0,1
Usar (Eq.59) para j =0,1,2,…,9.
x0 =0, y0 =2, a =0, b =1, m =
j =0:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Solução numérica de equações diferenciais ordinárias
7-86
j =1:
TABELA:
j
xj
0
1
2
3
4
5
6
7
8
9
10
0
yj
y (xj)
2
2
2,004837
2,018731
2,040818
2,07032
2,106531
2,148812
2,196585
2,249329
2,30657
2,367879
y j − y ( x j )= e j
Na pratica, não se dispõe da solução exata y ( x j ) do PVI. Daí a necessidade de se
determinar uma expressão matemática para o erro. Usa-se a fórmula de Taylor para
desenvolver y ( x ), solução teórica do PVI, em torno de x0 :
(Eq.61)
x − x0 ,
( x − x0 ) 2 ,,
( x − x0 ) 3 ,,,
y ( x )= y ( x0 )+
y ( x0 )+
y ( x0 )+
y ( x0 )+…
2!
3!
1!
Fazendo x = x1 e lembrando que y ( x0 )= y0 , x1 − x0 = h , y , ( x0 )= f ( x0 , y ( x0 )) e
y1 = y ( x1 ), toma-se os dois primeiros termos da (Eq.61):
y1 = y0 + h ⋅ f ( x0 , y0 ). Generalizando-se, tem-se (Eq.59).
7.2.2.3 Erro local de truncamento - ELT
O erro no método de Euler quando se calcula y1 é obtido a partir do resto da fórmula
( x − x0 ) 2 ,,
h2 , ,
de Taylor, que é:
y (ξ), x0 <ξ< x1 , ou e1 =
y (ξ), para h = x1 − x0 .
2!
2!
Numa etapa j dos cálculos, tem-se:
(Eq.62)
ej =
h2 , ,
y (ξ), x j −1 <ξ< x j ,
2!
que é o ERRO LOCAL DE TRUNCAMENTO – ELT.
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Solução numérica de equações diferenciais ordinárias 7-87
Na prática, procura-se estabelecer COTAS ou ESTIMATIVAS para que se possa
conduzir o cálculo do erro com segurança.
Toma-se k = y ,, (ξ), constante, e h suficientemente pequeno para ser tomado como
parâmetro do ELT. Diz-se que ELT é da ordem de h 2 e se escreve ( h 2 ).
7.2.3 Métodos de Runge-Kutta
7.2.3.1 Métodos de passo simples
Um método para resolver o PVI (Eq.56) é de passo simples se a aproximação y j +1
depende apenas do resultado y j da etapa anterior.
Forma geral para métodos de passo simples:
(Eq.63)
y j +1 = y j + h φ( x j , y j ; h ), para j =0,1,2,…, m −1.
Onde φ é a função incremento e h o comprimento do passo.
Para o método de Euler, a função incremento é φ( x j , y j ; h )= f ( x j , y j ). Um
caso especial de Runge-Kutta.
OBS. 24:
7.2.3.2 Métodos com Derivadas
O método de Euler possui ordem um pois, foi obtido da fórmula de Taylor com
desenvolvimento até o termo em h .
Ao fazer o mesmo desenvolvimento até o termo em h 2 , obtém-se o método de passo
simples e ordem dois.
(Eq.64)
h2 , ,
y j +1 = y j + h y ( x j )+
y ( x j ), para j =0,1,2,…, m −1.
2!
,
7.2.3.3 ELT – Erro local de truncamento
(Eq.65)
h3 , , ,
e j +1 =
y (ξ), x j <ξ< x j +1 .
3!
Em (Eq.64), y , ( x j )= f ( x j , y j ).
OBS. 25:
y ,, ( x j )=? Regra da cadeia de f em relação a x j :
∂f
∂f
∂x
∂f
∂y
(xj, yj )= ( x j , yj ) (xj, y j) + ( x j , yj ) (xj, yj)
∂x4243 ∂x
∂x4243 ∂y
∂x4243
1
1
1
=1
= y ,, ( x j )
y ,, ( x j )=
= f (x j ,y j )
∂f
∂f
( x j , y j )+ ( x j , y j ) f ( x j , y j )
∂x
∂y
 y, = x − y + 2
Achar aproximações para a solução do PVI 
na malha [0,1] com
 y( 0) = 2
h =0,1 usando o método da (Eq.64).
Exercício 98
Resolução:
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Solução numérica de equações diferenciais ordinárias
1− 0
x0 =0, y0 =2, a =0, b =1, m =
→ m =10.
0,1
Usar (Eq.64) para j =0,1,…,9.
7-88
j =0:
j =1:
TABELA:
j
xj
0
1
2
3
4
5
6
7
8
9
10
0
yj
y (xj)
2
2
2,004837
2,018731
2,040818
2,07032
2,106531
2,148812
2,196585
2,249329
2,30657
2,367879
y j − y ( x j )= e j
7.2.4 Método de Euler Aprimorado (Método de RungeKutta de Segunda Ordem)
Retomando a (Eq.62): y j +1 = y j + h φ( x j , y j ; h ), para j =0,1,2,…, m −1.
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Solução numérica de equações diferenciais ordinárias
1
Fazendo-se φ( x j , y j ; h )= ( k1 + k 2 ) e substituindo na equação, tem-se:
2
(Eq.66)
y j +1 = y j +
7-89
h
( k1 + k 2 ), para j =0,1,2,…, m −1
2
onde k1 = f ( x j , y j ) e k 2 = f ( x j + h , y j + h k1 ).
 dy
 = − xy
Exercício 99
Achar aproximações para a solução do PVI  dx
na malha [0,1] com
 y( 0) = 1
h =0,5 usando o método de Euler Aprimorado.
Resolução:
j
xj
yj
0
1
2
0
1
k1
k2
y ( x j )= e − x
2
/2
| y j − y ( x j )|
7.2.5 Fórmulas de Runge-Kutta de Quarta Ordem
Estas fórmulas são normalmente as mais utilizadas.
(Eq.67)
y j +1 = y j +
h
( k1 +2 k 2 +2 k3 + k 4 ), para j =0,1,2,…, m −1
6
onde k1 = f ( x j , y j ),
h
h
k 2 = f ( x j + , y j + k1 ),
2
2
k3 = f ( x j +
h
h
, y j + k2 ) e
2
2
k 4 = f ( x j + h , y j + h k3 ).
7.2.5.1 Erro local de truncamento: ETL
(Eq.68)
ej =
h5 (5 )
y (ξ), x j −1 <ξ< x j .
5!
 dy
 = − xy
Exercício 100 Calcular a solução do PVI  dx
com h =0,1, no interior do intervalo
 y( 0) = 1
[0,1], pelo método de Runge-Kutta de quarta ordem.
Resolução:
y j +1 = y j +
h
( k1 +2 k 2 +2 k3 + k 4 ), para j =0,1,2,…,9.
6
k1 = .........................................................................................
k 2 = .........................................................................................
k3 = .........................................................................................
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES
Cálculo Numérico
Solução numérica de equações diferenciais ordinárias
k 4 = .........................................................................................
j xj
yj
k1
k2
k3
k4
0 0
1
2
3
4
5
6
7
8
9
10
7-90
1
 y, = x − y + 2
Achar aproximação para a solução do PVI 
na malha [0,1] com
y
(
0
)
=
2

h =0,1 usando o método de Runge-Kutta de segunda ordem (Euler aprimorado).
Exercício 101
Resolução:
x0 =0, y0 =2, a =0, b =1, m =
y j +1 = y j +
1− 0
→ m =10.
0,1
0,1
( k1 + k 2 ), para j =0,1,2,…,9.
2
k1 = ......................................................................................... e k 2 = .........................................................................................
j
xj
yj
0
1
2
3
4
5
6
7
8
9
10
0
2
k1
k2
Universidade Tecnológica Federal do Paraná (UTFPR)
LAURO / NUNES