Introdução aos Métodos Numéricos...Fórmula de Lagrange Observação importante: deixe a expressão deste polinômio neste formato. É trabalho inútil e sujeito a erros tentar

Introdução aos Métodos Numéricos

Instituto de Computação UFFDepartamento de Ciência da Computação

Otton Teixeira da Silveira Filho

Conteúdo temático

● Interpolação

Conteúdo específico

● Fórmula de Lagrange

● Método de Newton-Gregory

● Custo computacional

Fórmula de Lagrange


Se não necessitamos do polinômio na forma canônica podemos usar a fórmula de Lagrange

Parece simples... Mas...

pn(x)=L0 (x) y0+L1(x) y1+L2(x) y2+⋯+Ln(x) yn=∑i=0

n

Li(x) y i


Se não necessitamos do polinômio na forma canônica podemos usar a fórmula de Lagrange

Parece simples... Mas...

pn(x)=L0 (x) y0+L1(x) y1+L2(x) y2+⋯+Ln(x) y n=∑i=0

n

Li(x) y i

Li(x)=(x−x0)(x−x1)(x−x2)⋯(x−x i−1)(x−x i+1)⋯(x−xn)

(xi−x0)(xi−x1)(x i−x2)⋯(xi−xi−1)(xi−xi+1)⋯(x i−xn)


Mas é simples mesmo! Para facilitar o entendimento...

L0(x)=(x−x1)(x−x2)⋯(x−xn−1)(x−xn)

(x0−x1)(x0−x2)⋯(x0−xn−1)(x0−xn)L1(x )=

(x−x0)(x−x2)⋯(x−xn−1)(x−xn)

(x1−x0)(x1−x2)⋯(x1−xn−1)(x1−xn)

L2(x )=(x−x0)(x−x1)⋯(x−xn−1)(x−xn)

(x2−x0)(x2−x1)⋯(x2−xn−1)(x1−xn)

⋮

Ln(x)=(x−x0)(x−x1)(x−x2)⋯(x−xn−1)

(xn−x0)(xn−x1)(xn−x2)⋯(xn−xn−1)


Não confunda simples com fácil e fácil não é o oposto de trabalhoso


Um exemplo

Já sabemos que teremos um polinômio de grau 2

Escrevamos o caso geral para este grau

x f(x)

0 0

π/2 1

π/4 √2/2


Um exemplo

Vamos substituir os valores da tabela

L0(x)=(x−x1)(x−x2)

(x0−x1)(x0−x2); L1(x)=

(x−x0)(x−x2)

(x1−x0)(x1−x2); L2(x )=

(x−x0)(x−x1)

(x2−x0)(x2−x1)

p2(x )=L0(x ) y0+L1(x ) y1+L2(x ) y2


Um exemplo

ou

L0(x)=(x−π/2)(x−π/4)

(0−π/4)(0−π/2); L1(x)=

(x−0)(x−π/4)

(π /2−0)(π/2−π/4); L2(x )=

(x−0)(x−π/2)

(π /4−0)(π /4−π/2)

p2(x )=L0(x )0+L1(x )×1+L2(x)√22


Um exemplo

logo

L0(x)=8

π2(x−π/2)(x−π/4); L1(x )=

8

π2(x−0)(x−π/4); L2(x)=−

16

π2(x−0)(x−π/2)


Um exemplo

logo

L0(x)=8

π2(x−π/2)(x−π/4); L1(x )=

8

π2(x−0)(x−π/4); L2(x)=−

16

π2(x−0)(x−π/2)

p2(x )=8

π2 x (x−π/4)−

16

π2 x (x−π/2)

√22

=8

π2 [ x (x−π/4)−x (x−π/2)√2 ]


Calculemos pela fórmula x = π/6 e x = π/3

p2(x )=8

π2 [ x (x−π/4)−x (x−π/2)√2 ]



p2(x )=8

π2 [ x (x−π/4)−x (x−π/2)√2 ]

