Página web del curso

Interpolación

En ocasiones necesitamos calcular el valor de una función en varios puntos pero

  • Puede que no tengamos el valor de la función para todos los puntos sino solo en algunos porque, por ejemplo, hemos obtenido los valores que tenemos a partir de algún experimento.
  • O puede que la función sea muy costosa de calcular desde el punto de vista computacional.
  • O puede que queramos una función, que siendo aproximada a la nuestra, sea más fácil de derivar o integrar o realizar otro tipo de operación.

Es en estos casos cuando la técnica de interpolación se hace útil. En interpolación utilizamos los datos (valor, derivada) de nuestra función en varios puntos de la función a aproximar y construímos una nueva función. Ejemplos de tipos de funciones que se emplean para construir estas nuevas funciones aproximadas son los polinomios y las funciones trigonométricas.

Vamos a ver:

  • Interpolación polinómica de Lagrange.
    • Utilizando los polinomios fundamentales de Lagrange.
    • Utilizando la forma de Newton con diferencias divididas.
  • Interpolación polinómica a trozos.

Interpolación polinómica de Lagrange

En este tipo de interpolación disponemos como datos de:

  • Nodos de interpolación: son $n+1$ puntos $x_0$, $x_1$, $x_2,\ldots,x_n$
  • Valores de la función $f$ en los nodos: $f(x_0)$, $f(x_1)$, $f(x_2),\ldots,f(x_n)$

Y buscamos un polinomio que pase por todos estos puntos. Esta puede ser toda la información que tengamos de $f$ y no necesariamente existe una expresión analítica de $f.$

Por ejemplo, sean los puntos

\begin{array}{|c|ccccc|} \hline k & 0 & 1 & 2 & 3 & 4 \\ \hline x & 1 & 2 & 3 & 4 & 5 \\ y & 1 & 2 & 4 & 3 & 5 \\ \hline \end{array}

Buscamos un polinomio que pase por todos estos puntos

Existencia y unicidad del polinomio de Lagrange

En el ejemplo del gráfico anterior, como tenemos $5$ puntos podemos plantear $5$ ecuaciones y necesitamos $5$ incógnitas que serán los coeficientes del polinomio. Por lo tanto:

$$P_4(x)=a_0+a_1x+a_2x^2+a_3x^3+a_4x^4$$

es decir $a_0,$ $a_1,$ $a_2,$ $a_3$ y $a_4$ son las $5$ incógnitas y el polinomio es, en principio, de grado $4$. En general, para los puntos $(x_0,y_0),(x_1,y_1),\ldots,(x_n,y_n)$ el polinomio interpolante será de grado $\leq n$.

Las ecuaciones serían:

$$ \begin{array}{lll} P_{4}(x_{0})=y_{0} & \quad & a_{0}+a_{1}x_{0}+a_{2}x_{0}^{2}+a_{3}x_{0}^{3}+a_{4}x_{0}^{4}=y_{0}\\ P_{4}(x_{1})=y_{1} & \quad & a_{0}+a_{1}x_{1}+a_{2}x_{1}^{2}+a_{3}x_{1}^{3}+a_{4}x_{1}^{4}=y_{1}\\ P_{4}(x_{2})=y_{2} & \quad & a_{0}+a_{1}x_{2}+a_{2}x_{2}^{2}+a_{3}x_{2}^{3}+a_{4}x_{2}^{4}=y_{2}\\ P_{4}(x_{3})=y_{3} & \quad & a_{0}+a_{1}x_{3}+a_{2}x_{3}^{2}+a_{3}x_{3}^{3}+a_{4}x_{3}^{4}=y_{3}\\ P_{4}(x_{4})=y_{4} & \quad & a_{0}+a_{1}x_{4}+a_{2}x_{4}^{2}+a_{3}x_{4}^{3}+a_{4}x_{4}^{4}=y_{4} \end{array} $$

Y el sistema a resolver es:

$$ \left(\begin{array}{ccccc} 1 & x_{0} & x_{0}^{2} & x_{0}^{3} & x_{0}^{4}\\ 1 & x_{1} & x_{1}^{2} & x_{1}^{3} & x_{1}^{4}\\ 1 & x_{2} & x_{2}^{2} & x_{2}^{3} & x_{2}^{4}\\ 1 & x_{3} & x_{3}^{2} & x_{3}^{3} & x_{3}^{4}\\ 1 & x_{4} & x_{4}^{2} & x_{4}^{3} & x_{4}^{4} \end{array}\right)\left(\begin{array}{c} a_{0}\\ a_{1}\\ a_{2}\\ a_{3}\\ a_{4} \end{array}\right)=\left(\begin{array}{c} y_{0}\\ y_{1}\\ y_{2}\\ y_{3}\\ y_{4} \end{array}\right) $$

