matlab实现Newton、Lagrange、分段线性插值多项式

坐忘1

于 2024-08-15 17:01:22 发布

阅读量327

点赞数 5

文章标签： matlab 开发语言线性代数算法

本文链接：https://blog.csdn.net/2301_81884045/article/details/141226279

版权

例题

Newton插值

clc
% 实验数据 x y 初始赋值，画出原始数据折线图
%%Newton插值、分段线性插值步长
fi=@(i) i.*0.2-1;
fx=@(x) 1./(1+9*x.^2);
i=0:1:10;            %设置i点
x=fi(i);             %插值节点
y=fx(x);             %各插值点函数值
x1=-1:0.01:1;
y1=fx(x1);
scatter(x,y,'k');   %插值节点图
hold on
plot(x1,y1,'r')     %被插函数图
hold on

%%Newton插值
syms t;             %系统变量
n=length(x);        %x长度
A=zeros(n,n);       %差商表
A(:,1)=y(:);        %初始化第一列y
    for  j=2:n
       for i=j:n
           A(i,j)=(A(i,j-1)- A(i-1,j-1))/(x(i)-x(i-j+1));  %逐步构建差商表
       end
    end
y0=A(1,1);                 %插值结果初始化
    for k=2:n
        p=1.0;
        for j=1:k-1
            p=p*(t-x(j));  %计算插值基函数
        end
        y0=y0+A(k,k)*p;    %加入新的基函数   
    end
y0 = collect(y0);               %合并同类项
y0 = vpa(y0,6);                 %截断表达
disp("Newton插值多项式结果：");  
disp(y0);                               %打印最后的Newton插值多项式

% 画图
f0 = subs(y0,'t',x1);
plot(x1,f0,'b');                 %画图
hold on

Lagrange插值

clear
clc
%%Lagrange插值
fi=@(i) cos(((i.*2+1)*pi)/22);    %求插值节点
fx=@(x) 1./(1+9*x.^2);            %求插值节点对应的y值
i=0:1:10;              %设置i点
x=fi(i);               %插值节点
y=fx(x);
x1=-1:0.01:1;
y1=fx(x1);
scatter(x,y,'k');   %画出插值节点图
hold on
plot(x1,y1,'r')     %被插函数图
hold on

%% lagrange插值过程
syms t;                               %系统变量
n = length(x);                        %x长度
f = 0.0;                              %lagrange插值多项式初始化
for i = 1:n                           %循环计算基函数
    l = y(i);
    for j = 1:i-1                     %计算基函数某个项前面一截，i之前的
        l = l*(t-x(j))/(x(i)-x(j));
    end
    for j = i+1:n                     %计算基函数某个项后面一截，i之后的
        l = l*(t-x(j))/(x(i)-x(j));
    end
    f = f + l;                        %得出最后的lagrange插值多项式
    simplify(f);                      %化简多项式
    if(i==n)
        f = collect(f);               %合并同类项
        f = vpa(f,6);                 %截断表达
    end
end
disp("lagrange插值多项式结果：");  
disp(f);                               %打印最后的lagrange插值多项式

%% 画图
f0 = subs(f,'t',x1);
plot(x1,f0,'b');                       %画图
hold on
legend("插值节点","被插函数","lagrange插值");
title("lagrange插值插值结果图")%标题

分段线性插值

%%分段线性插值
syms t;             %系统变量
x2=0.1+x;
n=length(x);        %x长度
  disp("分段线性插值多项式结果："); 
    for i=1:n-1            %输出分段函数表达式
        if x(i)<x2(i) && x2(i)<x(i+1)
            L1=(t-x(i+1))/(x(i)-x(i+1));
            L2=(t-x(i))/(x(i+1)-x(i));
            y1=L1*y(i)+L2*y(i+1);
            y1 = collect(y1);               %合并同类项
            y1 = vpa(y1,6);                 %截断表达
            fprintf('在区间[%.1f,%.1f],多项式为： ',x(i),x(i+1))
            disp(y1);
        elseif x(i)==x1
             y1=y(i);      
        end
    end
y2= interp1(x,y,x1,'linear');
plot(x1,y2,'g');                 %画图
hold on
legend("插值节点","被插函数","Newton插值","分段线性插值");
title("Newton插值、分段线性插值结果图")%标题