天天看點

分段線性插值

N=10;            %%定義插值節點的個數
x=-5:10/N:5;     %%依據課本xi=-5+i*h(i=0,1,...N),h=10/N
y=x./(1+x.^4);   %%依據課本(1)插值公式
xi=-5:0.25:5;    %%依據課本xk=-5+0.25k(k=0,1,...40)
m=length(x);     %%求出插值節點向量長度
n=length(y);
s=zeros(1,length(xi));
if m ~= n , error('向量x與y的長度必須一緻');
    %%這裡肯定一緻,隻是為了消除直接選取x,y資料的時候出錯而設定的
end;
%%這裡定義一個長度為計算節點個數的機關向量,可以加快運作速度
for m=1:length(xi)
    for b=1:n
        if (xi(m)>x(b))&&(xi(m)<=x(b+1))
            %%尋找取值區間
            s(m)=(xi(m)-x(b))*(y(b)-y(b+1))/(x(b)-x(b+1))+y(b);
            %%兩點式直線公式
        else
            continue
        end
    end
end




fprintf('\t\t\t插值點X\t\t插值點Y\n');


for i=1:length(x)
    fprintf('I(%2.0f)\t\t%5.3f\t\t%5.3f\n',i,x(i),y(i));
end
fprintf('\t\t\t計算點X\t\t計算點Y\n');
for i=1:length(s)
    fprintf('I(%2.0f)\t\t%5.3f\t\t%5.3f\n',i,xi(i),s(i));
    
end
hold on
%以下是設定繪圖區域的相關指令


plot(x,y,'r*');
ezplot('x./(1+x.^4)',[-5,5]);
plot(xi,s,'k-');
plot(xi,s,'go');


title('分段線性插值');
legend('插值點','原始函數','分段線性插值函數','計算點');
xlabel('x的取值');
ylabel('y的取值');
axis([-5,5,min([min(y) min(f)])-0.2,max([max(y) max(f)])+0.2]);
set(gca,'xtick',[-5:0.25:5]);
set(gca,'ytick',[-10:0.1:10]);
grid on
hold off
           

繼續閱讀