A interpolación con Matlab/Octave

O comando polyfit

Os coeficientes do polinomio de interpolación pódense calcular co comando polyfit. Por exemplo, se tomamos os nodos 0, 1 e 2, os coeficientes do polinomio que interpola a función exponencial calcúlanse con:
>> c = polyfit(0:2, exp(0:2), 2)

O resultado é un vector:
c =
1.47625   0.24204   1.00000

o que significa que o polinomio de interpolación é 1.47625 x²+0.24204 x+1. Obsérvese que o terceiro argumento de polyfit é o grado do polinomio que se desexa obter: o polinomio interpola á función cando o seu grao é unha unidade menor que o número de nodos. Se solicitamos calcular un polinomio de grao menor ao de interpolación, o comando polyfit realiza un axuste, minimizando a distancia por minimos cadrados entre os valores da función e o polinomio. Por exemplo, no caso anterior, a recta que axusta á exponencial nos tres nodos dados é 3.19453 x+0.50792, que se obtén con:
>> c = polyfit(0:2, exp(0:2), 1)

O comando polyval

Se queremos saber canto vale o polinomio de interpolación nun punto, debemos usar o comando polyval:
>> c = polyfit(0:2, exp(0:2), 2);
>> x = 0:0.01:2; y = polyval(c,x) %valores en x do polinomio de coeficientes c

A gráfica anterior obtívose cos comandos:
>> x = 0:0.01:2; plot(x,exp(x), x,polyval(polyfit(0:2,exp(0:2),2),x),
x,polyval(polyfit(0:2,exp(0:2),1),x), 0:2,exp(0:2),'ok')
>> legend('exponencial','pol. interp.','recta axuste','nodos')
>> set(gca(), 'keypos',2)
>> print -dpng polyfit.png

O comando interp1

Para a interpolación “a cachos” úsase o comando interp1:
>> y1 = interp1(0:2, exp(0:2), 0:0.01:2, 'linear')
>> y2 = interp1(0:2, exp(0:2), 0:0.01:2, 'spline')

O primeiro argumento son os nodos de interpolación; o segundo, os valores da función nos nodos; o terceiro, as abscisas onde avaliar o polinomio; o cuarto, o tipo de interpolación (‘linear’ no caso de interpolación por rectas e ‘spline’, no caso de splines cúbicos). O comando interp1 calcula o valor do poliniomio directamente, non é necesario o paso previo do cálculo dos seus coeficientes.

A gráfica anterior obtívose cos comandos:
>> x = 0:0.01:2; plot(x,exp(x), x,interp1(0:2,exp(0:2),x,'linear'),
x,interp1(0:2, exp(0:2),x,'spline'), 0:2,exp(0:2),'ok')
>> legend('exponencial','linear','spline','nodos')
>> set(gca(), 'keypos',2)
>> print -dpng interp1.png

Exercicio: O administrador dun cluster descobre que un programa que procesa ficheiros tarda 4 s para ficheiros de 1 MB, 6 s para ficheiros de 2 MB e 10 s para ficheiros de 3 MB. Para estimar canto tempo tardará en procesar ficheiros de 1,5 MB, calcula o polinomio de interpolación de grao 2 que pasa polos tres puntos anteriores. Usa tamén o polinomio linear a cachos e o spline cúbico para facer esa estimación.

Advertisements
Esta entrada foi publicada en Cálculo numérico, Interpolación polinomial, Matlab/Octave. Ligazón permanente.

Deixar unha resposta

introduce os teu datos ou preme nunha das iconas:

Logotipo de WordPress.com

Estás a comentar desde a túa conta de WordPress.com. Sair / Cambiar )

Twitter picture

Estás a comentar desde a túa conta de Twitter. Sair / Cambiar )

Facebook photo

Estás a comentar desde a túa conta de Facebook. Sair / Cambiar )

Google+ photo

Estás a comentar desde a túa conta de Google+. Sair / Cambiar )

Conectando a %s