p2(π /6)=8π

2 [π/6(π /6−π/4)−π/6(π /6−π/2)√2 ]=8 [ 16 ( 1

6−

14 )−1

6 ( 16−

12 )√2 ]=0,517428

p2(π /3)=8π

2 [π/3(π/3−π/4)−π/3(π/3−π/2)√2] 8 [ 13 ( 1

3−

14 )− 1

3 ( 13−

12 )√2 ]=0,850761



Obviamente com os mesmos resultados que os obtidos anteriormente

p2(x )=8

π2 [ x (x−π/4)−x (x−π/2)√2 ]

p2(π /6)=8π

2 [π/6(π /6−π/4)−π/6(π /6−π/2)√2 ]=8 [ 16 ( 1

6−

14 )−1

6 ( 16−

12 )√2 ]=0,517428

p2(π /3)=8π

2 [π/3(π/3−π/4)−π/3(π/3−π/2)√2] 8 [ 13 ( 1

3−

14 )− 1

3 ( 13−

12 )√2 ]=0,850761


Observação importante:

deixe a expressão deste polinômio neste formato. É trabalho inútil e sujeito a erros tentar “desenvolver“ o polinômio até ficar na forma canônica.


Observação importante:

deixe a expressão deste polinômio neste formato. É trabalho inútil e sujeito a erros tentar “desenvolver“ o polinômio até ficar na forma canônica.

Principalmente com o estresse que surge durante a prova.


Também não precisamos calcular o polinômio como fizemos

Podemos substituir os valores dos pontos mais o valor do ponto no qual queremos calcular o polinômio

(x i , y i)


A Fórmula de Lagrange pode ser obtida teoricamente usando a regra de Cramer


Qual a vantagem desta fórmula?



● Não resolvemos um sistema




Qual o problema?




Qual o problema?

● O cálculo do polinômio não é tão simples e barato como no caso anterior


Qual o custo computacional de calcular um ponto?

É fácil: observe que temos (n+1) Li(x) e em cada um fazemos n multiplicações, 2(n-1) somas e uma divisão

● Temos um algoritmo O(n2)

Interpolação

Quais são os problemas dos dois algoritmos que conhecemos até agora?

Vários, mas estudaremos um:

● Não há um procedimento simples se quisermos adicionar mais pontos

Fórmula de Newton-Gregory


Construamos o método

Comecemos com um ponto (x0, y0)


Um ponto

Qual o polinômio que passa por este ponto?

● p0(x0) = y0

● É interpolador somente para este ponto

● Mesmo assim, queremos preservar esta informação

(x0, y0)


Dois pontos

Qual o polinômio que passaria por estes pontos?

Sugestão:

(x0, y0) ,(x1, y1)

p1(x)=p0(x)+A (x−x0)


Dois pontos

Qual o polinômio que passaria por estes pontos?

Sugestão:

● Preserva a informação obtida anteriormente e permite adequarmos este polinômio à interpolar o outro ponto

(x0, y0) ,(x1, y1)

p1(x)=p0(x)+A (x−x0)


Dois pontos

Ele passa pelo primeiro ponto, obriguemos que passe pelo outro

(x0, y0) ,(x1, y1)

p1(x1)=p0(x1)+A ( x1−x0)= y1


Dois pontos


(x0, y0) ,(x1, y1)

p1(x1)=p0(x1)+A ( x1−x0)= y1⇒ y0+A (x1−x0)= y1⇒ A=y1− y0

x1−x0


Dois pontos


Logo

(x0, y0) ,(x1, y1)

p1(x1)=p0(x1)+A ( x1−x0)= y1⇒ y0+A (x1−x0)= y1⇒ A=y1− y0

x1−x0

p1(x)=p0(x)+y1− y0

x1−x0

(x−x0) ; p0(x)= y0


Dois pontos

Observe que o polinômio obtido interpola estes dois pontos

