El uso de polinomios de interpolación de alto grado puede producir
errores grandes debido al alto grado de oscilación de este tipo de
polinomios. Para evitar este problema se aproxima la función desconocida
en intervalos pequeños usando polinomios de grado bajo. El caso más
común de la interpolación por pedazos es usar polinomios
cúbicos.
Polinomios Cúbicos de Lagrange: Suponemos que n es divisible por tres y buscamos un polinomio cúbico que interpole a la función en los puntos . Este polinomio esta dado en su representación de Lagrange por:
para y
se asume constante
para toda i.
Interpolación de Hermite: Aqui buscamos un polinomio por pedazos Hn(x) que sea cúbico en cada subintervalo , y que interpole a f(x) y f'(x) en los puntos . La función Hn(x) queda determinada en forma única por estas condiciones y su cálculo requiere de la solución de n sistemas lineales de tamaño 4x4 cada uno. La desventaja de la interpolación de Hermite es que requiere de la disponibilidad de los lo cual no es el caso en muchas aplicaciones.
Interpolación usando Splines: Los dos tipos de polinomios por pedazos que hemos discutidos hasta ahora tienen la desventaja de que su segunda derivada no es continua en los puntos de interpolación. Se ha observado que en aplicaciones gráficas, el ojo humano es capaz de detectar discontinuidades en la segundas derivadas de una función, haciendo que los gráficos con este tipo de funciones no luscan uniformes. Esto motiva el uso de los splines que son funciones s(x) continuas por pedazos con las siguientes propiedades:
Si escribimos , entonces tenemos un total de 4n desconocidas. Las condiciones 2) y 4) nos dan 3(n-1) ecuaciones mientras que de 3) obtenemos n+1 para un total de 4n-3(n-1)-(n+1)=2 grados de libertad. Estos grados de libertad se fijan imponiendo condiciones de frontera adicionales en s(x).
Defina . Como s(x) es cúbico en , entonces s"(x) es lineal en y esta dado por:
Integrando esta ecuación dos veces, obtenemos que
Las condiciones
implican que
Sustituyendo esto en la formula anterior de s(x) obtenemos que
Note que s(x) y s"(x) son continuas por construcción y además tenemos que la condición de interpolación se cumple. Falta la condición de que s'(x) sea continua, i.e.,
En , tenemos que
Usando una expresión equivalente en , se obtiene que la condición de continuidad expresada por los limites de derecha e izquierda arriba, implica que
donde
Tenemos aqui n-1 ecuaciones para las n+1 desconocidas .
Condición Natural de los Splines: Para fijar los dos grados de libertad en el spline, requerimos que s(x) sea lineal en los intervalos lo cual es equivalente a las condiciones . Tenemos ahora el siguiente sistema de ecuaciones Am=d para las restantes desconocidas donde:
Note que la matriz de coeficientes de este sistema es tridiagonal y simétrica lo que hace que el spline s(x) pueda ser calculado en forma eficiente. El siguiente programa en MATLAB ensambla la matriz y lado derecho según definidos arriba y resuelve el sistema para determinar los M's:
%
% Los datos estan dados por los vectores x=[x(1)
x(n)] , y=[y(1)
y(n)]
%
n=length(x)
dx=x(2:n)-x(1:n-1);
yp=(y(2:n)-y(1:n-1))./dx;
a=sparse([1:n-2],[1:n-2],(dx(1:n-2)+dx(2:n-1))/3,n-2,n-2);
udiag=sparse([1:n-3],[2,n-2],dx(2:n-2)/6,n-2,n-2);
a=udiag'+a+udiag;
d=yp(2:n-1)-yp(1:n-2);
m=a\d;
La función spline de MATLAB se utiliza para calcular el spline s(x) directamente. En el siguiente ejemplo los datos se obtienen dividiendo el intervalo [-5,5] en seis subintervalos y evaluamos la función atan (tangente inversa) en los puntos de la partición. Luego construimos el spline que interpola en estos puntos y lo graficamos:
%
% Divide el intervalo [-5,5] en cinco pedazos generando
asi seis puntos
%
x=linspace(-5,5,6);
%
% Evalua la función atan en los puntos de la partición
%
y=atan(x);
%
% Calcula la representación del spline que interpola a los datos
%
pp=spline(x,y);
%
% Calcula 100 puntos en el intervalo [-5,5] para las graficas
%
z=linspace(-5,5,100);
%
% Evalua el spline y la función atan en los 100 puntos
%
sval=ppval(pp,z);
y1=atan(z);
%
% Grafica el spline, atan, y los puntos de interpolación en un
mismo
% sistema de coordenadas
%
plot(z,sval,z,y1,x,y,'+')
xlabel('x');
ylabel('y');
title('atan(x) en violeta y s(x) en amarillo')