使用切比雪夫点进行插值

Question

以 10 为增量在 Chebyshev 点对示例 10.6 的 Runge 函数进行插值，n 从 10 到 170。计算均匀评估网格上的最大插值误差 x = -1:.001:1 并绘制误差与多项式次数的关系，如下所示图 10.8 使用符号学。观察光谱精度。

梯级函数由下式给出： f(x) = 1 / (1 + 25x^2)

到目前为止我的代码：

x = -1:0.001:1;
n = 170;
i = 10:10:170;
cx = cos(((2*i + 1)/(2*(n+1)))*pi); %chebyshev pts

y = 1 ./ (1 + 25*x.^2); %true fct
%chebyshev polynomial, don't know how to construct using matlab

yc = polyval(c, x); %graph of approx polynomial fct
plot(x, yc);
mErr =  (1 / ((2.^n).*(n+1)!))*%n+1 derivative of f evaluated at max x in [-1,1], not sure how to do this
%plotting stuff

我对matlab知之甚少，所以我正在努力创建插值多项式。我做了一些谷歌工作，但我对当前的函数感到困惑，因为我没有找到一个只是简单地接受点和多项式进行插值的函数。在这种情况下，我也有点困惑我是否应该做i = 0:1:n，n=10:10:170或者是否n在这里修复。任何帮助表示赞赏，谢谢

Answer 1

由于您对MATLAB知之甚少，我将尝试逐步解释一切：

首先，要可视化 Runge 函数，您可以键入：

f = @(x) 1./(1+25*x.^2); % Runge function

% plot Runge function over [-1,1];
x = -1:1e-3:1;
y = f(x);
figure;
plot(x,y); title('Runge function)'); xlabel('x');ylabel('y');

@(x)代码的一部分是一个函数 handle，这是 MATLAB 一个非常有用的特性。请注意，该函数已正确vecotrized，因此它可以接收变量或数组作为参数。plot 函数很简单。

要理解龙格现象，请考虑一个包含 10 个元素的线性间隔向量，[-1,1]并使用这些点来获得插值（拉格朗日）多项式。您将获得以下信息：

% 10 linearly spaced points
xc = linspace(-1,1,10);
yc = f(xc);
p = polyfit(xc,yc,9); % gives the coefficients of the polynomial of degree 10
hold on; plot(xc,yc,'o',x,polyval(p,x)); 

该polyfit函数做了一个多项式拟合-它获得的插值多项式，给出的poins的系数x,y和多项式的程度n。您可以使用polyval函数轻松计算其他点的多项式。

请注意，靠近端域，您会得到一个振荡多项式，并且插值不是函数的良好近似。事实上，您可以绘制绝对误差，比较函数的值f(x)和插值多项式p(x)：

plot(x,abs(y-polyval(p,x))); xlabel('x');ylabel('|f(x)-p(x)|');title('Error');

如果不使用线性空间向量，而是使用其他点进行插值，则可以减少此误差。一个不错的选择是使用Chebyshev 节点，这应该会减少错误。事实上，请注意：

% find 10 Chebyshev nodes and mark them on the plot
n = 10;
k = 1:10; % iterator
xc = cos((2*k-1)/2/n*pi); % Chebyshev nodes
yc = f(xc); % function evaluated at Chebyshev nodes
hold on;
plot(xc,yc,'o')

% find polynomial to interpolate data using the Chebyshev nodes
p = polyfit(xc,yc,n-1); % gives the coefficients of the polynomial of degree 10
plot(x,polyval(p,x),'--'); % plot polynomial
legend('Runge function','Chebyshev nodes','interpolating polynomial','location','best')

请注意接近端域的误差是如何减少的。你现在没有得到插值多项式的高振荡行为。如果绘制错误，您将观察到：

plot(x,abs(y-polyval(p,x))); xlabel('x');ylabel('|f(x)-p(x)|');title('Error');

现在，如果您更改 Chebyshev 节点的数量，您将获得更好的近似值。对代码稍加修改，您就可以针对不同数量的节点再次运行它。您可以存储最大误差并将其绘制为节点数的函数：

n=1:20; % number of nodes

% pre-allocation for speed
e_ln = zeros(1,length(n)); % error for the linearly spaced interpolation
e_cn = zeros(1,length(n)); % error for the chebyshev nodes interpolation

for ii=1:length(n)
    % linearly spaced vector
    x_ln = linspace(-1,1,n(ii)); y_ln = f(x_ln);
    p_ln = polyfit(x_ln,y_ln,n(ii)-1);
    e_ln(ii) = max( abs( y-polyval(p_ln,x) ) );

    % Chebyshev nodes
    k = 1:n(ii); x_cn = cos((2*k-1)/2/n(ii)*pi); y_cn = f(x_cn);
    p_cn = polyfit(x_cn,y_cn,n(ii)-1);
    e_cn(ii) = max( abs( y-polyval(p_cn,x) ) );
end
figure
plot(n,e_ln,n,e_cn);
xlabel('no of points'); ylabel('maximum absolute error');
legend('linearly space','chebyshev nodes','location','best')