插值法的程序实现
一实验目的
1.熟悉Matlab编程；
2.学习插值方法及程序设计算法
二实验题目
分别用拉格朗日插值、牛顿插值、自然样条函数对1910、1965、2002的人口进行估算。
三实验原理与理论基础
1.拉格朗日插值算法设计
①利用已知条件得到xi，yi，i=0，1，2，…
②由Lk(x)=((x-x0)*…*(x-x(k-1))*(x-x(k+1))…(x-xn))/ ((xk-x0)*…*(xk-x(k-1))*(xk-x(k+1))…(xk-xn))得出Li(x);
③由Y=y1* L1(x)+…+yn*Ln(x)得出Y关于x的表达式。
④带值计算即可。
2.牛顿插值算法设计
①利用已知条件得到xi，yi，i=0，1，2，…
②利用差商公式
f[x0,…xk]=(f[x0,…,x(k-2),xk]-f[x0,…,x(k-1)])/(xk-x(k-1))各阶差商。
③利用牛顿插值公式
f(x)=f(x0)-f[x0,x1]*(x-x0)+…f[x0,x1,…xn]*(x-x0)*…(x-x(n-1)).
④带值计算即可。
3.自然样条函数算法设计
①利用已知条件得到xi，yi，i=0，1，2，…
②利用公式求出h(i) u(i) k(i);di;
③利用h(i) u(i) k(i);di解出向量M;
④将相关变量带入自然样条表达式中即可。
⑤带值计算即可。
四实验内容
（一）问题重述：
下面给出美国从1920年到1970年的人口表：
2002年的人口。在1910年的实际人口约为91772000，请判断插值计算得到的1965年和2002
年的人口数据准确性是多少？
2用牛顿插值估计：
（1）1965年的人口数；
（2）2002年的人口数。
3用自然样条函数估计在1910，1965和2002年的人口数。
请比较以上三种方法所求值的效果。那一种方法最优？
（二）实验代码：
1.①用matlab编写的拉格朗日插值M文件如下：
function [m,y]=cz5(x)
y1=105711;y2=123203;y3=131669;y4=150697;y5=179323;y6=203212;
x1=1920; x2=1930; x3=1940; x4=1950; x5=1960; x6=1970;
a1=((x-x2)*(x-x3)*(x-x4)*(x-x5)*(x-x6))/((x1-x2)*(x1-x3)*(x1-x4)*(x1-x5)*(x1-x6)); a2=((x-x1)*(x-x3)*(x-x4)*(x-x5)*(x-x6))/((x2-x1)*(x2-x3)*(x2-x4)*(x2-x5)*(x2-x6)); a3=((x-x1)*(x-x2)*(x-x4)*(x-x5)*(x-x6))/((x3-x1)*(x3-x2)*(x3-x4)*(x3-x5)*(x3-x6)); a4=((x-x1)*(x-x2)*(x-x3)*(x-x5)*(x-x6))/((x4-x1)*(x4-x2)*(x4-x3)*(x4-x5)*(x4-x6)); a5=((x-x2)*(x-x3)*(x-x4)*(x-x1)*(x-x6))/((x5-x1)*(x5-x3)*(x5-x2)*(x5-x4)*(x5-x6));
a6=((x-x2)*(x-x3)*(x-x4)*(x-x5)*(x-x1))/((x6-x1)*(x6-x2)*(x6-x3)*(x6-x4)*(x6-x5));
'实际值大约为：'
y=y1*a1+y2*a2+y3*a3+y4*a4+y5*a5+y6*a6;
if x==1910
m=91772;
'估计值大约为：'
else
m=(91772-31872)/31872*y+y;
'估计值大约为：'
End
②用c++编写程序如下：
#include
#include
#include
#include
double H(double x,double a,double p[6]) //求得插值基函数
{
double h=1.0,f=1.0;
for(int i=0;i<6;i++)
{
if(fabs(a-p[i])<0.000001)
continue;
h*=(x-p[i]); //利用累乘求得插值基函数的分子部分
f*=(a-p[i]); //利用累乘求得插值基函数的分母部分
}
return h/f;
}
double L(double x,double p[6],double q[6])//求得拉格朗日插值多项式代入x取值时的函数值
{
double L=0;
for(int i=0;i<6;i++)
{
if(fabs(x-p[i])<0.000001)
return p[i];
L+=q[i]*H(x,p[i],p);// 利用累加实现L=∑(yi*∏[(x-xj)/(xi-xj)])
}
return L;
}
void N(double q[6],double f[6][5])
{
double l=1.0,h=0.0;
for(int i=0;i<6;i++)
{
f[i][0]=q[i];
}
for( i=1;i<6;i++)
{
for(int j=1;j<6;j++)
for(int k=0;k<=j;k++)
;
}
}
double Lagrange(double x,double p[6],double q[6])
{
return L(x,p,q);
}
void main(){
int i;
double x[6],y[6],l=1.0,h=1.0;
double f[6][5];
cout<<"请输入x的取值:"< for(i=0;i<6;i++) cin>>x[i]; cout<<"请依次输入x对应的y的取值:"< for(i=0;i<6;i++) cin>>y[i]; cout<<"Lagrange多插值项式:"< cout<<"L=∑(yi*∏[(x-xj)/(xi-xj)]),其中(i=0,1,2,3,4,5;j=0,1,2,3,4,5;j≠ i)"< cout< cout<<"Lagrange插值多项式在1910，1965，2002处的函数值:"< cout<<"当x值为1910时"< cout<<"当x值为1965时"< cout<<"当x值为2002时"< cout< } 2.用matlab编写的牛顿插值M文件如下： function [y]=niuden5(x) y1=[105711 123203 131669 150697 179323 203212]; x0=[1920 1930 1940 1950 1960 1970]; for i=2:6 y2(i)=(y1(1)-y1(i))/(x0(1)-x0(i)); end for i=3:6 y3(i)=(y2(2)-y2(i))/(x0(2)-x0(i)); end for i=4:6 y4(i)=(y3(3)-y3(i))/(x0(3)-x0(i)); end for i=5:6 y5(i)=(y4(4)-y4(i))/(x0(4)-x0(i)); end y6(6)=(y5(5)-y5(6))/(x0(5)-x0(6)); y=y1(1)+y2(2)*(x-x0(1))+y3(3)*(x-x0(1))*(x-x0(2))+y4(4)*(x-x0(1))*(x-x0( 2))*(x-x0(3))+y5(5)*(x-x0(1))*(x-x0(2))*(x-x0(3))*(x-x0(4))+y6(6)*(x-x0( 1))*(x-x0(2))*(x-x0(3))*(x-x0(4))*(x-x0(5)); 3.用matlab编写的自然样条函数M文件如下： function []=yangtiaochazhi(x) y1=[105711 123203 131669 150697 179323 203212]; x0=[1920 1930 1940 1950 1960 1970]; for i=1:5 h(i)=x0(i+1)-x0(i); end for i=1:4 u(i)=h(i)/(h(i)+h(i+1)); end u(5)=0; k(1)=0; for i=2:5 k(i)=h(i)/(h(i-1)+h(i)); end d(1)=0;d(6)=0; for i=2:3 y21(i)=(y1(1)-y1(i))/(x0(1)-x0(i)); end for i=3:4 y22(i)=(y1(2)-y1(i))/(x0(2)-x0(i)); end for i=4:5 y23(i)=(y1(3)-y1(i))/(x0(3)-x0(i)); end for i=5:6 y24(i)=(y1(4)-y1(i))/(x0(4)-x0(i)); end d(1)=0;d(6)=0; d(2)=(y21(2)-y21(3))/(x0(2)-x0(3)); d(3)=(y22(3)-y22(4))/(x0(3)-x0(4)); d(4)=(y23(4)-y23(5))/(x0(4)-x0(5)); d(5)=(y24(5)-y24(6))/(x0(5)-x0(6)); A=[2 k(1) 0 0 0 0;u(1) 2 k(2) 0 0 0;0 u(2) 2 k(3) 0 0;0 0 u(3) 2 k(4) 0;0 0 0 u(4) 2 k(5);0 0 0 0 u(5) 2]; B=[d(1);d(2);d(3);d(4);d(5);d(6)]; M=[]; M=A\B; echo if x<=1930 y=M(1)*((x0(2)-x)*(x0(2)-x)*(x0(2)-x))/(6*h(1))+M(2)*((x-x0(1))*(x-x0(1))*(x-x0(1)) )/(6*h(1))+(y1(1)-(M(1)*h(1)*h(1))/6)*((x0(2)-x)/h(1))+(y1(2)-(M(2)*h(1)*h(1))/6)*( (x-x0(1))/h(1)); elseif x<=1940 y=M(2)*((x0(3)-x)*(x0(3)-x)*(x0(3)-x))/(6*h(2))+M(3)*((x-x0(2))*(x-x0(2))*(x-x0(2)) )/(6*h(2))+(y1(2)-(M(2)*h(2)*h(2))/6)*((x0(3)-x)/h(2))+(y1(3)-(M(3)*h(2)*h(2))/6)*( (x-x0(2))/h(2)); elseif x<=1950 y=M(3)*((x0(4)-x)*(x0(4)-x)*(x0(4)-x))/(6*h(3))+M(4)*((x-x0(3))*(x-x0(3))*(x-x0(3)) )/(6*h(3))+(y1(3)-(M(3)*h(3)*h(3))/6)*((x0(4)-x)/h(3))+(y1(4)-(M(4)*h(3)*h(3))/6)*( (x-x0(3))/h(3)); elseif x<=1960 y=M(4)*((x0(5)-x)*(x0(5)-x)*(x0(5)-x))/(6*h(4))+M(5)*((x-x0(4))*(x-x0(4))*(x-x0(4)) )/(6*h(4))+(y1(4)-(M(4)*h(4)*h(4))/6)*((x0(5)-x)/h(4))+(y1(5)-(M(5)*h(4)*h(4))/6)*( (x-x0(4))/h(4)); else y=M(5)*((x0(6)-x)*(x0(6)-x)*(x0(6)-x))/(6*h(5))+M(6)*((x-x0(5))*(x-x0(5))*(x-x0(5)) )/(6*h(5))+(y1(5)-(M(5)*h(5)*h(5))/6)*((x0(6)-x)/h(5))+(y1(6)-(M(6)*h(5)*h(5))/6)*( (x-x0(5))/h(5)); end y 五实验结果 1. ①5次拉格朗日插值在matlab命令窗口中输入如下代码： >>[m,y]=cz5(1910) >>[m,y]=cz5(1965) >>[m,y]=cz5(2002) 结果为： ans = 实际值大约为： ans = 估计值大约为： m = 91772 y = 31872 ans = 实际值大约为： ans = 估计值大约为： m = 5.559574703016167e+005 y = 1.930815117187500e+005 ans = 实际值大约为： ans = 估计值大约为： m = 7.526371798548073e+004 y = 2.613874841600098e+004 ②c++中的运行： 请输入x的取值: 1920 1930 1940 1950 1960 1970 请依次输入x对应的y的取值: 105711 123203 131669 150697 179323 203212 Lagrange多插值项式: L=∑(yi*∏[(x-xj)/(xi-xj)]),其中(i=0,1,2,3,4,5;j=0,1,2,3,4,5;j≠i) Lagrange插值多项式在1910，1965，2002处的函数值: 当x值为1910时31872 当x值为1965时193082 当x值为2002时26138.7 2. 牛顿插值在matlab命令窗口中输入如下代码： >>niuden5(1965),niuden5(2002) ----------------------------------------------------------------------------------------------------------------- 结果为： ans = 1.930815117187500e+005 ans = 2.613874841600133e+004 3. 样条插值在matlab命令窗口中输入如下代码： >> yangtiaochazhi(1965),yangtiaochazhi(2002) 结果为： y = 1.913745478468899e+005 y = 2.880973086315789e+005 六实验结果分析与小结 1. 5次拉格朗日插值与牛顿插值在估计人口时比较接近，但与真实值比起来误差较大，而用样条插值在区间内部估值比较准确，但在区间外估值误差相对 较大哦。 2.在matlab中数组名不能再用来作变量了哦。忽略此处.. 2 三次样条插值程序 三次样条插值利用方案二（求解固支样条或压紧样条） 按照要求要起点和终点的一阶导数值已知， 可得关于01,,.....,n M M M 的严格对角占优势的三对角方程组 然后利用三对角法（追赶法）解此线性方程组。 （1）编写M 文件，并保存文件名scfit.m % x,y 分别为n 个节点的横坐标和纵坐标值组成的向量 % dx0和dxn 分别为S 的导数在x0和xn 处的值，即m 0和m n n=length(x)-1; h=diff(x); d=diff(y)./h; a=h(2:n-1); b=2*(h(1:n-1)+h(2:n)); c=h(2:n); u=6*diff(d); b(1)=b(1)-h(1)/2; u(1)=u(1)-3*(d(1)-dx0); b(n-1)=b(n-1)-h(n)/2; u(n-1)=u(n-1)-3*(dxn-d(n)); %追赶法部分 for k=2:n-1 temp=a(k-1)/b(k-1); b(k)=b(k)-temp*c(k-1); u(k)=u(k)-temp*u(k-1); end m(n)=u(n-1)/b(n-1); for k=n-2:-1:1 m(k+1)=(u(k)-c(k)*m(k+2))/b(k); end %求S K1，S K2，S K3，S K4 m(1)=3*(d(1)-dx0)/h(1)-m(2)/2; m(n+1)=3*(dxn-d(n))/h(n)-m(n)/2; for k=0:n-1 00 ()S x m '=()n n S x m '=0011111111212212n n n n n n M d M d M d M d μλμλ----??????????????????????=?????????????????????????? 牛顿插值法 插值法是利用函数f (x)在某区间中若干点的函数值，作出适当的特定函数，在这些点上取已知值，在区间的其他点上用这特定函数的值作为函数f (x)的近似值。如果这特定函数是多项式，就称它为插值多项式。当插值节点增减时全部插值基函数均要随之变化，这在实际计算中很不方便。为了克服这一缺点，提出了牛顿插值。牛顿插值通过求各阶差商，递推得到的一个公式： f(x)=f[x0]+f[x0,x1](x-x0)+f[x0,x1,x2](x-x0)(x-x1)+...f[x0,...xn](x-x0 )...(x-xn-1)+Rn(x)。 插值函数 插值函数的概念及相关性质[1] 定义：设连续函数y-f(x) 在区间[a,b]上有定义，已知在n+1个互异的点 x0,x1,…xn上取值分别为y0,y1,…yn （设a≤ x1≤x2……≤xn≤b)。若在函数类中存在以简单函数P(x) ，使得P(xi)=yi,则称P(x) 为f(x)的插值函数. 称x1,x2,…xn 为插值节点，称[a,b]为插值区间。 定理：n次代数插值问题的解存在且唯一。 牛顿插值法C程序 程序框图#include Several Typical Interpolation in Matlab Lagrange Interpolation Supposing： If x=175, while y=? Solution： Lagrange Interpolation in Matlab: function y=lagrange(x0,y0,x); n=length(x0);m=length(x); for i=1:m z=x(i); s=0.0; for k=1:n p=1.0; for j=1:n if j~=k p=p*(z-x0(j))/(x0(k)-x0(j)); end end s=p*y0(k)+s; end y(i)=s; end input： x0=[144 169 225] y0=[12 13 15] y=lagrange(x0,y0,175) obtain the answer： x0 = 144 169 225 y0 = 12 13 15 y = 13.2302 Spline Interpolation Solution ： Input x=[ 1 4 9 6]；y=[ 1 4 9 6]；x=[ 1 4 9 6]；pp=spline(x,y) pp = form: 'pp' breaks: [1 4 6 9] coefs: [3x4 double] pieces: 3 order: 4 dim: 1 output : pp.coefs ans = -0.0500 0.5333 -0.8167 1.0000 -0.0500 0.0833 1.0333 2.0000 -0.0500 -0.2167 0.7667 4.0000 It shows the coefficients of cubic spline polynomial , so: S （x ）=, 169,3)9(1484.0)9(0063.0)9(0008.0,94,2)4(2714.0)4(0183.0)4(0008 .0, 41,1)1(4024.0)1(0254.0)1(0008.0232 3 23≥≤+-+---≥≤+-+---≥≤+-+---x x x x x x x x x x x x Newton’s Interpolation Resolve 65 Solution: Newton’s Interpolation in matlab ： function yi=newint(x,y,xi); n=length(x); ny=length(y); if n~=ny error end Y=zeros(n);Y(:,1)=y'; 沈阳航空航天大学 数学软件课程设计 （设计程序） 题目三次样条插值函数 班级 / 学号 学生姓名 指导教师 沈阳航空航天大学 课程设计任务书 课程名称数学软件课程设计 院（系）理学院专业信息与计算科学 班级学号姓名 课程设计题目三次样条插值函数 课程设计时间: 2010 年12月20日至2010 年12月31日 课程设计的内容及要求： 1.三次样条插值函数 给出函数在互异点处的值分别为。 (1)掌握求三次样条插值函数的基本原理； (2)编写程序求在第一边界条件下函数的三次样条插值函数； (3)在区间上取n=10，20，分别用等距节点对函数 作三次样条插值函数，利用(1)的结果画出插值函数的图形，并在该图形界面中同时画出的图形。 [要求] 1.学习态度要认真，要积极参与课程设计，锻炼独立思考能力； 2.严格遵守上机时间安排； 3.按照MATLAB编程训练的任务要求来编写程序； 4.根据任务书来完成课程设计论文； 5.报告书写格式要求按照沈阳航空航天大学“课程设计报告撰写规范”； 6.报告上交时间：课程设计结束时上交报告； 7.严谨抄袭行为。 指导教师年月日负责教师年月日学生签字年月日 沈阳航空航天大学 课程设计成绩评定单 课程名称数学软件课程设计 院（系）理学院专业信息与计算科学课程设计题目三次样条插值函数 学号姓名 指导教师评语： 课程设计成绩 指导教师签字 年月日 目录 一正文 (1) 1问题分析 (1) 1.1 题目 (1) 1.2 分析 (1) 2 研究方法原理 (1) 2.1 求三次样条插值多项式，算法组织 (1) 3 算例结果 (3) 二总结 (7) 参考文献 (8) 附录 (9) 源程序: (9) 程序1 (9) 程序2 (10) 程序3 (12) 程序 4 (12) 探索实验5 插值法 一、 实验目的 了解插值问题及其适用的场合，理解并掌握常用的插值算法的构造和计算，了解差商概念、Runge 现象及样条插值方法，学习用计算机求近似函数的一些科学计算方法和简单的编程技术。 二、概念与结论 1. 插值问题与插值函数： 由实验或测量的方法得到所求函数 y=f(x) 在互异点x 0 , x 1, ... , x n 处的值 y 0 , y 1 , … , y n ,构造一个简单函数 ?(x) 作为函数 y=f(x) 的近似表达式 y= f(x) ≈ ?(x) 使 ?(x 0)=y 0 , ?(x 1)=y 1 , ?, ?(x n )=y n , (1) 这类问题称为插值问题。 f(x) 称为被插值函数，?(x) 称为插值函数， x 0 , x 1, ... , x n 称为插值节点。(1)式称为插值条件。 常用的插值函数是多项式函数。且当n=1时是称为线性插值，n=2时称为Simpson 插值或抛物线插值。 2．插值定理： 假设x 0 ,x 1,…,x n 是n+1个互异节点,函数f(x)在这组节点的值f(x k )(k=0,1,…,n)是给定的，那么存在唯一的n 次次多项式p n (x)满足 p n (x k )=f(x k ), k=0,1,…,n 3.插值的截断误差 设?n (x)是过点x 0 ，x 1 ，x 2 ，…x n 的 n 次插值多项式， f (n+1)(x)在（a ，b ）上存在，其中[a ，b]是包含点x 0 ，x 1 ，x 2 ，…,x n 的任一区间，则对任意给定的x ∈[a ，b]，总存在一点ξ∈（a ，b ）（依赖于x ）使 其中 ωn+1 (x)=(x –x 0) (x - x 1)…(x -x n ) ，f (n+1)(ξ) 是f(x)的n+1阶微商在 ξ 的值。 4. 差商： 给定一个函数表 x | x 0 x 1 ….... x n -------- --------------------------------------------------------- y | y 0 ,y 1 ……. y n 其中当i ≠j 时 ,x i ≠x j 记 f[x i ]=f(x i ) ,定义f(x)关于x i ,x j 的一 阶差商 一般的, f(x)关于x i ,x i+1,…,x i+k 的k 阶差商定义为: ) ()! 1() ()()()(1) 1(x n x x f x R n n n n f +++= -=ωξ?j i j i j i x x x f x f x x f --= ][][],[ Matlab实现数值分析插值及积分 摘要： 数值分析(numerical analysis)是研究分析用计算机求解数学计算问题的数值计算方法及其理论的学科，是数学的一个分支，它以数字计算机求解数学问题的理论和方法为研究对象。在实际生产实践中，常常将实际问题转化为数学模型来解决，这个过程就是数学建模。学习数值分析这门课程可以让我们学到很多的数学建模方法。 分别运用matlab数学软件编程来解决插值问题和数值积分问题。题目中的要求是计算差值和积分，对于问题一，可以分别利用朗格朗日插值公式，牛顿插值公式，埃特金逐次线性插值公式来进行编程求解，具体matlab代码见正文。编程求解出来的结果为：=+。 其中Aitken插值计算的结果图如下： 对于问题二，可以分别利用复化梯形公式，复化的辛卜生公式，复化的柯特斯公式编写程序来进行求解，具体matlab代码见正文。编程求解出来的结果为： 0.6932 其中复化梯形公式计算的结果图如下： 问题重述 问题一：已知列表函数 表格 1 分别用拉格朗日，牛顿，埃特金插值方法计算。 问题二：用复化的梯形公式，复化的辛卜生公式，复化的柯特斯公式计算积分，使精度小于5。 问题解决 问题一：插值方法 对于问题一，用三种差值方法：拉格朗日，牛顿，埃特金差值方法来解决。 一、拉格朗日插值法： 拉格朗日插值多项式如下： 首先构造1+n 个插值节点n x x x ,,,10 上的n 插值基函数，对任一点i x 所对应的插值基函数 )(x l i ，由于在所有),,1,1,,1,0(n i i j x j +-=取零值，因此)(x l i 有因子 )())(()(110n i i x x x x x x x x ----+- 。又因)(x l i 是一个次数不超过n 的多项式，所以只 可能相差一个常数因子，固)(x l i 可表示成： )())(()()(110n i i i x x x x x x x x A x l ----=+- 利用1)(=i i x l 得： 插值法和拟合实验报告 一、 实验目的 1.通过进行不同类型的插值，比较各种插值的效果，明确各种插值的优越性； 2.通过比较不同次数的多项式拟合效果，了解多项式拟合的原理； 3.利用matlab 编程，学会matlab 命令； 4.掌握拉格朗日插值法； 5.掌握多项式拟合的特点和方法。 二、 实验题目 1.、插值法实验 将区间[-5,5]10等分，对下列函数分别计算插值节点 k x 的值，进行不同类型 的插值，作出插值函数的图形并与)(x f y =的图形进行比较： ;11)(2x x f += ;a r c t a n )(x x f = .1)(42 x x x f += (1) 做拉格朗日插值； (2) 做分段线性插值； (3) 做三次样条插值. 2、拟合实验 给定数据点如下表所示： 分别对上述数据作三次多项式和五次多项式拟合，并求平方误差，作出离散函数 ),(i i y x 和拟合函数的图形。 三、 实验原理 1.、插值法实验 ∏∑∏∏∏∑∑≠==≠=≠=≠=+-==--= =-= ==-=-=----==++==j i j j i i i i i n i i n n j i j j n j i j j i i n j i j j n i i i n i i n n n o i n i i n x x x x x y x l x L x x c n i x x c x x x c x x x x x x x x c y x l x L y x l y x l y x l x L ,00 ,0,0,01100 00 )(l )()() (1 ,1,0, 1)()(l ) ()())(()()()()()()()(， 故， 得 再由，设 2、拟合实验 拉格朗日插值的调用函数 function y=lagrange(x0,y0,x) n=length(x0);m=length(x); for i=1:m z=x(i); L=0.0; for j=1:n T=1.0; for k=1:n if k~=j T=T*(z-x0(k))/(x0(j)-x0(k)); end end L=T*y0(j)+L; end y(i)=L; end 四个图在一起： x=[-1:0.05:1]; y=1./(1+25*x.^2); x0=[-1:0.4:1]; y0=1./(1+25*x0.^2); y1=lagrange(x0,y0,x); x0=[-1:0.2:1]; y0=1./(1+25*x0.^2); y2=lagrange(x0,y0,x); x0=[-1:0.1:1]; y0=1./(1+25*x0.^2); y3= lagrange(x0,y0,x); plot(x,y,'-r') hold on plot(x,y1,'-b',x,y2,'-r',x,y3,'-r') l5和fx在一起： x=[-1:0.05:1]; y=1./(1+25*x.^2); x0=[-1:0.4:1]; y0=1./(1+25*x0.^2); y1=lagrange(x0,y0,x); plot(x,y,'-r') hold on plot(x,y1,'-b') l10和fx在一起： x=[-1:0.05:1]; y=1./(1+25*x.^2); x0=[-1:0.2:1]; y0=1./(1+25*x0.^2); y2= lagrange(x0,y0,x); plot(x,y,'-r') hold on plot(x,y2,'-b') l20和fx在一起： x=[-1:0.05:1]; y=1./(1+25*x.^2); x0=[-1:0.1:1]; y0=1./(1+25*x0.^2); y3= lagrange(x0,y0,x); plot(x,y,'-r') hold on plot(x,y3,'-b') 数值计算方法作业 实验4.3 三次样条差值函数 实验目的： 掌握三次样条插值函数的三弯矩方法。 实验函数: dt e x f x t ? ∞ -- = 2 221)(π 实验内容： （1） 编程实现求三次样条插值函数的算法，分别考虑不同的边界条件； （2） 计算各插值节点的弯矩值； （3） 在同一坐标系中绘制函数f(x)，插值多项式，三次样条插值多项式的曲线 比较插值结果。 实验4.5 三次样条差值函数的收敛性 实验目的： 多项式插值不一定是收敛的，即插值的节点多，效果不一定好。对三次样条插值函数如何呢？理论上证明三次样条插值函数的收敛性是比较困难的，通过本实验可以证明这一理论结果。 实验内容： 按照一定的规则分别选择等距或非等距的插值节点，并不断增加插值节点的个数。 实验要求： （1） 随着节点个数的增加，比较被逼近函数和三样条插值函数的误差变化情 况，分析所得结果并与拉格朗日插值多项式比较； （2） 三次样条插值函数的思想最早产生于工业部门。作为工业应用的例子，考 虑如下例子：某汽车制造商根据三次样条插值函数设计车门曲线，其中一 算法描述： 拉格朗日插值： 错误!未找到引用源。 其中错误!未找到引用源。是拉格朗日基函数，其表达式为：() ∏ ≠=--=n i j j j i j i x x x x x l 0) ()( 牛顿插值： ) )...()(](,...,,[.... ))(0](,,[)0](,[)()(1102101210100----++--+-+=n n n x x x x x x x x x x f x x x x x x x f x x x x f x f x N 其中????? ?? ?? ?????? --=--= --= -)/(]),...,[],...,[(]...,[..],[],[],,[)()(],[01102110x x x x x f x x x f x x x f x x x x f x x f x x x f x x x f x f x x f n n n n i k j i k j k j i j i j i j i 三样条插值： 所谓三次样条插值多项式Sn(x)是一种分段函数，它在节点Xi(a 佛山科学技术学院 实 验 报 告 课程名称 数值分析 实验项目 插值法与数据拟合 专业班级 机械工程 姓 名 余红杰 学 号 10 指导教师 陈剑 成 绩 日 期 月 日 一、实验目的 1、学会Lagrange 插值、牛顿插值和三次样条插值等基本插值方法； 2、讨论插值的Runge 现象 3、学会Matlab 提供的插值函数的使用方法，会用这些函数解决实际问题。 二、实验原理 1、拉格朗日插值多项式 2、牛顿插值多项式 3、三次样条插值 三、实验步骤 1、用MATLAB 编写独立的拉格朗日插值多项式函数 2、用MATLAB 编写独立的牛顿插值多项式函数 3、用MATLAB 编写独立的三次样条函数（边界条件为第一、二种情形） 4、已知函数在下列各点的值为： 根据步骤1,2,3编好的程序，试分别用4次拉格朗日多项式4()L x 、牛顿插值多项式4()P x 以及三次样条函数()S x (自然边界条件)对数据进行插值，并用图给出 ｛(,),0.20.08,0,1,2, ,10i i i x y x i i =+=｝，4()L x 、4()P x 和()S x 。 5、在区间[-1,1]上分别取10,20n =用两组等距节点对龙格函数 2 1 (),(11)125f x x x = -≤≤+作多项式插值，对不同n 值，分别画出插值函数及()f x 的图形。 6、下列数据点的插值 可以得到平方根函数的近似，在区间[0,64]上作图。 （1）用这9个点作8次多项式插值8()L x 。 （2）用三次样条（第一边界条件）程序求()S x 。 7、对于给函数2 1 ()125f x x = +在区间[-1,1]上取10.2(0,1, ,10)i x i i =-+=，试求3次 曲线拟合，试画出拟合曲线并打印出方程，与第5题的结果比较。 四、实验过程与结果： 1、Lagrange 插值多项式源代码： function ya=lag(x,y,xa) %x 所有已知插值点 %y 插值点对应函数值 %xa 所求点，自变量 %ya 所求点插值估计量 ya=0; mu=1; %初始化 %循环方式求L 系数，并求和： for i = 1:length(y) for j = 1:length(x) if i ~= j mu = mu * (xa - x(j) ) / ( x(i) - x(j) ); else continue end end ya = ya + y(i) * mu ; mu = 1; end 2、Newton 源代码： function ya = newton(x,y,xa) %x 所有已知插值点 %y 插值点对应函数值 %xa 所求点，自变量 %ya 所求点插值估计量 %建立系数零矩阵D 及初始化： 数值分析作业 姓名王建忠 学号132080202006 学院能源与动力工程 专业机械电子工程 2013年12月16日 https://www.doczj.com/doc/a013035293.html,grange插值多项式程序 function f=nalagr(x,y,xx) %x为节点向量 %y为节点函数值 %xx是插值点 syms s if(length(x)==length(y)) n=length(x); else disp('x和y的维数不相等!'); return; end f=0.0; for(i=1:n) l=y(i); for(j=1:i-1) l=l*(s-x(j))/(x(i)-x(j)); end; for(j=i+1:n) l=l*(s-x(j))/(x(i)-x(j));%计算拉格朗日基函数end; f=f+l;%计算拉格朗日插值函数 simplify(f); if(i==n) if(nargin==3) f=subs(f,'s');%计算插值点的函数值else f=collect(f);%将插值多项式展开 f=vpa(f,6);%将插值多项式的系数化成6位精度的小数 end end end >>x=[-2,-1,0,1];%已知节点向量y=[3,1,1,6];%节点函数值向量 f=nalagr(x,y) f= 0.5*s^3+ 2.5*s^2+ 2.0*s+ 1.0 >>f=nalagr(x,y,0) f=1 >> 2.牛顿插值多项式程序 function[p2,z]=newTon(x,y,t) %输入参数中x,y为元素个数相等的向量，t为待估计的点，可以为数字或向量。%输出参数中p2为所求得的牛顿插值多项式，z为利用多项式所得的t的函数值。 n=length(x); chaS(1)=y(1); for i=2:n x1=x;y1=y; x1(i+1:n)=[]; y1(i+1:n)=[]; n1=length(x1); s1=0; for j=1:n1 t1=1; for k=1:n1 if k==j continue; else t1=t1*(x1(j)-x1(k)); end end s1=s1+y1(j)/t1; end chaS(i)=s1; end b(1,:)=[zeros(1,n-1)chaS(1)]; cl=cell(1,n-1); for i=2:n u1=1; for j=1:i-1 u1=conv(u1,[1-x(j)]); cl{i-1}=u1; end cl{i-1}=chaS(i)*cl{i-1}; b(i,:)=[zeros(1,n-i),cl{i-1}]; end p2=b(1,:); for j=2:n p2=p2+b(j,:); end if length(t)==1 rm=0; 课题一：拉格朗日插值法 1.实验目的 1．学习和掌握拉格朗日插值多项式。 2.运用拉格朗日插值多项式进行计算。 2.实验过程 作出插值点（1.00，0.00），（-1.00，-3.00），（2.00，4.00）二、算法步骤 已知：某些点的坐标以及点数。 输入：条件点数以及这些点的坐标。 输出：根据给定的点求出其对应的拉格朗日插值多项式的值。 3.程序流程： （1）输入已知点的个数； （2）分别输入已知点的X坐标； （3）分别输入已知点的Y坐标； 程序如下： #include 插值算法*/ { int i,j; float *a,yy=0.0; /*a作为临时变量，记录拉格朗日插值多项*/ a=(float*)malloc(n*sizeof(float)); for(i=0;i<=n-1;i++) { a[i]=y[i]; for(j=0;j<=n-1;j++) if(j!=i) a[i]*=(xx-x[j])/(x[i]-x[j]); yy+=a[i]; } free(a); return yy; } int main() { int i; int n; float x[20],y[20],xx,yy; printf("Input n:"); scanf("%d",&n); if(n<=0) { printf("Error! The value of n must in (0,20)."); getch();return 1; } for(i=0;i<=n-1;i++) { printf("x[%d]:",i); scanf("%f",&x[i]); } printf("\n"); for(i=0;i<=n-1;i++) { printf("y[%d]:",i);scanf("%f",&y[i]); } printf("\n"); printf("Input xx:"); scanf("%f",&xx); yy=lagrange(x,y,xx,n); printf("x=%f,y=%f\n",xx,yy); getch(); } 举例如下：已知当x=1，-1,2时f(x)=0，-3,4，求f(1.5)的值。 Matlab中插值函数汇总和使用说明 MATLAB中的插值函数为interp1，其调用格式 为： yi= interp1(x,y,xi,'method') 其中x，y为插值点，yi为在被插值点xi处的插值结果；x,y为向量， 'method'表示采用的插值方法，MATLAB提供的插值方法有几种： 'method'是最邻近插值， 'linear'线性插值； 'spline'三次样条插值； 'cubic'立方插值．缺省时表示线性插值 注意：所有的插值方法都要求x是单调的，并且xi不能够超过x的范围。 例如：在一天24小时内，从零点开始每间隔2小时测得的环境温度数据分别为 12，9，9，10，18 ，24，28，27，25，20，18，15，13， 推测中午12点（即13点）时的温度． x=0:2:24; y=[12 9 9 10 18 24 28 27 25 20 18 15 13]; a=13; y1=interp1(x,y,a,'spline') 结果为： 27.8725 若要得到一天24小时的温度曲线，则： xi=0:1/3600:24; yi=interp1(x,y,xi, 'spline'); plot(x,y,'o' ,xi,yi) 命令1 interp1 功能一维数据插值（表格查找）。该命令对数据点之间计算内插值。它找出一元函数f(x)在中间点的数值。其中函数f(x)由所给数据决定。 x：原始数据点 Y：原始数据点 xi：插值点 Yi：插值点 格式 (1)yi = interp1(x,Y,xi) 返回插值向量yi，每一元素对应于参量xi，同时由向量x 与Y 的内插值决定。参量x 指定数据Y 的点。 若Y 为一矩阵，则按Y 的每列计算。yi是阶数为length(xi)*size(Y,2)的输出矩阵。(2)yi = interp1(Y,xi) 假定x=1:N，其中N 为向量Y 的长度，或者为矩阵Y 的行数。 (3)yi = interp1(x,Y,xi,method) 用指定的算法计算插值： ’nearest’：最近邻点插值，直接完成计算； ’linear’：线性插值（缺省方式），直接完成计算； ’spline’：三次样条函数插值。对于该方法，命令interp1 调用函数spline、ppval、mkpp、umkpp。这些命令生成一系列用于分段多项式操作的函数。命令spline 用它们执行三次样条函数插值； ’pchip’：分段三次Hermite 插值。对于该方法，命令interp1 调用函数pchip，用于对 实验二插值法 1、实验目的： 1、掌握直接利用拉格郎日插值多项式计算函数在已知点的函数值；观察拉格郎日插值的龙格现象。 2、了解Hermite插值法、三次样条插值法原理，结合计算公式，确定函数值。 2、实验要求: 1)认真分析题目的条件和要求，复习相关的理论知识，选择适当的解决方案和算法； 2)编写上机实验程序，作好上机前的准备工作； 3)上机调试程序，并试算各种方案，记录计算的结果（包括必要的中间结果）； 4)分析和解释计算结果； 5)按照要求书写实验报告； 3、实验内容： 1) 用拉格郎日插值公式确定函数值；对函数f(x)进行拉格郎日插值，并对f(x)与插值多项式的曲线作比较。 已知函数表：（0.56160,0.82741）、（0.56280,0.82659）、（0.56401,0.82577）、（0.56521,0.82495）用三次拉格朗日插值多项式求x=0.5635时函数近似值。 2) 求满足插值条件的插值多项式及余项 1) 4、题目：插值法 5、原理： 拉格郎日插值原理： n次拉格朗日插值多项式为：L n (x)=y l (x)+y 1 l 1 (x)+y 2 l 2 (x)+…+y n l n (x) n=1时，称为线性插值， L 1(x)=y (x-x 1 )/(x -x 1 )+y 1 (x-x )/(x 1 -x )=y +(y 1 -x )(x-x )/(x 1 -x ) n=2时，称为二次插值或抛物线插值， L 2(x)=y (x-x 1 )(x-x 2 )/(x -x 1 )/(x -x 2 )+y 1 (x-x )(x-x 2 )/(x 1 -x )/(x 1 -x 2 )+y 2 (x -x 0)(x-x 1 )/(x 2 -x )/(x 2 -x 1 ) n=i时， Li= （X-X0）……(X-X i-1)(x-x i+1) ……(x-x n) （X-X0）……(X-X i-1)(x-x i+1) ……(x-x n) 6、设计思想： 拉格朗日插值法是根据n + 1个点x0, x1, ... x n(x0 < x1 < ... x n)的函数值f (x0), f (x1) , ... , f (x n)推出n次多項式p(x)，然后n次多項式p (x)求出任意的点x对应的函数值f (x)的算法。 7、对应程序： 1 ) 三次拉格朗日插值多项式求x=0.5635时函数近似值 #include"stdio.h" #define n 5 void main() { int i,j; float x[n],y[n]; float x1; float a=1; float b=1; float lx=0; printf("\n请输入想要求解的X：\n x="); scanf("%f",&x1); printf("请输入所有点的横纵坐标：\n"); for(i=1;i 题目一：多项式插值 某气象观测站在8：00（AM ）开始每隔10分钟对天气作如下观测，用三次多项式插值函数（Newton ）逼近如下曲线，插值节点数据如上表，并求出9点30分该地区的温度（x=10）。 二、数学原理 假设有n+1个不同的节点及函数在节点上的值（x 0，y 0），……（x n ，y n ），插值多项式有如下形式： ）（） ）（（）（）（）（n 10n 102010n x -x )(x -x x -x x P x x x x x x -??-+??+-++=αααα （1） 其中系数i α（i=0,1,2……n ）为特定系数，可由插值样条i i n y x P =） （（i=0,1,2……n ）确定。 根据均差的定义，把x 看成[a,b]上的一点,可得 f(x)= f （0x ）+f[10x x ，]（0x -x ） f[x, 0x ]= f[10x x ，]+f[x,10x x ，] （1x -x ） …… f[x, 0x ，…x 1-n ]= f[x, 0x ，…x n ]+ f[x, 0x ，…x n ](x-x n ) 综合以上式子，把后一式代入前一式，可得到： f(x)= f[0x ]+f[10x x ，]（0x -x ）+ f[210x x x ，，]（0x -x ）（1x -x ）+ …+ f[x, 0x ，…x n ]（0x -x ）…(x-x 1-n )+ f[x, 0x ，…x n ,x ]） （x 1n +ω= N n （x ）+） （x n R 其中 N n （x ）= f[0x ]+f[10x x ，]（0x -x ）+ f[210x x x ，，]（0x -x ）（1x -x ）+ …+ f[x, 0x ，…x n ]（0x -x ）…(x-x 1-n ) （2） ）（x n R = f(x)- N n （x ）= f[x, 0x ， (x) n ,x ]） （x 1n +ω （3） ） （x 1n +ω=（0x -x ）…(x-x n ) Newton 插值的系数i α（i=0,1,2……n ）可以用差商表示。一般有 f k =α[k 10x x x ??，] （k=0,1,2，……，n ） （4） 把（4）代入（1）得到满足插值条件N ）（） （i i n x f x =（i=0,1,2，……n ）的n 次Newton 插值多项式 N n （x ）=f （0x ）+f[10x x ，]（1x -x ）+f[210x x x ，，]（1x -x ）（2x -x ）+……+f[n 10x x x ??，]（1x -x ）（2x -x ）…（1-n x -x ）. 其中插值余项为： ） （）！ （） （）（）（）（x 1n f x N -x f x R 1n 1 n n +++==ωξ ξ介于k 10x x x ??，之间。 三、程序设计 function [y,A,C,L]=newdscg(X,Y,x,M) % y 为对应x 的值，A 为差商表，C 为多项式系数，L 为多项式 % X 为给定节点，Y 为节点值，x 为待求节点 n=length(X); m=length(x); % n 为X 的长度 for t=1:m MATLAB 程序设计期中考查 在许多问题中，通常根据实验、观测或经验得到的函数表或离散点上的信息，去研究分析函数的有关特性。其中插值法是一种最基本的方法，以下给出最基本的插值问题——三次样条插值的基本提法： 对插值区间[]b a ,进行划分：b x x x a n ≤ 插值法数值上机实验报告 实验题目： 利用下列条件做插值逼近，并与R (x) 的图像比较 考虑函数：R x y=1 1+x2 （1）用等距节点X i=?5+i,i=0,1,...,10.给出它的10次Newton插值多项式的图像； π),i=0,1,...,20.给出它的20次Lagrange插值多项式（2）用节点X i=5cos(2i+1 42 的图像； （3）用等距节点X i=?5+i,i=0,1,...,10.给出它的分段线性插值函数的图像；（4）用等距节点X i=?5+i,i=0,1,...,10.给出它的三次自然样条插值函数的图像； （5）用等距节点X i=?5+i,i=0,1,...,10.给出它的分段三次Hermite插值函数的图像； 实验图像结果： 实验结果分析： 1.为了验证Range现象，我还特意做了10次牛顿插值多项式和20次牛顿插值多项式的对比图像，结果如下图（图对称，只截取一半） 可以看出，Range现象在高次时变得更加明显。这也是由于高次多项式在端点处的最值随次数的变大很明显。可以料定高次多项式在两侧端点处剧烈震荡，在更小的间距内急剧上升然后下降，Range现象非常明显。 2.分析实验（2）的结果，我们会惊讶地发现，由于取21个点逼近，原本预料的Range现象会很明显，但这里却和f(x)拟合的很好。（即下图中Lagrange p(x)的图像）。可是上图中取均匀节点的20次牛顿多项式逼近的效果在端点处却很差。料想是由于节点X i=5cos2i+1 42 π ,i=0,1,...,20 取得很好。由书上第五章的 知识，对于函数y=1 1+x ，y 1 2对应的cherbyshev多项式的根恰好为X i= 5cos2i+1 42 π ,i=0,1,...,20 。由于所学限制，未能深入分析。 （3）比较三次样条插值图像和Hermit插值图像对原函数图像的逼近情形。见下图： MATLAB程序设计期中作业 ——编程实现牛顿插值 成员:刘川（P091712797）签名＿＿＿＿＿ 汤意（P091712817）签名＿＿＿＿＿ 王功贺（P091712799）签名＿＿＿＿＿ 班级:2009信息与计算科学 学院:数学与计算机科学学院 日期:2012年05月02日 牛顿插值的算法描述及程序实现 一：问题说明 在我们的实际应用中，通常需要解决这样的问题，通过一些已知的点及其对应的值，去估算另外一些点的值，这些数据之间近似服从一定的规律，于是，这就引入了插值法的思想。 插值法是利用函数f (x)在某区间中若干点的函数值，作出适当的特定函数，在这些点上取已知值，在区间的其他点上用这特定函数的值作为函数f (x)的近似值。如果这特定函数是多项式，就称它为插值多项式。利用插值基函数很容易得到拉格朗日插值多项式，公式结构紧凑，在理论分析中甚为方便，但当插值节点增减时全部插值基函数均要随之变化，整个公式也将发生变化，这在实际计算中是很不方便的，为了克服这一缺点，提出了牛顿插值。 二：算法分析 newton 插值多项式的表达式如下： 010011()()()()()n n n N x c c x x c x x x x x x -=+-+???+--???- 其中每一项的系数c i 的表达式如下： 12011010 [,,,][,,,] [,,,]i i i i i f x x x f x x x c f x x x x x -???-???=???= - 即为f (x)在点01,,,i x x x ???处的i 阶差商，（[]()i i f x f x =，1,2,,i n = ），由差商01[,,,]i f x x x ???的性质可知： () 010 1 [,,,]()i i i j j k j k k j f x x x f x x x ==≠???=-∑∏ 牛顿插值的程序实现方法： 第一步：计算[][][][]001012012,,,,,,,n f x f x x f x x x f x x x x 、、、 、。 第二步：计算牛顿插值多项式中01[,,,]i f x x x ???011()()()i x x x x x x ---???-，1,2,,i n = ，得到n 个多项式。三次样条插值代码
牛顿插值法原理及应用
matlab插值法实例
三次样条插值函数
实验名称 插值法
matlab实现数值分析报告插值及积分
插值法和拟合实验报告(数值计算)
拉格朗日插值matlab程序
数值分析作业-三次样条插值
MATLAB数值实验一(数据的插值运算及其应用完整版)
matlab插值法,迭代法程序
数值分析拉格朗日插值法上机实验报告
matlab插值(详细 全面)
插值法实验报告
matlab牛顿插值法例题与程序
(精选)三次样条插值的MATLAB实现
插值法数值上机实验报告
牛顿插值MATLAB算法
相关主题
文本预览