Borrar filtros
Borrar filtros

Compute linear interpolant in pp form and find derivative

1 visualización (últimos 30 días)
I would like to compute a linear interpolant in pp-form and then use the function fnder to compute its derivative. [The reason to do this is that I am solving a functional equation using collocation and I approximate the policy function using splines, but this is not so relevant here].
I used to do the following with interp1 and 'pp' (this is a minimum working example), but Matlab says that it will be removed in a later release. So my question is how to update this code.
Note: if I use cubic spline or other splines like pchip or makima I know how to do it, e.g.
pp = spline(x_grid,f_grid);
der_pp = fnder(pp,1);
y = ppval(der_pp,x)
but I want to do it for linear splines
Any help is greatly appreciated!
clear;clc;close all
n = 10;
x_grid = linspace(-5,5,n)';
x_grid_fine = linspace(-5,5,10*n)';
f_grid = x_grid.^2;
% Linear interpolant
pp = interp1(x_grid,f_grid,'linear','pp'); % MATLAB does not like it :(
% Derivative of linear interpolant
der_pp = fnder(pp,1);
fun_f_x = @(x) ppval(pp,x);
fun_df_dx = @(x) ppval(der_pp,x);
figure
plot(x_grid,fun_f_x(x_grid),'-o')
hold on
plot(x_grid_fine,fun_f_x(x_grid_fine))
legend('Raw data','Fitted')
figure
plot(x_grid,2*x_grid,'-o')
hold on
plot(x_grid_fine,fun_df_dx(x_grid_fine))
legend('Raw data','Fitted')

Respuesta aceptada

John D'Errico
John D'Errico el 30 de Dic. de 2023
Editada: John D'Errico el 30 de Dic. de 2023
n = 10;
x_grid = linspace(-5,5,n)';
f_grid = x_grid.^2; % A simple quadratic polynomial
fn = spapi(2,x_grid,f_grid)
fn = struct with fields:
form: 'B-' knots: [-5 -5 -3.8889 -2.7778 -1.6667 -0.5556 0.5556 1.6667 2.7778 3.8889 5 5] coefs: [25 15.1235 7.7160 2.7778 0.3086 0.3086 2.7778 7.7160 15.1235 25] number: 10 order: 2 dim: 1
Convert to a pp-form. I might be able to do this directly in the call to spapi, but why bother?
pp = fn2fm(fn,'pp')
pp = struct with fields:
form: 'pp' breaks: [-5 -3.8889 -2.7778 -1.6667 -0.5556 0.5556 1.6667 2.7778 3.8889 5] coefs: [9×2 double] pieces: 9 order: 2 dim: 1
fnplt(pp)
So clearly a piecewise linear spline. Now differentiate using fnder.
dfdx = fnder(pp)
dfdx = struct with fields:
form: 'pp' breaks: [-5 -3.8889 -2.7778 -1.6667 -0.5556 0.5556 1.6667 2.7778 3.8889 5] coefs: [9×1 double] pieces: 9 order: 1 dim: 1
fnplt(dfdx)
So easy enough, especially if the interp1 usage is going to turn into a pumpkin one day.
  1 comentario
Alessandro Maria Marco
Alessandro Maria Marco el 30 de Dic. de 2023
Thank you so much for your answer!
I wanted to add a comment that might be useful for a general audience. Since the routine spapi returns the interpolant in B-form, one might wonder why I want instead the pp-form (which, after all, requires calling another routine for conversion, fn2fm). The reason is that the pp-form is much faster to evaluate.
This code
V_interp = spapi(2,k_grid,V0);
V_interp = fn2fm(V_interp,'pp');
% omitted
% - Evaluate linear spline in pp-form
V_interp_val = ppval(V_interp,kp_val);
is about 4 times faster on my laptop than this code that does the same:
V_interp = spapi(2,k_grid,V0);
% omitted
% - Evaluate linear spline in B-form
V_interp_val = fnval(V_interp,kp_val);

Iniciar sesión para comentar.

Más respuestas (0)

Categorías

Más información sobre Spline Postprocessing en Help Center y File Exchange.

Etiquetas

Productos


Versión

R2023b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by