方法一:Lagrange插值

function [y]=lagrange(x0,y0,x)
    n=length(x0);y=zeros(1,length(x));
    for k=1:n
        t=ones(1,length(x));
        for i=1:n
            if i~=k
                t=t.*(x-x0(i))/(x0(k)-x0(i));
            end
        end
        y=y+t*y0(k);
    end
end

方法二:Newton插值

function [y]=newton(x0,y0,x)
    n=length(x0);f=zeros(n,n);
    for k=1:n      
        f(k,1)=y0(k);
    end
    for i=2:n
        for k=i:n
            f(k,i)=(f(k,i-1)-f(k-1,i-1))/(x0(k)-x0(k+1-i));  
        end
    end
    y=f(1,1);          
    for k=2:n
        t=ones(1,length(x));
        for j=1:k-1
            t=t.*(x-x0(j)); 
        end
        y=f(k,k)*t+y;
    end
end

方法二*:牛顿前插公式

function [y]=newton_forward(x0,h,xn,Y,x)
    n=(xn-x0)/h+1;f=zeros(n,n);
    for k=1:n      
        f(k,1)=Y(k);
    end
    for i=2:n
        for k=i:n
            f(k,i)=f(k,i-1)-f(k-1,i-1);  
        end
    end
    t=(x-x0)/h; y=f(1,1);    
    for k=2:n
        p=ones(1,length(x));
        for j=1:k-1
            p=p.*(t-j+1)/j;
        end
        y=f(k,k)*p+y;
    end
end

方法三:分段线性插值

function [y]=linear(x0,y0,x)
    n=length(x0);m=length(x);y=zeros(1,m);
    for t=1:m
        for i=1:n-1
            if (x(t)<=x0(i+1))&&(x(t)>=x0(i))
              y(t)=y0(i)*(x(t)-x0(i+1))/(x0(i)-x0(i+1))+y0(i+1)*(x(t)-x0(i))/(x0(i+1)-x0(i));
            end
        end
    end
end

方法四:Hermite插值

function [y]=hermite(x0,y0,m,x)
    n=length(x0);h=length(x);y=zeros(1,h);
    for t=1:h
        for i=1:n-1
            if (x(t)<=x0(i+1))&&(x(t)>=x0(i))
              y(t)=y0(i)*(1+2*((x(t)-x0(i))/(x0(i+1)-x0(i))))*((x(t)-x0(i+1))/(x0(i)-x0(i+1)))^2 ...
                  +y0(i+1)*(1+2*((x(t)-x0(i+1))/(x0(i)-x0(i+1))))*((x(t)-x0(i))/(x0(i+1)-x0(i)))^2 ...
                  +m(i)*(x(t)-x0(i))*((x(t)-x0(i+1))/(x0(i)-x0(i+1)))^2 ...
                  +m(i+1)*(x(t)-x0(i+1))*((x(t)-x0(i))/(x0(i+1)-x0(i)))^2;
            end
        end
    end
end

方法五:三次样条插值(自然边界条件)

function [y]=threesimple(X,Y,M0,x)
     n=length(X); h=zeros(1,n-1);
     for i=1:n-1
         h(i)=X(i+1)-X(i);
     end
     u=zeros(1,n-1);v=zeros(1,n-1);d=zeros(1,n-1);
     for i=2:n-1
         u(i)=h(i-1)/(h(i-1)+h(i)); 
         v(i)=h(i)/(h(i-1)+h(i));
         d(i)=6*(((Y(i+1)-Y(i))/h(i)-(Y(i)-Y(i-1))/h(i-1)))/(h(i-1)+h(i));
     end
     A=zeros(n-2,n-2);b=zeros(n-2,1);M=zeros(1,n);
     for i=1:n-2
         A(i,i)=2; b(i)=d(i+1);
     end
     b(1)=b(1)-u(2)*M0; b(n-2)=b(n-2)-v(n-1)*M0;
     for i=2:n-2
         A(i,i-1)=u(i+1); A(i-1,i)=v(i);
     end
     M(1)=M0;M(n)=M0;M(2:n-1)=A\b;    
     m=length(x); y=zeros(1,m);   
     for t=1:m
        for i=1:n-1
           if (x(t)<=X(i+1))&&(x(t)>=X(i))
              p1=M(i)*(X(i+1)-x(t))^3/(6*h(i));
              p2=M(i+1)*(x(t)-X(i))^3/(6*h(i));
              p3=(Y(i)-M(i)/6*(h(i))^2)*(X(i+1)-x(t))/h(i);
              p4=(Y(i+1)-M(i+1)/6*(h(i))^2)*(x(t)-X(i))/h(i);
              y(t)=p1+p2+p3+p4; 
           end
        end
     end
