《数值分析上机作业1-1解析.pdf》由会员分享,可在线阅读,更多相关《数值分析上机作业1-1解析.pdf(9页珍藏版)》请在三一文库上搜索。
1、1 数值计算方法上机题目1 1、实验 1. 病态问题 实验目的 : 算法有“优”与“劣”之分,问题也有“好”和“坏”之别。所谓坏问题就是问题本身 的解对数据变化的比较敏感,反之属于好问题。 希望读者通过本实验对此有一个初步的体会。 数值分析的大部分研究课题中,如线性代数方程组、矩阵特征值问题、非线性方程及方 程组等都存在病态的问题。病态问题要通过研究和构造特殊的算法来解决,当然一般要付出 一些代价(如耗用更多的机器时间、占用更多的存储空间等)。 问题提出 : 考虑一个高次的代数多项式 20 1 )()20).(2)(1()( k kxxxxxp(E1-1) 显然该多项式的全部根为l,2, 20
2、,共计 20 个,且每个根都是单重的(也称为简 单的)。现考虑该多项式方程的一个扰动 0)( 19 xxp(E1-2) 其中是一个非常小的数。这相当于是对(E1-1)中 19 x的系数作一个小的扰动。我们希望 比较( E1-1)和( E1-2)根的差别,从而分析方程(E1-1)的解对扰动的敏感性。 实验内容 : 为了实现方便,我们先介绍两个Matlab 函数:“ roots”和“ poly”,输入函数 uroots(a) 其中若变量a存储1n维的向量,则该函数的输出u为一个n维的向量。设a 的元素依次 为 121 ,., n aaa,则输出u 的各分量是多项式方程 0. 1 1 21nn nn
3、 axaxaxa 的全部根,而函数 b=poly(v) 的输出 b 是一个 n1 维变量,它是以 n 维变量 v 的各分量为根的多项式的系数。可见 “roots” 和“ Poly”是两个互逆的运算函数. ve=zeros(1,21); ve(2)=ess; roots(poly(1:20)+ve) 上述简单的Matlab 程序便得到 (E1-2)的全部根, 程序中的 “ess ”即是( E1-2)中的。 实验要求 : (1)选择充分小的ess,反复进行上述实验,记录结果的变化并分析它们。如果扰动项 的系数很小,我们自然感觉(E1-1)和 (E1-2)的解应当相差很小。计算中你有什么出乎意 料的
4、发现?表明有些解关于如此的扰动敏感性如何? (2)将方程( E1-2)中的扰动项改成 18 x或其他形式,实验中又有怎样的现象出现? 2 实验步骤: (1)程序 function t_charpt1_1 clc result=inputdlg( 请输入扰动项: 在0 20之间的整数 :, charpt 1_1,1,19); Numb=str2num(char(result); if(Numb20)|(Numb0)errordlg( 请输入正确的扰动项:0 20之间的整 数!); return; end result=inputdlg( 请输入 (0 1)之间的扰动常数:,charpt 1_1,
5、1,0.00001); ess=str2num(char(result); ve=zeros(1,21); ve(21-Numb)=ess; root=roots(poly(1:20)+ve); x0=real(root); y0=imag(root); plot(x0,y0, *); disp( 对扰动项 ,num2str(Numb), 加扰动 ,num2str(ess), 得到的全部根 为:); disp(num2str(root); 二、实验结果分析 ess分别为 1e-6,1e-8.1e-10,1e-12. 对扰动项19 加扰动 1e-006 得到的全部根为: 21.3025+1.56
6、717i 21.3025-1.56717i 18.5028+3.6004i 18.5028-3.6004i 15.1651+3.76125i 15.1651-3.76125i 12.4866+2.88278i 12.4866-2.88278i 10.5225+1.71959i 10.5225-1.71959i 9.04485+0.594589i 9.04485-0.594589i 7.9489+0i 7.00247+0i 5.99995+0i 5+0i 4+0i 3+0i 2+0i 1+0i 对扰动项19 加扰动 1e-010 得到的全部根为: 19.9953+0i 19.0323+0i 17
7、.8696+0i 17.2186+0i 15.4988+0.0211828i 15.4988-0.0211828i 13.7707+0i 13.1598+0i 11.9343+0i 11.029+0i 9.99073+0i 9.00247+0i 7.99952+0i 7.00007+0i 5.99999+0i 5+0i 4+0i 3+0i 2+0i 1+0i ess分别为 1e-6,1e-8.1e-10,1e-12 的图像如下: 3 从实验的图形中可以看出,当 ess充分小时, 方程 E.1.1 和方程 E.1.2 的解相差很小, 当 ess逐渐增大时,方程的解就出现了病态解,这些解都呈现复共
8、轭性质。 (2) 将扰动项加到x 18 上后, ess=1e-009 时方程的解都比较准确,没有出现复共轭现象。 ess=1e-008时误差与x 19(ess=1e-009)时相当,即扰动加到 x 18 上比加到x 19 小一个数量级。对 x 8 的扰动 ess=1000时没有出现复共轭, 误差很小;对 x 的扰动 ess=10e10时没有出现复共轭, 误差很小。因此,扰动作用到xn上时, n 越小,扰动引起的误差越小。 2、实验 2。多项式插值的振荡现象,即插值的龙格(Runge)现象 问题提出 : 考虑在一个固定的区间上用插值逼近一个函数。显然,拉格朗日插值中使用的节点越多, 插值多项式的
9、次数就越高、自然关心插值多项式的次数增加时,)(xLn是否也更加靠近被逼 近的函数。龙格给出的一个例子是极著名并富有启发性的。设区间 1 , 1上函数 2 251 1 )( x xf 实验内容 : 考虑区间1 , 1的一个等距划分,分点为 ni n i xi,.,2, 1 ,0, 2 1 则拉格朗日插值多项式为 n i i i n xl x xL 0 2 )( 251 1 )( 其中的)(xli,ni,.,2, 1 ,0是 n 次拉格朗日插值基函数。 实验要求 : (l) 选择不断增大的分点数目,.3 , 2n, 画出原函数)(xf及插值多项式函数)(xLn在 1 , 1上的图像,比较并分析实
10、验结果。 (2)选择其他的函数,例如定义在区间-5,5上的函数 4 xxg x x xharctan)(, 1 )( 4 重复上述的实验看其结果如何。 (3)区间,ba上切比雪夫点的定义为 1,.,2, 1, ) 1(2( )12( cos 22 nk n kabab xk 以 121 ,., n xxx为插值节点构造上述各函数的拉格朗日插值多项式,比较其结果。 实验步骤: (1) 试验程序: function y=Lagrange(x0, y0, x); % Lagrange插值 n= length(x0); m=length(x); for i=1:m z=x(i); s=0.0; for
11、 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 = s + p*y0(k); end y(i) = s; end function t_charpt2 promps = 请选择实验函数,若选f(x),请输入 f,若选 h(x),请输入 h, 若选 g(x),请输 入 g:; titles = charpt_2; result = inputdlg(promps,charpt 2,1,f); Nb_f = char(result); if(Nb_f = f return; end resul
12、t = inputdlg( 请输入插值结点数N:,charpt_2,1,10); Nd = str2num(char(result); if(Nd 1)errordlg( 结点输入错误! ); return; end switch Nb_f 5 casef f=inline(1./(1+25*x.2); a = -1;b = 1; caseh f=inline(x./(1+x.4); a = -5; b = 5; caseg f=inline(atan(x); a = -5; b= 5; end x0 = linspace(a, b, Nd+1); y0 = feval(f, x0); x =
13、 a:0.1:b; y = Lagrange(x0, y0, x); fplot(f, a b, co); hold on ; plot(x, y, b-); xlabel(x); ylabel(y = f(x) o and y = Ln(x)-); 增大分点 n=2,3,时,拉格朗日插值函数曲线如图所示。 n=3 n=6 n=7 n=8 从图中可以看出,随着n 的增大,拉格朗日插值函数在x=0 附近较好地逼近了原来的函 数 f(x) ,但是却在两端x= -1 和 x=1 处出现了很大的振荡现象。通过分析图形,可以看出, 当 n 为奇数时,虽然有振荡,但振荡的幅度不算太大,n 为偶数时,其振荡
14、幅度变得很大。 (2) 将原来的 f(x) 换为其他函数如h(x)、g(x),结果如图所示。 其中 h(x), g(x) 均定义在 -5, 5区间上, h(x)=x/(1+x 4),g(x)=arctan x 。 6 h(x), n=7 h(x), n=8 h(x), n=9 h(x), n=10 g(x), n=8 g(x), n=9 g(x), n=12 g(x), n=13 分析两个函数的插值图形,可以看出: 随着 n 的增大, 拉格朗日插值函数在x=0 附近较好地逼近了原来的函数f(x) ,但是却在 两端 x= -5 和 x=5 处出现了很大的振荡现象。通过图形可以看出,当n 为偶数时
15、,虽然有振 荡,但振荡的幅度不算太大,n 为奇数时,其振荡幅度变得很大。原因和上面f(x) 的插值类 似, h(x)、 g(x)本身是奇函数,如果n 为偶数,那么Lagrange 插值函数Ln(x)的最高次项 x n-1 是奇次幂,比较符合h(x)、g(x)本身是奇函数的性质;如果n 为奇数,那么Lagrange 插值函 数 Ln(x)的最高次项 x n-1 是偶次幂,与h(x)、g(x)本身是奇函数的性质相反,因此振荡可能更 7 剧烈。 3、实验 3。 样条插值的收敛性 问题提出 : 一般的多项式插值不能保证收敛性,即插值的节点多,效果不一定就好。对样条函数插 值又如何呢?理论上证明样条插值
16、的收敛性是比较困难的,也超出了本课程的内容。通过本 实验可以验证这一理论结果。 实验内容 : 请按一定的规则分别选择等距或者非等距的插值节点,并不断增加插值节点的个数。考 虑实验2.中的函数或选择其它你有兴趣的函数,可以用Matab 的函数“sp line” 作此函数的 三次样条插值。在较新版本的Matlab 中,还提供有spline 工具箱,你可以找到极丰富的样 条工具,包括B- 样条。 实验要求 : (1)随节点个数增加,比较被逼近函数和样条插值函数误差变化情况。分析所得结果 并与拉格朗目多项式插值比较。 (2)样条插值的思想最早产生于工业部门。作为工业应用的例子,考虑如下问题:某 汽车制
17、造商用三次样条插值设计车门的曲线,其中一段的数据如下: k x0 1 2 3 4 5 6 7 8 9 10 k y0.0 0.79 1.53 2.19 2.71 3.03 3.27 2.89 3.06 3.19 3.29 k y0.8 0.2 实验步骤: (1)程序: function t_charpt2 promps = 请选择实验函数,若选f(x),请输入 f,若选 h(x),请输入 h, 若选 g(x),请输 入 g:; titles = charpt_2; result = inputdlg(promps,charpt 2,1,f); Nb_f = char(result); if(N
18、b_f = f return; end result = inputdlg( 请输入插值结点数N:,charpt_2,1,10); Nd = str2num(char(result); if(Nd 1)errordlg( 结点输入错误! ); return; end switch Nb_f casef f=inline(1./(1+25*x.2); a = -1;b = 1; caseh f=inline(x./(1+x.4); a = -5; b = 5; caseg 8 f=inline(atan(x); a = -5; b= 5; end x0 = linspace(a, b, Nd+1
19、); y0 = feval(f, x0); x = a:0.1:b; cs = spline(x0, y0); y = ppval(cs, x); plot(x0, y0, o); hold on ; plot(x, y, k-); xlabel(x); ylabel(y = f(x) o and y = Spline(x)-); 实验结果: 如图所示。 f(x), n=5 n=10 n=20 h(x), n=5 h(x), n=10 n=20 g(x), n=5 n=10 n=20 图中可以看出, 由于其采用了分段三次多项式拟合的方法,随着三次样条插值的插值结 点的增加,并没有出现振荡现象。 (2) 程序: x0=0:10; y0=0.0 0.79 1.53 2.19 2.71 3.03 3.27 2.89 3.06 3.19 3.29; x=0:0.1:10; pp=csape(x0,y0,complete,0.8 0.2); y = ppval(pp, x); plot(x0, y0, o); hold on ; plot(x, y, k-); xlabel(x); ylabel(y = f(x) o and y = Spline(x)-); 9 车门的曲线如下图:
链接地址:https://www.31doc.com/p-4954668.html