Interpolação polinomial com Python

Um pouco de teoria

Prof. Doherty Andrade — www.metodosnumericos.com.br

Sabemos da teoria de interpolação polinomial que o polinômio de grau mámximo $n$ que interpola os pontos $(x_k,y_k), k=0,1,…n$ deve satisfazer ao sistema $Ax=b$, \begin{equation} \underbrace{\left[ \begin{array}{ccccc} 1 & x_0 & x_0^2 & \ldots & x_0^n \\ 1 & x_1 & x_1^2 & \ldots & x_1^n \\ \vdots & \vdots & \vdots & \cdots & \vdots \\ 1 & x_n & x_n^2 & \ldots & x_n^n \\ \end{array} \right]}_{=A}\cdot \left[ \begin{array}{c} a_0 \\ a_1\\ \vdots \\ a_n \\ \end{array} \right]= \left[ \begin{array}{c} f(x_0) \\ f(x_1)\\ \vdots \\ f(x_n) \\ \end{array} \right]. \end{equation}

Notemos que $\det A \neq 0$, pois os pontos $x_0,x_1,\ldots, x_n$ TEM que serem distintos e $A$ é uma matriz de Vandermonde. Assim, o sistema linear admite uma única solução $a_0, a_1,\ldots, a_n$. Isto é, existe um único polinômio $p_n(x)= a_nx^n+ a_{n-1}x^{n-1}+ \cdots + a_1x+a_0$ de grau máximo $n$ que interpola os pontos.

Note que a solução é $$a= A^{-1}b.$$ In [1]:

Exemplos

Vamos ativar alguns pacotes.

from __future__ import division
import numpy as np
import matplotlib.pyplot as plt
import scipy as sci
from scipy import optimize

Exemplo 1:

Determinar o polinomio de grau 3 que interpola os pontos (0,1), (1,6), (2,5), (3,-8).

xi = np.array([0,1,2,3], dtype='double')
yi = np.array([1,6,5,-8], dtype='double')
A = np.array([xi**3,xi**2,xi**1,xi**0]).transpose()
a = np.linalg.inv(A).dot(yi) # resolvendo o sistema

xx = np.linspace(-0.5,3.25);
plt.plot(xi,yi,'ro',xx,np.polyval(a,xx),'b-')
plt.grid();plt.show(),
#coeficientes do polinomio
a

Out[2]:

array([-1.00000000e+00, -1.77635684e-15,  6.00000000e+00,  1.00000000e+00])

Exemplo 2:

Determinar o polinomio de grau 4 que interpola os pontos (-1,2), (0,3), (1,5), (4,-2) , (5,8).

xi = np.array([-1,0,1,4,5], dtype='double')
yi = np.array([2,3,5,-2,8], dtype='double')
A = np.array([xi**4, xi**3,xi**2,xi**1,xi**0]).transpose()
a = np.linalg.inv(A).dot(yi) # resolvendo o sistema

xx = np.linspace(-1.5,6);
plt.plot(xi,yi,'ro',xx,np.polyval(a,xx),'b-')
plt.grid();plt.show(),
#coeficientes do polinomio
print('os coeficientes do polinônio são:', a)
os coeficientes do polinônio são: [ 0.19166667 -1.08333333  0.30833333  2.58333333  3.        ]

In [ ]:

Tags :

Compartilhe:

Deixe um comentário

O seu endereço de e-mail não será publicado. Campos obrigatórios são marcados com *