La matriz de coeficientes del sistema se llama Matriz de Vandermonde y su determinante es

$$\det\left(A\right)=\prod_{0\leq l\lt k\leq n}\left(x_{k}-x_{l}\right)\neq0\quad \mathrm{si}\quad x_k \ne x_l$$

es decir

$$\det(A) = (x_1-x_0)(x_2-x_0)(x_3-x_0)(x_4-x_0)(x_2-x_1)(x_3-x_1)(x_4-x_1)(x_3-x_2)(x_4-x_2)(x_4-x_3)$$

que es distinto de cero si no hay dos $x_i$ iguales. Y en este caso, la solución del sistema existe y es única porque es un sistema determinado.

En el ejemplo anterior los puntos $(x_k,y_k),$ con $y_k = f(x_k)$ son

\begin{array}{|c|ccccc|} \hline k & 0 & 1 & 2 & 3 & 4 \\ \hline x & 1 & 2 & 3 & 4 & 5 \\ y & 1 & 2 & 4 & 3 & 5 \\ \hline \end{array}

El sistema a resolver es:

$$ \left(\begin{array}{ccccc} 1 & x_{0} & x_{0}^{2} & x_{0}^{3} & x_{0}^{4}\\ 1 & x_{1} & x_{1}^{2} & x_{1}^{3} & x_{1}^{4}\\ 1 & x_{2} & x_{2}^{2} & x_{2}^{3} & x_{2}^{4}\\ 1 & x_{3} & x_{3}^{2} & x_{3}^{3} & x_{3}^{4}\\ 1 & x_{4} & x_{4}^{2} & x_{4}^{3} & x_{4}^{4} \end{array}\right)\left(\begin{array}{c} a_{0}\\ a_{1}\\ a_{2}\\ a_{3}\\ a_{4} \end{array}\right)=\left(\begin{array}{c} y_{0}\\ y_{1}\\ y_{2}\\ y_{3}\\ y_{4} \end{array}\right) $$$$ \left(\begin{array}{ccccc} 1 & 1 & 1 & 1 & 1\\ 1 & 2 & 4 & 8 & 16\\ 1 & 3 & 9 & 27 & 81\\ 1 & 4 & 16 & 64 & 256\\ 1 & 5 & 25 & 125 & 625 \end{array}\right)\left(\begin{array}{c} a_{0}\\ a_{1}\\ a_{2}\\ a_{3}\\ a_{4} \end{array}\right)=\left(\begin{array}{c} 1\\ 2\\ 4\\ 3\\ 5 \end{array}\right) $$

Como tenemos 5 nodos distintos el determinante no se anula

$$\det\left(A\right)=(2-1)(3-1)(4-1)(5-1)(3-2)(4-2)(5-2)(4-3)(5-3)(5-4)=288$$

Y, resolviendo el sistema, el polinomio interpolador de Lagrange es

$$P_4(x)=15-28.66666667\,x+19.08333333\,x^2-4.83333333\,x^3+0.41666667\,x^4$$

Una forma de obtener el polinomio interpolador de Lagrange es resolver el sistema anterior. Pero este sistema está mal condicionado (pequeñas errores en los datos pueden producir grandes errores en los resultados) y es comparativamente costoso.

Vamos a ver dos formas alternativas de calcular el polinomio interpolador de Lagrange:

  • Mediante los polinomios fundamentales de Lagrange.
  • Mediante las diferencias divididas de Newton.

Ejercicio

Dados los nodos $x_{0}=-1,$ $x_{1}=1,$ $x_{2}=3$ y $x_{3}=5$ y la función $f\left(x\right)=\mathrm{sen}\left(\dfrac{\pi}{6}x\right)$

  • Calcular los polinomios fundamentales de Lagrange y dibujarlos.
  • Calcular el polinomio interpolante por el método de Lagrange.
  • Aproximar el valor en el punto $x=2,$ calcular la cota de error y compararla con el error.

Calcular los polinomios fundamentales de Lagrange y dibujarlos

Los polinomios fundamentales de lagrange tienen un valor 0 en todos los nodos excepto en uno de ellos, donde vale 1. Hay tantos polinomios fundamentales de Lagrange como nodos y los llamaremos $L_0(x),$ $L_1(x),$ $L_2(x)$ y $L_3(x).$