(x0, y0) ,(x1, y1)

p1(x)=p0(x)+y1− y0

x1−x0

(x−x0) ; p0(x)= y0


Dois pontos

Observe que o polinômio obtido interpola estes dois pontos

Embora as coisas estejam indo bem, faremos uma simplificação

(x0, y0) ,(x1, y1)

p1(x)=p0(x)+y1− y0

x1−x0

(x−x0) ; p0(x)= y0


Antes de adicionar mais pontos, imporemos as seguintes restrições:

● Os pontos estarão ordenados em relação a x




● Os pontos estarão igualmente espaçados, ou seja,




● Os pontos estarão igualmente espaçados, ou seja,

xi+1−x i=h;∀ i


Dois pontos

Assim o polinômio obtido será escrito como

(x0, y0) ,(x1, y1)

p1(x)=p0(x)+y1− y0

h(x−x0) ; p0(x)= y0


Três pontos

Preservando o que obtivemos, vamos experimentar o polinômio

que funciona como interpolador dos dois primeiros pontos. Obriguemos que ele interpole o terceiro ponto.

(x0, y0) ,(x1, y1),( x2, y2)

p2(x)=p1(x)+B (x−x0)(x−x1)


Três pontos (x0, y0) ,(x1, y1),( x2, y2)

p2(x2)=p1(x2)+B (x2−x0)(x2−x1)= y2⇒ p0(x2)+y1− y0

h(x2−x0)+B (2h)(h)= y2


Três pontos (x0, y0) ,(x1, y1),( x2, y2)

p2(x2)=p1(x2)+B (x2−x0)(x2−x1)= y2⇒ p0(x2)+y1− y0

h(x2−x0)+B (2h)(h)= y2

y0+y1− y0

h2h+2h2 B= y2⇒B=

y2−2 y1+ y0

2h2


Três pontos

Então

(x0, y0) ,(x1, y1),( x2, y2)

p2(x2)=p1(x2)+B (x2−x0)(x2−x1)= y2⇒ p0(x2)+y1− y0

h(x2−x0)+B (2h)(h)= y2

y0+y1− y0

h2h+2h2 B= y2⇒B=

y2−2 y1+ y0

2h2

p2(x)=p1(x)+y2−2 y1+ y0

2h2 (x−x0)(x−x1)


Quatro pontos

Já parece natural experimentar o polinômio

que funciona como interpolador dos três primeiros pontos. Obriguemos que ele interpole o quarto ponto.

(x0, y0) ,(x1, y1),(x2, y2) ,(x3, y3)

p3(x)=p2(x)+C (x−x0)(x−x1)(x−x2)


Quatro pontos

logo

(x0, y0) ,(x1, y1),(x2, y2) ,(x3, y3)

p3(x3)=p2(x3)+C (x3−x0)(x3−x1)(x3−x2)= y3

p1(x3)+y2−2 y1+ y0

2h2(x3−x0)(x3−x1)+C (3h)(2h)(h)= y3

p0(x3)+y1− y0

h(x3−x0)+

y2−2 y1+ y0

2h26h2+C 6h3= y3

y0+3( y1− y0)+3( y2−2 y1+ y0)+6Ch3= y3⇒C=

y3−3 y2+3 y1− y0

6h3


Quatro pontos

que interpola os pontos dados

(x0, y0) ,(x1, y1),(x2, y2) ,(x3, y3)

p3(x)= p2( x)+y3−3 y2+3 y1− y0

6h3(x−x0)(x−x1)(x−x2)


Observe que temos padrões surgindo

p3(x)= p2( x)+y3−3 y2+3 y1− y0

6h3(x−x0)(x−x1)(x−x2)

p2(x )=p1(x)+y2−2 y1+ y0

2h2 (x−x0)(x−x1)

p1(x )= p0(x )+y1− y0

h(x−x0)

p0( x)= y0



Temos algo como

O “?“ é algo que parece exigir mais análise

