Problem in finding the correct Chebyshev interpolation line
2 visualizaciones (últimos 30 días)
Mostrar comentarios más antiguos
Ignatius Tommy Pratama
el 27 de Mzo. de 2024
Comentada: Mathieu NOE
el 2 de Abr. de 2024
Hi, I have problem in finding the correct Chebyshev interpolation line. I have inputted the calculation procedure in accordance with the Chebyshev theorem. However the regression line is apparently does not agree with the trend of the provided data point. Please let me know I missed something or how to improve my code. Your help is highly appreciated.
%This program is for governing a polinomial equation using Chebyshev
%interpolation formula
clc; clear all; close all;
t=tic;
%% Data input
x_data = [1830; 1920; 1940; 1975; 2000; 2020]; %x values
y_data = [0.089; 3.52; 4.86; 17.05; 25.5; 35.01]; %y values
%% Chebyshev interpolating formula
% Normalize x-values
x_min = min(x_data);
x_max = max(x_data);
x_data_norm = 2 * (x_data - x_min) / (x_max - x_min) - 1;
% Number of data points
n = length(x_data);
% Compute Chebyshev nodes
x_nodes = cos((2*(1:n) - 1) * pi / (2*n));
% Compute Chebyshev coefficients
a = zeros(n, 1);
for i = 1:n
T = cos((i-1) * acos(x_nodes));
a(i) = (2/n) * sum(y_data .* T');
end
% Define the Chebyshev polynomial function
T = @(k, x) cos((k-1) * acos(x));
% Evaluate the interpolating polynomial
x_query = linspace(-1, 1, 100); % Query points
y_interp = zeros(size(x_query));
for i = 1:n
y_interp = y_interp + a(i) * T(i, x_query);
end
%% Plotting the data and show polynomial equation
% Plot the given data points and the interpolating polynomial
% Construct the interpolation equation
interpolation_eqn = '';
for i = 1:n
if i == 1
interpolation_eqn = [interpolation_eqn num2str(a(i))];
else
interpolation_eqn = [interpolation_eqn ' + ' num2str(a(i)) ' * T' num2str(i-1) '(x)'];
end
end
disp(['Interpolation Equation:']);
disp(['P(x) = ' interpolation_eqn]);
% Construct the interpolation equation in terms of powers of x
expanded_eqn = '';
for i = 1:n
if i == 1
expanded_eqn = [expanded_eqn num2str(a(i))];
else
expanded_eqn = [expanded_eqn ' + ' num2str(a(i)) ' * (2 * cos(' num2str(i-1) '* acos(x)) - 1)^' num2str(i-1)];
end
end
disp(['Expanded Interpolation Equation:']);
disp(['P(x) = ' expanded_eqn]);
% Plot the data
figure;
plot(x_data, y_data, 'ro', 'MarkerSize', 10); % plot data points
hold on;
plot(x_query*(x_max-x_min)/2+(x_max+x_min)/2, y_interp, 'b-', 'LineWidth', 2); % plot interpolating polynomial
xlabel('Year');
ylabel('CO2 Concentration');
title('Chebyshev Interpolation');
legend('Given Data', 'Interpolating Polynomial', 'Location', 'best');
grid on;
%% Estimate P(x) for specific x values
% Evaluate the interpolating polynomial at specific x-values
x_values = [1915, 2010, 2050];
y_values = zeros(size(x_values));
% Evaluate the interpolating polynomial at the given x-values
for i = 1:length(x_values)
x_norm = 2 * (x_values(i) - x_min) / (x_max - x_min) - 1;
for j = 1:n
y_values(i) = y_values(i) + a(j) * T(j, x_norm);
end
end
disp(['For x = 1915, y = ', num2str(y_values(1))]);
disp(['For x = 2010, y = ', num2str(y_values(2))]);
disp(['For x = 2050, y = ', num2str(y_values(3))]);
toc(t)
0 comentarios
Respuesta aceptada
Mathieu NOE
el 27 de Mzo. de 2024
hello
maybe this ?
%% Data input
x_data = [1830; 1920; 1940; 1975; 2000; 2020]; %x values
y_data = [0.089; 3.52; 4.86; 17.05; 25.5; 35.01]; %y values
% Chebyshev polynomial interpolation
x = linspace(x_data(1),x_data(end),25);
% Barycentric interpolation
m = numel(x_data);
ts = -1+1/m:2/m:1-1/m;
c = (-1).^(0:m-1).*sin(pi*(ts+1)/2);
numer = zeros(size(x));
denom = zeros(size(x));
exact = zeros(size(x));
for j = 1:m
xdiff = x-x_data(j);
temp = c(j)./xdiff;
numer = numer + temp*y_data(j);
denom = denom + temp;
exact(xdiff==0) = j;
end
y = numer./denom;
jj = find(exact);
y(jj) = y_data(exact(jj));
% plot
plot(x_data,y_data,'*b',x,y,'r');
3 comentarios
Más respuestas (0)
Ver también
Categorías
Más información sobre Interpolation en Help Center y File Exchange.
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!