4.2 Interpolación polinómica de Lagrange con Python

[ Interpola Lagrange ] [ Ejercicio ] [ Analítico ] [ Algoritmo ]
..


1. Interpolación polinómica de Lagrange

Referencia: Chapra 18.2 p516, Burden 3.1 p80, Rodríguez 6.2 p195

El polinomio de interpolación de Lagrange reformula el polinomio de interpolación de Newton evitando el cálculo de la tabla de diferencias divididas. El método de Lagrange tolera las diferencias entre distancias x de los puntos de muestra.

El polinomio de Lagrange se construye a partir de las fórmulas:

f_{n} (x) = \sum_{i=0}^{n} L_{i} (x) f(x_{i}) L_{i} (x) = \prod_{j=0, j \neq i}^{n} \frac{x-x_j}{x_i - x_j}

Donde una vez que se han seleccionado los puntos a usar que generan la misma cantidad de términos que puntos.

[ Interpola Lagrange ] [ Ejercicio ] [ Analítico ] [ Algoritmo ]
..


2. Ejercicio

Dados los 4 puntos en la tabla se requiere un polinomio de grado 3.

p(x) = a_0 x^3 + a_1 x^2 + a_2 x^1 + a_3
xi 0 0.2 0.3 0.4
fi 1 1.6 1.7 2.0

[ Interpola Lagrange ] [ Ejercicio ] [ Analítico ] [ Algoritmo ]

..


3. Desarrollo analítico

Se calculan cuatro términos de Lagrange según las fórmulas,

término 1

L_{0} (x) = \frac{(x-0.2)(x-0.3)(x-0.4)}{(0-0.2)(0-0.3)(0-0.4)}

término 2

L_{1} (x) = \frac{(x-0)(x-0.3)(x-0.4)}{(0.2-0)(0.2-0.3)(0.2-0.4)}

término 3

L_{2} (x) = \frac{(x-0)(x-0.2)(x-0.4)}{(0.3-0)(0.3-0.2)(0.3-0.4)}

término 4

L_{3} (x) = \frac{(x-0)(x-0.2)(x-0.3)}{(0.4-0)(0.4-0.2)(0.4-0.3)}

se construye el polinomio usando la fórmula para fn(x) para cada valor fi,

p_3(x) = 1 L_{0} (x) + 1.6 L_{1} (x) + 1.7 L_{2} (x) + 2 L_{3} (x) p_3(x) = 1 \frac{(x-0.2)(x-0.3)(x-0.4)}{(0-0.2)(0-0.3)(0-0.4)} + + 1.6 \frac{(x-0)(x-0.3)(x-0.4)}{(0.2-0)(0.2-0.3)(0.2-0.4)} + 1.7 \frac{(x-0)(x-0.2)(x-0.4)}{(0.3-0)(0.3-0.2)(0.3-0.4)} + 2 \frac{(x-0)(x-0.2)(x-0.3)}{(0.4-0)(0.4-0.2)(0.4-0.3)}

La simplificación de la expresión del polinomio se puede dejar como tarea o realizarla con Sympy con la instrucción sym.expand().

p_3(x) = 41.6666 x^3 - 27.5 x^2 + 6.8333 x + 1

Una forma de verificar que el polinomio es correcto es usar un punto original xi[i] y comprobar que la evaluación tiene como resultado fi[i].

[ Interpola Lagrange ] [ Ejercicio ] [ Analítico ] [ Algoritmo ]
..


3. Algoritmo en Python

En el algoritmo en Python, para construir las expresiones de cada término de Lagrange se usa la forma simbólica con Sympy, con los datos:

xi = [0, 0.2, 0.3, 0.4]
fi = [1, 1.6, 1.7, 2.0]

En cada término Li(x) se usan todos los elementos i , excepto el mismo elemento i, en el numerador y denominador de la expresión.

En polinomio se agrupan todos los términos multiplicados por su respectivo valor fi[i].

    valores de fi:  [1.  1.6 1.7 2. ]
divisores en L(i):  [-0.024  0.004 -0.003  0.008]

Polinomio de Lagrange, expresiones
400.0*x*(x - 0.4)*(x - 0.3) 
  - 566.666666666667*x*(x - 0.4)*(x - 0.2) 
  + 250.0*x*(x - 0.3)*(x - 0.2) 
  - 41.6666666666667*(x - 0.4)*(x - 0.3)*(x - 0.2)

Polinomio de Lagrange: 
41.666666666667*x**3 - 27.5*x**2 + 6.8333333333334*x + 1.0

Las expresiones del polinomio contiene los binomios en su forma básica, para resolver y simplificar las ecuaciones se usa polinomio.expand().

Para realizar la gráfica del polinomio es conveniente convertirlo a la forma lambda con Numpy, de esta forma se evalúa en una sola linea todos los puntos para el intervalo [a,b] en x.

# Interpolacion de Lagrange
# divisoresL solo para mostrar valores
import numpy as np
import sympy as sym
import matplotlib.pyplot as plt

# INGRESO , Datos de prueba
xi = [0, 0.2, 0.3, 0.4]
fi = [1, 1.6, 1.7, 2.0]

# PROCEDIMIENTO
xi = np.array(xi,dtype=float)
fi = np.array(fi,dtype=float)
# Polinomio de Lagrange
n = len(xi)
x = sym.Symbol('x')
polinomio = 0
divisorL = np.zeros(n, dtype = float)
for i in range(0,n,1):
    
    # Termino de Lagrange
    numerador = 1
    denominador = 1
    for j  in range(0,n,1):
        if (j!=i):
            numerador = numerador*(x-xi[j])
            denominador = denominador*(xi[i]-xi[j])
    terminoLi = numerador/denominador

    polinomio = polinomio + terminoLi*fi[i]
    divisorL[i] = denominador

# simplifica el polinomio
polisimple = polinomio.expand()

# para evaluación numérica
px = sym.lambdify(x,polisimple)

# Puntos para la gráfica
muestras = 101
a = np.min(xi)
b = np.max(xi)
pxi = np.linspace(a,b,muestras)
pfi = px(pxi)

# SALIDA
print('    valores de fi: ',fi)
print('divisores en L(i): ',divisorL)
print()
print('Polinomio de Lagrange, expresiones')
print(polinomio)
print()
print('Polinomio de Lagrange: ')
print(polisimple)

# Gráfica
plt.plot(xi,fi,'o', label = 'Puntos')
plt.plot(pxi,pfi, label = 'Polinomio')
plt.legend()
plt.xlabel('xi')
plt.ylabel('fi')
plt.title('Interpolación Lagrange')
plt.show()

[ Interpola Lagrange ] [ Ejercicio ] [ Analítico ] [ Algoritmo ]