第四章问题一一、问题综述在离地球表面高度为y处的重力加速度如下:计算高度y=55000m处的重力加速度值。
二、问题分析以高度y作为自变量,重力加速度的值为因变量。
得到以下信息:f(0)=9.8100;f(30000)=9.7487;f(60000)=9.6879;f(90000)=9.6278;f(120000)=9.5682;本题要求的就是f(55000)的值。
以下将采用课堂中学到的Lagrange插值多项式法、Newton插值多项式法、分段低次插值法和样条插值法求解该问题。
三、问题解决1. lagrange插值多项式法对某个多项式函数,已知有给定的k+ 1个取值点:其中对应着自变量的位置,而对应着函数在这个位置的取值。
假设任意两个不同的x j都互不相同,那么应用拉格朗日插值公式所得到的拉格朗日插值多项式为:其中每个为拉格朗日基本多项式(或称插值基函数),其表达式为:拉格朗日基本多项式的特点是在上取值为1,在其它的点上取值为0。
源程序lagrange.mfunction [c,f]=lagrange(x,y,a)% 输入:x是自变量的矩阵;y是因变量的矩阵;a是要计算的值的自变量;% 输出:c是插值多项式系数矩阵;f是所求自变量对应的因变量;m=length(x);l=zeros(m,m); % l是权矩阵f=0;for i=1:mv=1;for j=1:mif i~=jv=conv(v,poly(x(j)))/(x(i)-x(j)); % v是l_i(x)的系数矩阵endendl(i,:)=v; % l矩阵的每一行都是x从高次到低次的系数矩阵endc=vpa(y*l,10); % 对应阶次的系数相加,乘以y,显示10位有效数字for k=1:mf=f+c(k)*a^(m-k);end输入矩阵x=[0 30000 60000 90000 120000]y=[9.81 9.7487 9.6879 9.6278 9.5682]a=55000再运行源函数,可得:c =[ -2.057613169e-23, 4.938271605e-18, -3.703703702e-14, -0.000002046111111, 9.81]f =9.6979851723251649906109417384537即此时得到的插值函数为:f(x)=−2.06×10−23x4+4.94×10−18x3−3.70×10−14x2−2.05×10−6x+9.81所以可以得到本题的结果为:f(55000)=9.69802. Newton插值多项式法Newton插值多项式法通过各阶差商来得到插值函数。
牛顿基本插值多项式为:N n(x)=f(x0)+f[x0,x1](x−x0)+⋯+f[x0,…,x1](x−x0)…(x−x n−1)其中零阶差商 f[x i]=f(x i)一阶差商 f[x i,x j]=f(x i)−f(x j) x i−x j二阶差商 f[x i,x j,x k]=f[x i,x j]−f[x j,x k]x i−x k……k阶差商 f[x i,x i+1,…,x i+k]=f[x i+1,…,x i+k]−f[x i,…,x i+k−1]x i+k−x i源程序newton.m:function [c,d,f]=newton(x,y,a)% 输入:x是自变量矩阵;y是因变量矩阵;a是要计算的值的自变量;% 输出:c是插值多项式系数矩阵;d是差商矩阵;f是所求自变量对应的因变量;n=length(x);d=zeros(n,n);d(:,1)=y';f=0;for j=2:nfor i=j:nd(i,j)=(d(i,j-1)-d(i-1,j-1))/(x(i)-x(i-j+1)); % 构造差商矩阵 endendc=d(n,n);for k=(n-1):-1:1c=conv(c,poly(x(k)));m=length(c);c(m)=c(m)+d(k,k); % 又插值公式构造各阶次对应系数矩阵endc=vpa(c,10); % 显示10位有效数字for k=1:mf=f+c(k)*a^(m-k);endend输入矩阵x=[0 30000 60000 90000 120000]y=[9.81 9.7487 9.6879 9.6278 9.5682]a=55000再运行源函数,可得:c =[ -2.057613169e-23, 4.938271605e-18, -3.703703703e-14, -0.000002046111111, 9.81]d =9.8100 0 0 0 09.7487 -0.0000 0 0 09.6879 -0.0000 0.0000 0 09.6278 -0.0000 0.0000 0.0000 09.5682 -0.0000 0.0000 -0.0000 -0.0000f =9.6979851723251030491735014653435即此时得到的插值函数为:f(x)=−2.06×10−23x4+4.94×10−18x3−3.70×10−14x2−2.05×10−6x+9.81本题结果为:f(55000)=9.69803.分段低次插值法当插值点比较多的时候,拉格朗日插值多项式的次数可能会很高,因此具有数值不稳定的特点,也就是说尽管在已知的几个点取到给定的数值,但在附近却会和“实际上”的值之间有很大的偏差。
这类现象也被称为Runge现象,解决的办法是分段用较低次数的插值多项式。
Runge现象源程序seplow.m:function [c,f1,f2]=seplow(x,y)% 输入:x是自变量矩阵;y是因变量矩阵;% 输出:c的每一行是相邻三个点插值多项式系数矩阵;f1、f2分部是所求值落在两段内的插值函数值;m=length(x);c=zeros(m-2,3);for k=1:m-2 % 使用拉格朗日二次插值法对每一小段插值l=zeros(m,3); % l是权矩阵for i=1:3v=1;for j=1:3if i~=jv=conv(v,poly(x(k+j-1)))/(x(k+i-1)-x(k+j-1)); % v是l_i(x)的系数矩阵endendl(i,:)=v; % l矩阵的每一行都是x从高次到低次的系数矩阵endc(k,:)=y*l; % 对应阶次的系数相加,乘以yendc=vpa(c,10);f1=0;for i=1:3f1=f1+c(1,i)*55000^(3-i);endf2=0;for i=1:3f2=f2+c(2,i)*55000^(3-i);endend输入矩阵x=[0 30000 60000 90000 120000]y=[9.81 9.7487 9.6879 9.6278 9.5682]c =[ 2.777777778e-13, -0.000002051666667, 9.81][ 2.777777778e-13, -0.000002068333333, 9.8718][ 2.777777778e-13, -0.000002085, 9.9341]f1 =9.6979986111111076173436997590726f2 =9.7588819444444404252968307863624本题要求的f(55000)中自变量分别落在第一段和第二段中,将其分别带入其中可以得到:f(55000)=9.6980和f(55000)=9.7589取平均值可得:f(55000)=9.7284。
4.样条插值法使用多项式插值,对给定数据集进行插值的n阶多项式就被给定数据点所唯一地定义出来。
但是,对同样的数据进行插值的n阶样条并不是唯一的,为了构建一个唯一的样条插值式它还必须满足另外n-1个自由度。
三次样条插值法(简称Spline插值)是通过一系列基点的光滑曲线,S(x)在区间上存在二阶的连续导数,每个子区间上S(x)都是一个不高于3次的多项式。
通过内部节点函数值相等,内部节点一阶和二阶导数相等,再加上边界条件,即可得到三次样条插值函数。
Matlab中有三次样条插值的函数spline。
本解法直接利用该函数。
源程序ytinterpolation.m:function f=ytinterpolation (x,y,a)% 输入:x是自变量矩阵;y是因变量矩阵;a是要计算的值的自变量;% 输出:f是所求自变量对应的因变量;f=spline(x,y,a); % 求解本题xx = 0:1000:120000;yy = spline(x,y,xx);plot(x,y,'o',xx,yy); % 画三次样条图end输入矩阵x=[0 30000 60000 90000 120000]y=[9.81 9.7487 9.6879 9.6278 9.5682]f =9.6980三次样条法插值图即由三次样条插值法求得:f(55000)=9.6980。
四、总结1.本题分别采用了Lagrange插值多项式法、Newton插值多项式法、分段低次插值法和样条插值法四种方法。
从每一种解法的结果中可以看到Lagrange插值多项式法和Newton插值多项式法求得的结果一模一样,与三次样条插值法求得的结果也相同,而分段低次插值法其中一段的结果与之相同,而另外一组结果则有些不一致,平均计算之后,误差减小了一些。
总体而言,对于本题这种较普通,没有特别极限或者是苛刻条件的问题,上述四种方法均为可行方法。
2.相比于Newton插值法,Lagrange插值法的优点是不需要保存差商,编程简单;缺点是无承袭性,即增加一个节点,所有的基函数都要重新计算。
而Newton插值法则能通过有限差商的方式很好的解决这个问题,大大减小了计算量。
grange插值法和Newton插值法都存在一个问题,就是如果适当提高插值多项式的次数,可能提高计算结果的准确程度。
然而次数越高意味着参加插值节点数越多,这样就会导致计算量增大;而且还可能会出现Runge现象,即插值函数曲线在部分区间上(两端)发生激烈振荡。
所以高次插值效果并不一定比低次插值好。
解决此问题的方式就是采用后两种方法即分段低次插值或样条插值。
4.分段低次插值的优点就是局部性质好,能保证精度要求,而且改变一个节点只会对附近两个节点产生影响。
缺点就是节点处导数可能不连续,从而导致曲线不光滑。
而三次样条插值法分别考虑了一阶导数和二阶导数,所以形成的曲线是光滑的,一般来说可以满足工程上的需求。
问题二一、问题综述停车过程中需要行驶的距离包括大脑反应过程中行驶的距离和刹车过程中行驶的距离,两者都是速度的函数。