计算实习报告

1、计算实习报告题目1 已知函数，试用以下数据通过构造插值（1）求出的近似值并估计插值误差；（2）求出的近似值并画出图形；0.450.550.650.750.85-0.798508-0.597837-0.430783-0.287682-0.1625191.1基本原理：插值形式是，其中li(x)是对应于每个插值节点的基函数，ci为系数。对于给出n+1个插值节点的代数插值问题，称满足下列条件的函数li(x) (i=0,1,2,、，n) 为lagrange型的节点基函数。是n次多项式 li(xj)=ij (i,j=0,1,2,、，n)ij=1 （j=i）ij=0 (ji)即lagrange节点基函数是在

2、该节点上取值为1，而在其余节点取值为0的多项式。 构造方法：由，设，a待定由，li(xi)=1可求出a,故得节点基函数为 （i=0,1,2,n）由插值条件ln(xi)=yi直接得ci=yi (i=0,1,2,n)综上得lagrange插值函数是1.2 matlab源程序多项式插值%文件名:laginterp.mfunction yy=laginterp(x,y,xx)%功能:lagrange多项式插值%输入:x,y-原始数据（要求x,y同长度且x(i)互异）% xx-欲插值点数据（可为向量）% 输出:yy-返回插值点处的插值结果n=length(x); yy=zeros(size(xx);fo

3、r j=1:n t=ones(size(xx);for i=1:n; if i=j; t=t.*(xx-x(i)/(x(j)-x(i);end; endyy=yy+t*y(j);end 误差%文件名wuch1function wuch1n=5;xx=0.1:0.1:1.8;y0=log(xx);x=0.45 0.55 0.65 0.75 0.85;y=-0.798508 -0.597837 -0.430783 -0.287682 -0.162519; yy=laginterp(x,y,xx);er1=(y0-yy)*100;z=0*xx;plot(xx,z,k,xx,er1,b,xx,y0,r

4、,xx,yy,g);axis(0 2 -2 1);gtext(e(x)*100),gtext(ln(x),gtext(laginterp(x,y,xx) 1.3运行结果与误差（1）yy = -0.5342 = -0.5108 绝对误差e=0.0234 相对误差（2）的近似值xx0.10.20.30.40.50.60.70.80.9yy-1.9907-1.5436-1.1927-0.9155-0.6932-0.5108-0.3567-0.2231-0.1058er*100-31.1862-6.5832-1.1306-0.08230.0077-0.00280.0025-0.00520.0422xx

5、1.01.11.21.31.41.51.61.71.8yy-0.00420.07850.13550.15660.12750.0307-0.1553-0.4554-0.8983er*1000.42301.68544.679910.580320.895437.477462.527798.6006148.6057差值图形：误差图形：总结：由以上插值数据及图形可以看到在插值节点范围内及其附近（0.1-1.2），插值值与精确值较逼近，在较远的地方误差叫大,说明插值方法不稳定。1.4结果分析 lagrange插值函数的优点是容易建立，插值节点基函数简明。但当增加节点时原有插值函数不可利用，需重新计算。且形

6、式不易简化，计算量大，不易编程运算。实际应用时，并不完全使用全部给出的节点。因为高次插值计算量大，且高次插值精度不一定高。对某些被插函数用lagrange高次插值多项式还可能降低精度，甚至出现不收敛的情况。故常用低次线性插值和分段低次插值。题目2 利用函数在节点处的值，用newton插值求点处的函数值，并作图比较。2.1 基本原理插值函数形式是 ，是newton插值基函数，ci是系数。取 （i=1,2,n）则 nn(x)=c0+c1(x-x0)+c2(x-x0)(x-x1)+cn(x-x0)(x-x1)(x-xn-1)根据插值条件依次求出ci一般形式有是f(x)在点， 上的n阶差商。故newt

7、on插值函数用差商作为系数，一般形式是2.2 matlab源程序function yy,c,p=newtoninterp(x,y,xx)%功能:newton插值%输入:x,y-插值原始数据，要求x按升序排列,且x,y长度一致% xx-欲插值点,可为标量或向量%输出: yy-插值结果;c(可选)-newton插值系数;%p(可选)-插值多项式系数(基底为1,x,.,x(n-1)，按x的幂降序排列%step1：变量初始化n=length(x); q = zeros(n,n);%差商表p=zeros(1,n);%插值多项式系数(基底为1,x,.,x(n-1)，按x的幂降序排列 q(:,1)=y(:)

8、;%第一列为函数值% step2：构造差商表for i = 1:n-1 for j = 1:i;q(i+1,j+1) = (q(i+1,j) - q(i,j) / (x(i+1) - x(i-j+1);end end% step 3:利用九昭算法计算牛顿插值结果yy=q(n,n);for i=n-1:-1:1; yy=yy.*(xx-x(i)+q(i,i);end% step 4:可选输出if nargout1 c=q; if nargout2 s=zeros(1,n);p=zeros(1,n); for i=1:n t=1; for j=1:i-1;t=conv(t,1,-x(j);end

9、s(1,n-i+1:n)=t*q(i,i); p=p+s; end endend2.3运行结果点处函数值分别为： 91.4050 21.2417 12.7889 -6.1932作图比较2.4结果分析newton插值函数每增加一个新节点不需全部重新计算，只需计算新增加的一项。故newton插值便于计算和编程运行。题目3 对函数分别用等距节点和点插值并作图比较误差。3.1基本原理若给出函数y=f(x)在等距节点xi=x0+ih(i=0,1,2,n)上的值，由差商与差分关系得差分表示的是newton插值函数。对于chebyshev多项式，当取，（k=0,1,n）。tn（x）在-1,1轮流取最大值1与

10、最小值-1. （k=0,1,n）3.2 matlab源程序构造n阶chebyshev正交多项式function pn=cheby_pw(n) pbb=1;if n=0,pn=pbb; endpb=1 0;if n=1,pn=pb;endfor i=2:n pn=2*pb,0-0,0,pbb; pbb=pb;pb=pn; end 3.3 运行结果与误差（newton等距节点插值）y110.03850.10730.15700.19040.21010.21850.21780.20980.19640.17910.15930.13830.11710.09670.0777等距节点插值误差e100.0673

11、0.11540.14710.16500.17150.16860.15850.14270.12300.10050.07660.05230.0286.0.0061节点x00-0.9808 -0.8315-0.5556-0.1951 0.19510.5556 0.8315 0.9808(节点插值)y220.02340.04050.05210.05920.06270.06350.06230.05970.05630.05260.04910.04620.04420.04330.0439节点插值误差e20.01510.00050.01050.01590.01760.01650.01320.00830.002

12、6 0.00360.00970.01550.02060.02480.0278（精确）y0.03850.04000.04160.04330.04510.04710.04910.05130.05360.05610.05880.06170.06480.06810.07163.4结果分析 当插值区间限制在上时，若取为点作为插值节点，即的零点作为我们的插值点列，此时得到相应的插值误差界最小，也称极小化插值。当插值区间为一般区间上时，可以先做变量代换化为。题目4对函数分别利用分段线性和分段三次插值并作图比较误差。4.1 基本原理1、分段线性插值对于给出的n+1个节点可以用其中每两个点构造分段线性插值。分段

13、线性插值函数满足：整体上连续，局部上是线性函数，即在xi,xi+1上是线性函数，i=0,1,2n-1满足插值条件 用局部非零的分段线性插值基函数来构造，其中： 则分段线性插值函数2、分段三次插值利用基函数插值的方法可写出分段三次插值多项式的一般形式如下其中的形式如下： 当时可写为： 4.2 matlab源程序1. 分段线性插值% 文件名：piecelin.mfunction yy = piecelin(x,y,xx)% 功能：分段线性插值% 输入: x,y -原始数据对（要求xi互异且升序排列;xx -欲插值点% 输出: yy-插值结果delta = diff(y)./diff(x);%作插商

14、i = binsearch(x,xx); % 二分法寻找节点所在区间s = xx - x(i); yy = y(i) + s.*delta(i); % 计算插值结果end上述程序中用到二分法搜索插值点所在区间的子函数binsearch.m如下：%文件名:binsearch.mfunction ia = binsearch(x,xx)% 功能：二分法搜索插值点所在区间% 输入: x-数据向量（要求升序排列）% xx-欲搜索点数据（可为向量或标量，要求位于x所在区间）%输出: ia-使得x(i)= xx 1 im = fix(ia(k)+ib(k)/2); % 取中间值 if x(im) x=-2:0.1:2; y=1./(1+25*x.2); xx=-2:0.03:2;yy=piecelin(x,y,xx);plot(xx,yy,r);hold on;ezplot(1/(1+25*x2)clc; x=-2:0.1:2; y=1./(1+25*x.2); xx=-2:0.03:2; f1=diff(1/(1+25*x2);f1 =-(50*x)/(25*x