Veamos como se construyen. Primero hacemos que los polinomios valgan cero en todos los nodos de interpolación menos en uno y tendremos

Los nodos son

$$x_{0}=-1\quad x_{1}=1\quad x_{2}=3\quad x_{3}=5$$

Y un polinomio que es cero en los tres últimos los nodos es

$$l_0(x)=(x-x_1)(x-x_2)(x-x_3)=(x-1)(x-3)(x-5)$$

Si dividimos este polinomio por $l_0(x_0)$

$$L_0(x)=\frac{(x-x_1)(x-x_2)(x-x_3)}{(x_0-x_1)(x_0-x_2)(x_0-x_3)} =\frac{(x-1)(x-3)(x-5)}{(-1-1)(-1-3)(-1-5)}$$

Estamos omitiendo en el numerador y el denominador el término en $x_0,$ que es $(x-x_0)$ y

$$L_0(x_0)=1 \quad L_0(x_1)=0 \quad L_0(x_2)=0 \quad L_0(x_3)=0$$

es decir

$$L_k(x_i)=\begin{cases} 0 & \mathrm{si} & k\ne i\\ 1 & \mathrm{si} & k= i \end{cases}$$

Los nodos son

$$x_{0}=-1\quad x_{1}=1\quad x_{2}=3\quad x_{3}=5$$

Y un polinomio que es cero en los tres últimos los nodos es

$$l_1(x)=(x-x_0)(x-x_2)(x-x_3)=(x+1)(x-3)(x-5)$$

Si dividimos este polinomio por $l_1(x_1)$

$$L_1(x)=\frac{(x-x_0)(x-x_2)(x-x_3)}{(x_1-x_0)(x_1-x_2)(x_1-x_3)} =\frac{(x+1)(x-3)(x-5)}{(1+1)(1-3)(1-5)}$$

Estamos omitiendo en el numerador y el denominador el término en $x_1,$ que es $(x-x_1)$ y tenemos que

$$L_1(x_0)=0 \quad L_1(x_1)=1 \quad L_1(x_2)=0 \quad L_1(x_3)=0$$

es decir

$$L_k(x_i)=\begin{cases} 0 & \mathrm{si} & k\ne i\\ 1 & \mathrm{si} & k= i \end{cases}$$

Los nodos son

$$x_{0}=-1\quad x_{1}=1\quad x_{2}=3\quad x_{3}=5$$

Y un polinomio que es cero en todos los nodos menos el tercero es

$$l_2(x)=(x-x_0)(x-x_1)(x-x_3)=(x+1)(x-1)(x-5)$$

Si dividimos este polinomio por $l_2(x_2)$

$$L_2(x)=\frac{(x-x_0)(x-x_1)(x-x_3)}{(x_2-x_0)(x_2-x_1)(x_2-x_3)} =\frac{(x+1)(x-1)(x-5)}{(3+1)(3-1)(3-5)}$$

Estamos omitiendo en el numerador y el denominador el término en $x_2,$ que es $(x-x_2)$ y

$$L_2(x_0)=0 \quad L_2(x_1)=0 \quad L_2(x_2)=1 \quad L_2(x_3)=0$$

es decir

$$L_k(x_i)=\begin{cases} 0 & \mathrm{si} & k\ne i\\ 1 & \mathrm{si} & k= i \end{cases}$$

Los nodos son

$$x_{0}=-1\quad x_{1}=1\quad x_{2}=3\quad x_{3}=5$$

Y un polinomio que es cero en todos los nodos menos el último es

$$l_3(x)=(x-x_0)(x-x_1)(x-x_2)=(x+1)(x-1)(x-3)$$

Si dividimos este polinomio por $l_3(x_3)$

$$L_3(x)=\frac{(x-x_0)(x-x_1)(x-x_2)}{(x_3-x_0)(x_3-x_1)(x_3-x_2)} =\frac{(x+1)(x-1)(x-3)}{(5+1)(5-1)(5-3)}$$

Estamos omitiendo en el numerador y el denominador el término en $x_3,$ que es $(x-x_3)$ y

$$L_3(x_0)=0 \quad L_3(x_1)=0 \quad L_3(x_2)=0 \quad L_3(x_3)=1$$

es decir

$$L_k(x_i)=\begin{cases} 0 & \mathrm{si} & k\ne i\\ 1 & \mathrm{si} & k= i \end{cases}$$