end

使用范例:

clear
clc
%Q1
x0=[25 40 50 60];y0=[95 75 63 54];x=[70];
disp('answer for Q1');
y=lagrange(x0,y0,x);
y=roundn(y,-1);
disp('lagrange interpolation:y(70)=');disp(y);
%Q2
x0=[1994 1995 1996 1997 1998 1999 2000 2001 2002 2003];
y0=[67.052 68.008 69.803 72.024 73.400 72.063 74.669 74.487 74.065 76.777];
disp('answer for Q2');
X=1990:0.01:2005;
Y=newton(x0,y0,X);plot(x0,y0,'o',X,Y);saveas(gcf,'newton','png');
Y=newton_forward(1994,1,2003,y0,X);plot(x0,y0,'o',X,Y);saveas(gcf,'newton_forward','png');
X=1994:0.01:2003;
Y=linear(x0,y0,X);plot(x0,y0,'o',X,Y);saveas(gcf,'linear','png');
m=[1 2 -7 6 2 -2 2 3 4 -3];
Y=hermite(x0,y0,m,X);plot(x0,y0,'o',X,Y);saveas(gcf,'hermite','png');
x=2010;
y=newton(x0,y0,x);y=roundn(y,-3);
disp('lagrange interpolation:y(2010)=');disp(y);
%Q3
X=[0 1 2 3];Y=[3 5 4 1];M0=1;x=0:0.01:3;
[y_lagrange]=lagrange(X,Y,x);
[y_linear]=linear(X,Y,x);
[y_hermite]=hermite(X,Y,[2 1 1.5 0],x);
[y_threesimple]=threesimple(X,Y,M0,x);
plot(X,Y,'r.',x,y_lagrange,'k-',x,y_linear,'k--',x,y_hermite,'k:',x,y_threesimple,'k-.');
legend('original point','lagrange','linear','hermite','threesimple');
saveas(gcf,'plot of the function','png');