pi(x )= pi−1(x)+?

i!hi(x−x0)(x−x1)(x−x2)⋯(x−xi−1)

p0( x)= y0

Operador de diferenças progressivas

Vamos definir o operador

que é linear. Assim, observe que

Δ f (x)=f (x+δ x)−f (x)

Δ2 f (x)=ΔΔ f (x)=Δ [ f (x+δ x)−f (x)]=Δ f (x+δ x)−Δ f ( x)




Δ f (x)=f (x+δ x)−f (x)


Δ2 f (x)= f (x+δ x+δ x)−f (x+δ x)− [ f (x+δ x)− f (x) ]




Δ f (x)=f (x+δ x)−f (x)


Δ2 f (x)= f (x+δ x+δ x)−f (x+δ x)− [ f (x+δ x)− f (x) ]

Δ2 f (x)= f (x+2δ x)−2 f (x+δ x)+f (x)


Da mesma forma, podemos escrever

Δ3 f (x)=ΔΔ2 f (X )=Δ [ f (x+2δ x )−2 f (x+δ x )+ f (x )]Δ

3 f (x )=Δ f (x+2δ x)−2Δ f ( x+δ x)+Δ f (x )

Δ3 f (x)=f (x+3δ x )−f (x+2δ x )−2 [ f (x+2δ x)−f (x+δ x) ]+ f (x+δ x )−f ( x)


Da mesma forma, podemos escrever

Então

Δ3 f (x)=ΔΔ2 f (X )=Δ [ f (x+2δ x )−2 f (x+δ x )+ f (x )]Δ

3 f (x )=Δ f (x+2δ x)−2Δ f ( x+δ x)+Δ f (x )

Δ3 f (x)=f (x+3δ x )−f (x+2δ x )−2 [ f (x+2δ x)−f (x+δ x) ]+ f (x+δ x )−f ( x)

Δ3 f (x)=f (x+3δ x )−3 f (x+2δ x)+3 f (x+δ x)−f (x )


Observe que podemos continuar a calcular as aplicações múltiplas do operador de diferenças progressivas o quanto quisermos.

No entanto, pararemos por aqui pois já dá para perceber o que queremos


Vamos substituir nas expressões encontradas

Δ f (x0)=f (x0+h)−f (x0)=f (x1)−f (x0)= y1− y0

x⇒ x0 ;δ x⇒h



Δ2 f (x0)=f (x0+2h)−2 f (x0+h)+ f (x0)=f (x2)−2 f (x1)+ f (x0)= y2−2 y1+ y0

Δ f (x0)=f (x0+h)−f (x0)=f (x1)−f (x0)= y1− y0

x⇒ x0 ;δ x⇒h



onde usamos o fato de trabalharmos com um polinômio interpolador.

Δ3 f ( x0)=f (x0+3h)−3 f (x0+2h)+3 f (x0+h)−f (x0)=f (x3)−3 f (x2)+3 f (x1)+ f (x0)

Δ3 f (x0)= y3−3 y2+3 y1+ y0

Δ2 f (x0)=f (x0+2h)−2 f (x0+h)+ f (x0)=f (x2)−2 f (x1)+ f (x0)= y2−2 y1+ y0

Δ f (x0)=f (x0+h)−f (x0)=f (x1)−f (x0)= y1− y0

x⇒ x0 ;δ x⇒h


Observando o que obtemos

p3(x)=p2( x)+y3−3 y2+3 y1− y0

6h3(x−x0)(x−x1)(x−x2)

p2(x )=p1(x)+y2−2 y1+ y0

2h2 (x−x0)(x−x1)

p1(x )= p0(x )+y1− y0

h(x−x0)

p0( x)= y0

pi(x )= pi−1(x)+?

i!hi(x−x0)(x−x1)(x−x2)⋯(x−xi−1)

⋮


vemos que o fator “?“ é exatamente dado pelos valores que encontramos aqui com as diferenças progressivas.