Ya hemos visto los polinomios fundamentales de Lagrange para nuestros nodos (4 nodos, 4 polinomios fundamentales). En general, para $x_0,$ $x_1,\ldots, x_n$ nodos tenemos $n+1$ polinomios fundamentales de Lagrange que son

$$ L_{k}\left(x\right)= \frac{(x-x_{0})}{(x_{k}-x_{0})}\cdots\frac{(x-x_{k-1})}{(x_{k}-x_{k-1})}\frac{(x-x_{k+1})}{(x_{k}-x_{k+1})}\cdots\frac{(x-x_{n})}{(x_{k}-x_{n})} \quad k=0,1,\ldots,n $$

o expresado con productorios

$$ L_{k}\left(x\right)=\prod_{\begin{array}{c} j=0\\ j\neq k \end{array}}^{n}\dfrac{x-x_{j}}{x_{k}-x_{j}} \quad k=0,1,\ldots,n $$

Estos polinomios verifican

$$L_k(x_i)=\begin{cases} 0 & \mathrm{si} & k\ne i\\ 1 & \mathrm{si} & k= i \end{cases}$$

Calcular el polinomio interpolante por el método de Lagrange

El polinomio interpolante es de la forma

$$ P_{3}\left(x\right)=y_{0}L_{0}\left(x\right)+ y_{1}L_{1}\left(x\right)+y_{2}L_{2}\left(x\right)+y_{3}L_{3}\left(x\right) $$

Los polinomios fundamentales de Lagrange se construyen solo con los nodos. Pero ahora, para calcular el polinomio de interpolación $P_3(x)$ necesitamos los valores de la función en los nodos

$$f\left(x\right)=\mathrm{sen}\left(\dfrac{\pi}{6}x\right)\\[0.7cm] y_0=f\left(-1\right)=\mathrm{sen}\left(-\dfrac{\pi}{6}\right)=-\frac{1}{2}\quad y_1=f\left(1\right)=\mathrm{sen}\left(\dfrac{\pi}{6}\right)=\frac{1}{2}\\[0.7cm] y_2=f\left(3\right)=\mathrm{sen}\left(3\dfrac{\pi}{6}\right)=\mathrm{sen}\left(\dfrac{\pi}{2}\right)=1\quad y_3=f\left(5\right)=\mathrm{sen}\left(\dfrac{5\pi}{6}\right)=\frac{1}{2}$$
\begin{array}{|l|cccc|} \hline k & 0 & 1 & 2 & 3\\ \hline x & -1 & 1 & 3 & 5\\ \hline y &-\dfrac{1}{2} & \dfrac{1}{2} & 1 & \dfrac{1}{2}\\ \hline \end{array}
  • Si multiplicamos cada polinomio fundamental de Lagrange $L_i(x)$ por el valor de la función correspondiente $y_i$ el polinomio $y_i\,L_i(x)$ toma el valor $y_i$ en el nodo $i$ y cero en los otros nodos.
  • Así conseguimos 4 polinomios, cada uno de los cuales pasa por uno de los puntos a interpolar.
  • Si sumamos estos 4 polinomios, como la suma se hace punto a punto, el valor en los nodos será el valor de la función, $y_i,$ por lo que este polinomio suma, $P_3(x)$ es el polinomio de interpolación que estábamos buscando.
  • Todos los polinomios son de grado 3, por lo que al sumarlos $P_3(x)$ será de grado menor o igual que 3. En este caso el polinomio es de grado tres pero si, por ejemplo, los 4 puntos estuvieran alineados, el polinomio sería de grado 1. O si los 4 puntos estuvieran sobre una recta horizontal el polinomio sería de grado cero. O si estuvieran sobre una parábola, sería de grado 2.
\begin{array}{|l|cccc|} \hline k & 0 & 1 & 2 & 3\\ \hline x & -1 & 1 & 3 & 5\\ \hline y &-\dfrac{1}{2} & \dfrac{1}{2} & 1 & \dfrac{1}{2}\\ \hline \end{array}

Y teniendo en cuenta los polinomios fundamentales de Lagrange que calculamos arriba, el polinomio de interpolación es

$$ P_{3}\left(x\right)=-\frac{1}{2}\,L_{0}\left(x\right)+ \frac{1}{2}\,L_{1}\left(x\right)+(1)\,L_{2}\left(x\right)+\frac{1}{2}\,L_{3}\left(x\right) $$