用matlab进行函数插值的几种方法相关推荐

  1. Matlab中数组元素引用——三种方法

    Matlab中数组元素引用--三种方法 1.Matlab中数组元素引用有三种方法 1 2 3 1.下标法(subscripts) 2.索引法(index) 3.布尔法(Boolean) 注意:在使用这 ...

  2. matlab保存所有图,Matlab中图片保存的5种方法

    matlab的绘图和可视化能力是不用多说的,可以说在业内是家喻户晓的. Matlab提供了丰富的绘图函数,比如ez**系类的简易绘图函数,surf.mesh系类的数值绘图函数等几十个.另外其他专业工具 ...

  3. Matlab中计算程序运行时间的几种方法

    平常科研当中,当我们在看文献时,没看到一个优秀的算法时都有想要自己动手编程去实现的愿望,算法好坏可以用代码的运行时间来评估,在MATLAB中大致有以下几种方法来计算程序的运行时间: 1.tic和toc ...

  4. matlab读取cvs文件的几种方法

    matlab读取CVS文件的几种方法: 1,实用csvread()函数 csvread()函数有三种使用方法: 1.M = csvread('filename') 2.M = csvread('fil ...

  5. Matlab 计算均方误差MSE的三种方法

    Matlab 计算均方误差MSE的三种方法 数据说明: ytest 测试集y,真实的y值,是一维数组: ytest_fit 基于测试集 x 预测的y值,是一维数组: test_error 是预测误差. ...

  6. matlab求pi值的三种方法

    https://www.icourse163.org/learn/CSU-1002475002?tid=1450231442#/learn/content?type=detail&id=121 ...

  7. MATLAB求解非线性方程组的五种方法

    MATLAB求解非线性方程组的五种方法 求解线性方程分为两种方法–二分法和迭代法 常见的方法一共有5种 二分法 迭代法 牛顿法 割线法 拟牛顿法 Halley法 使用条件 二分法需要知道两个自变量,分 ...

  8. matlab创建三维数组的三种方法

    在Matlab中习惯性的会将二维数组中的第一维称为"行"第二维称为"列",而对于三维数组的第三位则是习惯性的称为"页".在Matlab中将三 ...

  9. 用matlab设计fir高阶滤波器,用matlab设计fir滤波器的三种方法.doc

    用matlab设计fir滤波器的三种方法.doc 用MATLAB信号处理工具箱进行FIR滤波器设计的三种方法摘要介绍了利用MATLAB信号处理工具箱进行FIR滤波器设计的三种方法程序设计法.FDATO ...

  10. React绑定事件处理函数this的几种方法

    在以类继承的方式定义的组件中,为了能方便地调用当前组件的其他成员方法或属性(如:this.state),通常需要将事件处理函数运行时的 this 指向当前组件实例. 绑定事件处理函数this的几种方法 ...

最新文章

  1. 如何理解李飞飞价值十亿美金的“人文AI”计划 ?
  2. 扫盲行动之九:Vi编辑器的基本使用方法!
  3. 即时通讯音视频开发(三):视频编解码之编码基础
  4. 网络游戏server编程,第一章笔记
  5. UML学习——类图(三)
  6. Hibernate Collection乐观锁定
  7. mysql 并发 锁表_MySQL中的锁(表锁、行锁) 并发控制锁
  8. 利用Serverless Kubernetes和Kaniko快速自动化构建容器镜像
  9. linux多个文件打包命令行,linux命令五十七之tar命令;linux多个文件压缩打包到一个压缩文件...
  10. 英特尔开源技术中心招收Linux内核高手一名
  11. w8fuckcdn 通过扫描全网绕过CDN获取网站IP地址
  12. AXURE中SVG矢量图标的使用方法,以及图标颜色的改变方法
  13. https://github.com/qiangqiang666/demo
  14. 关于Youtube榜单数据的探索,排名第一的视频播放次数已接近90亿次!
  15. hbase 使用lzo_装配HBase LZO
  16. PHP+Mysql 实现最简单的注册登录
  17. 飘雪的夜晚,在那短短的几分钟内,她是我的,然而转过那条街,我们各奔东西
  18. 数值法求解微分博弈问题(〇)——定义
  19. 【十三】景区人流量统计:python日志生成+logstash+kafka+storm+mysql+springBoot+高德地图
  20. 仿真系统rivz进入后不显示机器人,Global Status:Error

热门文章

  1. 基于Python的信用评分卡模型建立和分析,万字阐述,收藏
  2. 教务系统自动评教_FAFU教务管理系统
  3. 点击click触发两次事件解决办法
  4. idea企业破解版安装
  5. 【转载】COMSOL Multiphysics 5.3a 安装教程
  6. 1vcpu等于几核?vcpu是什么意思
  7. 2021年N1叉车司机考试资料及N1叉车司机模拟试题
  8. 【钢结构·技术】国内经典的钢结构建筑BIM应用
  9. 查找交换机IP笨方法
  10. IDEA打包jar包并运行