Assim, poderemos escrever



pi(x )=pi−1(x)+Δ

i f (x0)

i !hi(x−x0)(x−x1)(x−x2)⋯(x−xi−1)

p0( x)= y0


A fórmula parece simples mas o cálculo das diferenças parece um problema.


A fórmula parece simples mas o cálculo das diferenças parece um problema.

De fato, não é...


Observe a tabela abaixo

x y Δf(x) Δ2f(x) Δ3f(x) ... Δnf(x)

x0

y0

y1-y

0y

2-2y

1+y

0y

3-3y

2+3y

1 -y

0.

x1

y1

y2-y

1y

3-2y

2+y

1.

x2

y2

y3-y

2. .

x3

y3

. . .

. . . . .

xn-2

yn-2

yn-1

-yn-2

yn-2y

n-1+y

n-2

xn-1

yn-1

yn-y

n-1

xn

yn


● Colocamos duas colunas em uma tabela contendo os valores de x e y

● Criamos uma nova coluna subtraindo os valores de y progressivamente. Obtemos a primeira diferença

● Criamos outra coluna subtraindo os valores da primeira diferença progressivamente. Obtemos a segunda diferença

● Continuamos sucessivamente até os valores se reduzirem à última diferença


Para ficar mais claro, façamos um exemplo


Experimentando...

Três pontos: x f(x)

0 0

π/2 1

π/4 √2/2


Experimentando...

Três pontos:

Podemos usar os pontos como estão?

x f(x)

0 0

π/2 1

π/4 √2/2


Experimentando...

Três pontos:

Podemos usar os pontos como estão?

Certamente não: os pontos são igualmente espaçados mas não estão ordenados

x f(x)

0 0

π/2 1

π/4 √2/2


Ordenando: h = π/4

Tabela

x y Δf(x) Δ2f(x)

0 0

π/4 √2/2

π/2 1


Ordenando: h = π/4

Tabela

Inicialmente temos

x y Δf(x) Δ2f(x)

0 0 √2/2 1-√2

π/4 √2/2 1-√2/2

π/2 1

p0( x)= y 0=0


Ordenando: h = π/4

Tabela

Inicialmente temos

Daí obtemos

x y Δf(x) Δ2f(x)

0 0 √2/2 1-√2

π/4 √2/2 1-√2/2

π/2 1

p0( x)= y 0=0

p1(x )= p0(x )+Δ f (x0)

1! h( x−x0)


Ordenando: h = π/4

Tabela

Inicialmente temos

Daí obtemos

x y Δf(x) Δ2f(x)

0 0 √2/2 1-√2

π/4 √2/2 1-√2/2

π/2 1

p0( x)= y 0=0

p1(x )= p0(x )+Δ f (x0)

1! h( x−x0)=p0(x)+

√2/21×π/4

(x−0)= p0(x )+2π √2 x


Claramente é um polinômio que interpola os dois primeiros pontos

Calculemos em x=π/6

idêntico aos cálculos anteriores, como era de se esperar.

Adicionemos mais um ponto

p1(π/6)=p0(π/6)+2π √2 π

6=0+

√23

=0,471404


h = π/4

Tabela

Já temos

x y Δf(x) Δ2f(x)

0 0 √2/2 1-√2

π/4 √2/2 1-√2/2

π/2 1

p0(x)= y0=0; p1(x)=p0( x)+2π √2 x


h = π/4

Tabela

Já temos

E assim

x y Δf(x) Δ2f(x)

0 0 √2/2 1-√2

π/4 √2/2 1-√2/2

π/2 1

p0(x)= y0=0; p1(x)=p0( x)+2π √2 x

p2(x )=p1(x)+Δ2 f (x0)

2! h2 (x−x0)(x−x1)=p1(x )+1−√2

2×π2/16

(x−0)(x−π/4)

p2(x )=p1(x)+8

π2(1−√2) x (x−π/4)


