Polynomial Interpolation Interpolate(f,a,b,n,\'spacing) that returns an n-degree
ID: 2261198 • Letter: P
Question
Polynomial Interpolation Interpolate(f,a,b,n,'spacing) that returns an n-degree polynomial p that approximates the function f on the interval [a,b] polyfit. If the input argument spacing is omitted, then space the (n + 1) inter- polation nodes evenly (the function linspace may be useful). If spacing is equal to Your task is to write a Matlab function function p using chebyshev', then use the Chebyshev nodes given by 2k-1 Test your function on the polynomial f(x) = 25,7 on the interval [-1,1] with n- 1,2,...,15. Compare the valucs of your interpolating functions pa(x) and () against f(x) for x=-1 : 0.001 : 1, and print out a table with the information 10 1. The output polynomial p should be a function handle, not a vector representing a 2. When using Chebyshev nodes, try to define the nodes z, entirely using vector notation. 3. As a bonus, consider looking up the extra features of polyfit to make your function polynomial. Consider combining anonymous functions with polyval to define p. You should be able to do it without using a for loop. work well in cases where [a, b] is badly scaled or not centered around zeroExplanation / Answer
1.interpolation
% Define f(x)
f=@(x) 1./(1+25.*x.^2);
% Define a grid to evaluate f(x) and P(x) on
m=10000;
x_int=linspace(-1,1,m);
x_int=x_int';
fg=f(x_int);
pause('on'); % Switch pausing on
jp=1; % Stepsize for n
nmax=15; % Max n
figure(1);
for n=jp:jp:nmax
x= cos( (2* (0:n-1)+1) / (2*n) * pi );
% Define the set of equally-spaced interpolation points x
y=f(x); % Compute the data
x=x';
y=y';
% Compute polynomial interpolant P(x)
w=baryweights(x);
P=baryinterp(x,w,y,x_int);
errP=max(abs(fg-P));
% Plot f(x) and P(x)
plot(x_int,fg,x_int,P);
pause; % Pause the for loop until keystroke
end
Q=polyfit(x_int,fg,x_int,10)
end
baryinterp.m
function u= baryinterp(x,w,y,grid)
n=length(x);
m=length(grid);
u=zeros(m,1);
for i=1:m
diff=grid(i).*ones(n,1)-x;
l=sum(diff==0);
%disp([num2str(diff)]);
%disp(['l=' num2str(l)]);
if l==0
z=w./diff;
u(i)=(z'*y)/sum(z);
else
u(i)=y(diff==0);
end
%disp(['u=' num2str(u(i))]);
end
end
baryweights.m
function w = baryweights(x)
%
% Description
n=length(x);
w=zeros(n,1);
if i ==0
w = 1/2;
elseif i == 1:n-1
w = (-1)^i;
else
w = (1/2)*((-1)^n);
end
end