El polinomio de interpolación, para un caso general, donde tenemos $n+1$ nodos $x_0,$ $x_1,\ldots, x_n$ será de la forma

$$P_n(x)=y_0\,L_0(x)+y_1\,L_1(x)+\cdots+y_n\,L_n(x)=\sum_{k=0}^n y_k\,L_k(x)$$

donde $P_n(x)$ es un polinomio de grado menor o igual que $n.$

Aproximar el valor en el punto $x=2,$ calcular la cota de error y compararla con el error

$$ P_{3}\left(x\right)=-\frac{1}{2}\,L_{0}\left(x\right)+ \frac{1}{2}\,L_{1}\left(x\right)+(1)\,L_{2}\left(x\right)+\frac{1}{2}\,L_{3}\left(x\right) $$

Sustituyendo $x=2$ en las expresiones que calculamos para los polinomios fundamentales de Lagrange

$$ P_{3}\left(2\right)=-\frac{1}{2}\,L_{0}\left(2\right)+ \frac{1}{2}\,L_{1}\left(2\right)+(1)\,L_{2}\left(2\right)+\frac{1}{2}\,L_{3}\left(2\right) $$
$$ P_{3}\left(2\right)=-\frac{1}{2}\,(-0.0625)+ \frac{1}{2}(0.5625)+(1)\,(0.5625)+\frac{1}{2}\,(-0.0625) $$
$$\mathrm{sen}\left(2\frac{\pi}{6}\right)=\mathrm{sen}\left(\frac{\pi}{3}\right)\approx P_3(2)=0.84375$$

Error de interpolación

El error de interpolación viene dado por

$$ E(x)=f(x)-P_{n}(x)=f^{(n+1)}(c)\dfrac{(x-x_{0})\ldots(x-x_{n})}{(n+1)!}, $$

donde las $x_{i}$ son los puntos de interpolación y $c$ un punto del intervalo de interpolación. En este caso, como tenemos cuatro nodos de interpolación

$$ \left|E(x)\right|=\left|f(x)-P_{3}(x)\right|=\left|f^{(4)}(c)\right|\dfrac{\left|(x-x_{0})(x-x_{1})(x-x_{2})(x-x_{3})\right|}{4!}. $$

Como $c$ es desconocido, aunque sabemos que está en el intervalo de interpolación, en nuestro caso $c\in\left(-1,5\right)$, tenemos que encontrar una cota para ese valor

$$f(x)=\mathrm{sen}\left(\frac{\pi}{6} x\right) \quad f'(x)=\frac{\pi}{6}\,\cos\left(\frac{\pi}{6} x\right) \quad f''(x)=-\frac{\pi^2}{6^2}\,\mathrm{sen}\left(\frac{\pi}{6} x\right)$$
$$ f'''(x)=-\frac{\pi^3}{6^3}\,\cos\left(\frac{\pi}{6} x\right) \quad f^{(4)}(x)=\frac{\pi^4}{6^4}\,\mathrm{sen}\left(\frac{\pi}{6} x\right)$$

Y como $$\left|\mathrm{sen}\left(\dfrac{\pi}{6}x\right)\right|\le 1$$

Entonces $$\left|f^{(4)}(c)\right|\le \frac{\pi^4}{6^4}(1)=\frac{\pi^4}{6^4}$$

Y por lo tanto podemos dar como cota del error

$$ \left|E(2)\right|\le\frac{\pi^4}{6^4}\dfrac{\left|(2-x_{0})(2-x_{1})(2-x_{2})(2-x_{4})\right|}{4!} $$
$$ \left|E(2)\right|\le \frac{\pi^4}{6^4}\dfrac{\left|(2+1)(2-1)(2-3)(2-5)\right|}{4!}=0.0282 $$

Como $f(2)=\mathrm{sen}\left(\dfrac{\pi}{6}(2)\right)=\mathrm{sen}\left(\dfrac{\pi}{3}\right)=0.86603$

$$\mathrm{Error} = |f(2)-P_3(2)|=|0.86603-0.84375|=0.0223\lt 0.0282$$

¿Cómo se evalua el buen o mal condicionamiento de la matriz de coeficientes? Con el número de condición

$$ \mathrm{cond}(A)=\det(A)\det(A^{-1}) $$

El número de condición es siempre mayor que uno, pero cuanto más próximo a uno, mejor condicionada está la matriz y viceversa. En este ejemplo de la matriz de Vandermonde $\mathrm{cond}(A)=26170.$