É este polinômio que interpola os três pontos


Temos o resultado anterior para o polinômio do primeiro grau. Assim...

p2(π/6)=p1(π/6)+8

π2(1−√2) π

6( π

6−π

4)




Temos o resultado anterior para o polinômio do primeiro grau. Assim,

p2(π/6)=√23

+8(1−√2)16(

16−

14)=

√23

−1−√2

9=0,517428

p2(π/6)=p1(π/6)+8

π2(1−√2) π

6( π

6−π

4)




Temos o resultado anterior para o polinômio do primeiro grau. Assim,

novamente o resultado é idêntico aos anteriores.

p2(π/6)=√23

+8(1−√2)16(

16−

14)=

√23

−1−√2

9=0,517428

p2(π/6)=p1(π/6)+8

π2(1−√2) π

6( π

6−π

4)


O mesmo ocorrerá se calcularmos o polinômio em x=π/3

Mas e se queremos maior precisão?


Usando das simetrias do seno, podemos adicionar um ponto

x y Δf(x) Δ2f(x) Δ3f(x)

0 0 √2/2 1-√2

π/4 √2/2 1-√2/2

π/2 1

3π/4 √2/2



e construirmos uma tabela ampliada! Vejamos o polinômio de terceiro grau


0 0 √2/2 1-√2 2√2-3

π/4 √2/2 1-√2/2 √2-2

π/2 1 √2/2-1

3π/4 √2/2



e construirmos uma tabela ampliada! Vejamos o polinômio de terceiro grau


0 0 √2/2 1-√2 2√2-3

π/4 √2/2 1-√2/2 √2-2

π/2 1 √2/2-1

3π/4 √2/2

p3(x)=p2( x)+Δ

3 f (x0)

3 !h3(x−x0)(x−x1)(x−x2)=p2(x )+

2√2−3

6×π3/64(x−0)(x−π/4)(x−π/2)


Ou dando uma arrumada

Calculando em x=π/6

p3(x)= p2( x)+32

3π3(2√2−3) x (x−π/4)(x−π/2)

p3(π/6)=p2(π/6)+32

3 π3(2√2−3) π

6( π

6−π

4)( π

6− π

2)


Ou dando uma arrumada

Calculando em x=π/6

p3(x)= p2( x)+32

3π3(2√2−3) x (x−π/4)(x−π/2)

p3(π/6)=p2(π/6)+32

3 π3(2√2−3) π

6( π

6−π

4)( π

6− π

2)

p3(π/6)=√23

−1−√2

9+

323

(2√2−3)16(16−

14)(

16−

12)=√2

3−

1−√29

+481

(2√2−3)

p3(π/6)=0,508955


As vantagens desta fórmula se tornam evidentes:

● Adicionar um ponto novo é de custo baixo (O(n))

● Podemos aproveitar todos os cálculos anteriores

● Não necessitamos resolver sistemas


E o custo total disto?

Construir a tabela

Fazemos n-1 subtrações na primeira coluna, n-2 na segunda, n-3 na terceira, etc. Ou seja, (n-1)(n)/2 subtrações que nos dá custo O(n2)



Cálculo de um valor do polinômio

Observe que essencialmente o número de operações para o cálculo dos pontos interpolados é quadrático

p0(x)= y0 ; p1(x)=p0(x)+Δ f (x0)

h(x−x0); p2(x)=p1(x)+

Δ2 f (x0)

2h2(x−x0)( x−x1)

p3(x)=p2( x)+Δ

3 f (x0)

6h3(x−x0)(x−x1)(x−x2)



● Construir a tabela custa O(n2)

● Cada cálculo do polinômio custa O(n2)

● Só é válida para pontos ordenados e igualmente espaçados

Documents

Introdução aos Métodos Numéricos...Fórmula de Lagrange Observação importante: deixe a expressão deste polinômio neste formato. É trabalho inútil e sujeito a erros tentar