19. Polinômio Interpolador de Newton (Diferenças Divididas)#

19.1. Interpolação Linear#

Excluindo-se o caso de função constante, a forma mais simples de interpolação é ligar dois pontos dados com uma reta. Usando semelhança de triângulos entre nós e valores de função, obtemos

\frac{f_{1} (x) - f (x_{0})}{x - x_{0}} = \frac{f (x_{1}) - f (x_{0})}{x_{1} - x_{0}}

a qual pode ser reorganizada para fornecer

f_{1} (x) = f (x_{0}) + \frac{f (x_{1}) - f (x_{0})}{x_{1} - x_{0}} (x - x_{0}) (1)

A notação $f_{1} (x)$ indica que esse é um polinômio interpolador de primeiro grau. Observe que, além de representar a inclinação da reta ligando os pontos, o termo $[f (x_{1}) - f (x_{0})] / (x_{1} - x_{0})$ é uma aproximação por diferenças divididas da primeira derivada.

19.2. Interpolação Quadrática#

Com três pontos, a interpolação quadrática é obtida a partir de

f_{2} (x) = b_{0} + b_{1} (x - x_{0}) + b_{2} (x - x_{0}) (x - x_{1}) (2)

Um procedimento simples pode ser usado para determinar os valores dos coeficientes.

Para $b_{0}$ , a Equação (2) com $x = x_{0}$ pode ser usada para calcular

b_{0} = f (x_{0}) (3)

A Equação (3) pode ser substituída na Equação (2), a qual pode ser calculada em $x = x_{1}$ para

b_{1} = \frac{f (x_{1}) - f (x_{0})}{x_{1} - x_{0}} (4)

Finalmente, as Equações (3) e (4) podem ser substituídas na Equação (2), a qual pode ser calculada em $x = x_{2}$ e resolvida (depois de algumas manipulações algébricas) por

b_{2} = \frac{\frac{f (x_{2}) - f (x_{1})}{x_{2} - x_{1}} - \frac{f (x_{1}) - f (x_{0})}{x_{1} - x_{0}}}{x_{2} - x_{0}} (5)

19.3. Forma Geral dos Polinômios Interpoladores de Newton#

A análise anterior pode ser generalizada para ajustar um polinômio de grau $n$ a $n + 1$ pontos dados. O polinômio de grau $n$ é

f_{n} (x) = b_{0} + b_{1} (x - x_{0}) + \dots + b_{n} (x - x_{0}) (x - x_{1}) \dots (x - x_{n - 1}) (6)

Como foi feito anteriormente com as interpolações linear e quadrática, os pontos dados podem ser usados para calcular os coeficientes $b_{0}$ , $b_{1}$ , \dots , $b_{n}$ . Para um polinômio de grau $n$ , $n + 1$ pontos dados são necessários: $(x_{0}, f (x_{0})), (x_{1}, f (x_{1})), \dots, (x_{n}, f (x_{n}))$ . Usamos esses pontos dados e as seguintes equações para calcular os coeficientes

\begin{array}{r} b_{0} = f (x_{0}) (7) \\ b_{1} = f [x_{1}, x_{0}] (8) \\ b_{2} = f [x_{2}, x_{1}, x_{0}] (9) \\ ⋮ \\ b_{n} = f [x_{n}, x_{n - 1}, \dots, x_{1}, x_{0}] (10) \end{array}

onde a função com colchetes corresponde a diferenças divididas. Por exemplo, a primeira diferença dividida finita é representada em geral por

f [x_{i}, x_{j}] = \frac{f (x_{i}) - f (x_{j})}{x_{i} - x_{j}}

A segunda diferença dividida finita, que representa a diferença das duas primeiras diferenças divididas, é expressa em geral por

f [x_{i}, x_{j}, x_{k}] = \frac{f [x_{i}, x_{j}] - f [x_{j}, x_{k}]}{x_{i} - x_{k}}

Analogamente, a n-ésima diferença dividida é

f [x_{n}, x_{n - 1}, \dots, x_{1}, x_{0}] = \frac{f [x_{n}, x_{n - 1}, \dots, x_{1}] - f [x_{n - 1}, x_{n - 2}, \dots, x_{0}]}{x_{n} - x_{0}}

19.4. Funções de Base de Newton#

Código gerador de funções de base de Newton de grau $n$ por computação simbólica.

from sympy import Symbol

def symbolic_vector(n,var):
    """Cria uma lista com n variáveis simbólicas.
    
        entrada:
            n: numero de pontos
            var: uma string (ex. 'x')
        saida:
            V: ['var0','var1',...,'varn-1']
    """
    if not isinstance(var,str):
        raise TypeError("{0} must be a string.".format(var))                 
        
    V = [Symbol(var  + str(i)) for i in range(0,n)]
        
    return V

def N_nj(X,j):
    """ Calcula a função de base de Newton N_{n,j}(x). 
        
        entrada:
            X: um vetor contendo variáveis simbólicas
    """
    # pega a variavel base passada e converte para simbólica
    x = X[1]
    x = str(x)
    x = Symbol(x[0:-1])
    N = x/x;
    if j > 0: 
        for k in range(0,j):
            N *= (x - X[k])
            
    return N

import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline

# número de nós de interpolação: interpolação de (n-1)-ésimo grau
n = 8

# domínio de interpolação
x0,x1 = -1,1

# constroi vetor simbolico em x
X = symbolic_vector(n,'x')

# constroi pontos numericos 
xp = np.linspace(x0,x1,num=n,endpoint=True)

# cria malha numérica
xv = np.linspace(x0,x1)

# matriz das funções
Y = np.zeros((n,len(xv)))
for i in range(0,n):
    Y[i,] = np.zeros(np.shape(xv))

# montagem de dict para substituição: [xk,x0,x1,x2,...]
k = [str(i) for i in X]
k.insert(0,'x')

# preenche matriz
for i in range(0,Y.shape[0]):
    for j in range(0,np.size(xv)):
        v = list(np.concatenate([np.asarray([xv[j]]),xp]))
        d = dict(zip(k,v))
        Y[i,j] = N_nj(X,i).subs(d)

# plotagem das funções

# nós 
plt.scatter(xp,np.zeros(xp.shape),c= 'k')

leg = []
for i in range(0,Y.shape[0]):
    plt.plot(xv,Y[i,])
    s = '$N_{' + str(i) + '}(x)$'
    leg.append(s)
        
plt.grid()        
plt.legend(leg,loc='best')
plt.title('Funções de base de Newton até ordem ' + str(n-1) + ' em ['+str(x0)+','+str(x1)+']');

_images/aula-15-interpolacao-newton_3_0.png

19.4.1. Exemplo:#

Encontre o polinômio interpolador de Newton de ordem 1 $P_{1} (x)$ para a tabela abaixo

x	y
-1	4
0	1

Compute o valor de $P_{1} (0.35)$ .

import numpy as np
import matplotlib.pyplot as plt

# interpolação linear

# coeficiente a0 = y0
# coeficiente a1 = (y1-y0)/(x1-x0)

# pontos
x0,y0 = -1,4
x1,y1 =  0,1

# ordem 0
a0 = y0

# interpolador de Newton
a1 = (y1-y0)/(x1-x0)
P1 = lambda x: a0 + a1*(x-x0)

# ponto interpolado
xp = -0.35
yp = P1(xp)
yp

# plotagem

# nós
plt.plot([x0,x1],[0,0],'ok')

# valores nodais
plt.plot([x0,x1],[y0,y1],'ok')

# interpolador
x = np.linspace(x0,x1,30,endpoint=True)
plt.plot(x,P1(x),label='$P_1(x)$')

# ponto interpolado
plt.plot(xp,0,'sr')
plt.plot(xp,yp,'sr')
plt.axvline(xp,0,yp,c='r',ls='dashed')

plt.grid()
plt.legend(loc='best');

_images/aula-15-interpolacao-newton_5_0.png

P1(xp)

2.05

19.4.2. Exemplo:#

Encontre o polinômio interpolador de Newton de ordem 2 $P_{2} (x)$ para a tabela abaixo

x	y
-1	4
0	1
2	-1

Compute o valor de $P_{2} (0.35)$ .

# leitura de planilha / dataset
x,y = np.loadtxt('t.txt',unpack=True)

x = np.concatenate( (x, [0.5]) )
y = np.concatenate( (y, [-4.2]) )

for i in range(len(x)):
    exec(f'x{i} = {x[i]}')
    exec(f'y{i} = {y[i]}')


DD0 = y0 # f[x0]
DD1 = (y1 - y0)/(x1 - x0) # f[x0,x1]
DD2 =  ( (y2 - y1)/(x2 - x1) - (y1 - y0)/(x1 - x0) ) / ( x2 - x0 )
DD3 =  ( ( (y3 - y2)/(x3 - x2) - (y2 - y1)/(x2 - x1) ) / ( x3 - x1 ) - DD2 ) / (x3 - x0)

# interpolante
PIN_2 = lambda x: DD0 + DD1*(x - x0) + DD2*(x - x0)*(x - x1)
PIN_3 = lambda x: PIN_2(x) + DD3*(x - x0)*(x - x1)*(x - x2)

plt.plot(x,y*0,'vr',ms=10, alpha=1)
xx = np.linspace(x.min(),x.max(),180)
plt.plot(xx,PIN_3(xx),'sg',ms = 20, alpha=0.2)
plt.plot(x,y,'oy',ms=10, alpha=1)

[<matplotlib.lines.Line2D at 0x7fdf49e808b0>]

_images/aula-15-interpolacao-newton_8_1.png

# interpolação quadrática

# Usando tabela DD: 
# https://vnicius.github.io/numbiosis/interpolador-newton/index.html


# par adicional 
x2,y2 = 2.,-1.

# coeficiente a2 = f[x0,x1,x2] = ( f[x1,x2] - f[x0,x1] ) / (x2 - x0)
# a2 = ( (y2-y1)/(x2-x1) - (y1-y0)/(x1-x0) )/(x2-x0)*(xx-x0)*(xx-x1)

# interpolador de Newton
P2 = lambda xx: P1(xx) + ( ( (y2-y1)/(x2-x1) - (y1-y0)/(x1-x0) )/(x2-x0) )*(xx-x0)*(xx-x1)


# ponto interpolado
yp = P2(xp)
yp

# plotagem

# nós
plt.plot([x0,x1,x2],[0,0,0],'ok')

# valores nodais
plt.plot([x0,x1,x2],[y0,y1,y2],'ok')

# interpolador
x = np.linspace(x0,x2,30,endpoint=True)
plt.plot(x,P2(x),label='$P_2(x)$')

# ponto interpolado
plt.plot(xp,0,'sr')
plt.plot(xp,yp,'sr')
plt.axvline(xp,0,yp,c='r',ls='dashed')

plt.grid()
plt.legend(loc='best');

_images/aula-15-interpolacao-newton_9_0.png

P2(xp)

1.8983333333333332

19.4.3. Exemplo:#

Encontre o polinômio interpolador de Newton de ordem 3 $P_{3} (x)$ para a tabela abaixo

x	y
-1	4
0	1
2	-1
3	1

Compute o valor de $P_{3} (0.35)$ .

# interpolação quadrática

# Usando tabela DD: 
# https://vnicius.github.io/numbiosis/interpolador-newton/index.html


# par adicional 
x3,y3 = 3.,1.

# coeficiente a3 = f[x0,x1,x2,x3] 

# interpolador de Newton
P3 = lambda xxx: P2(xxx) + 1/12*(xxx-x0)*(xxx-x1)*(xxx-x2)


# ponto interpolado
yp = P3(xp)
yp

# plotagem

# nós
plt.plot([x0,x1,x2,x3],[0,0,0,0],'ok')

# valores nodais
plt.plot([x0,x1,x2,x3],[y0,y1,y2,y3],'ok')

# interpolador
x = np.linspace(x0,x3,30,endpoint=True)
plt.plot(x,P3(x),label='$P_3(x)$')

# ponto interpolado
plt.plot(xp,0,'sr')
plt.plot(xp,yp,'sr')
plt.axvline(xp,0,yp,c='r',ls='dashed')

plt.grid()
plt.legend(loc='best');

_images/aula-15-interpolacao-newton_12_0.png

P3(-0.35)

1.9428854166666665

# ponto interpolado
YP = [P1(xp),P2(xp),P3(xp)]

# plotagem
# nós
plt.plot([x0,x1,x2,x3],[0,0,0,0],'ok')

# valores nodais
plt.plot([x0,x1,x2,x3],[y0,y1,y2,y3],'ok')

# interpoladores
x = np.linspace(x0,x3,30,endpoint=True)
plt.plot(x,P1(x),'b',label='$P_1(x)$')
plt.plot(x,P2(x),'c',label='$P_2(x)$')
plt.plot(x,P3(x),'g',label='$P_3(x)$')

# ponto interpolado
plt.plot(xp,0,'sr')
plt.plot([xp,xp,xp],YP,'sr')
plt.axvline(xp,0,max(YP),c='r',ls='dashed')

plt.grid()
plt.legend(loc='best');

_images/aula-15-interpolacao-newton_14_0.png

19.5. Comparação (zoom)#

# interpoladores
x = np.linspace(-0.4,-0.3,30,endpoint=True)
plt.plot(x,P1(x),'b',label='$P_1(x)$')
plt.plot(x,P2(x),'c',label='$P_2(x)$')
plt.plot(x,P3(x),'g',label='$P_3(x)$')

# ponto interpolado
plt.plot([xp,xp,xp],YP,'sr')
plt.axvline(xp,0,max(YP),c='r',ls='dashed')

plt.grid()
plt.legend(loc='best');

_images/aula-15-interpolacao-newton_16_0.png

import sympy as sym
from sympy.abc import x

P3x = 512 + 63*(x-2015) - 2.34*(x-2015)*(x-2017) - 2.45*(x-2015)*(x-2017)*(x-2020)
P3xn = sym.lambdify(x, P3x) 

20092402659.2048

x = np.array([2015,2017,2020,2022])
y = np.array([512,638,792,700])

xn = np.linspace(2015,2022,num=100)

plt.plot(x,y,'o')
plt.plot(xn,P3xn(xn))

xx = 2021
yy = P3xn(xx)
plt.plot(xx,yy,'sg')
plt.axvline(x=xx)
print(f'Censo em {xx} = {yy}')

Censo em 2021 = 775.0400000000081

_images/aula-15-interpolacao-newton_18_1.png

Métodos Numéricos para Ciências Computacionais e Engenharia

Polinômio Interpolador de Newton (Diferenças Divididas)

Contents