当前位置:文档之家› 数值分析上机作业1-1

数值分析上机作业1-1

数值计算方法上机题目11、实验1. 病态问题实验目的:算法有“优”与“劣”之分,问题也有“好”和“坏”之别。

所谓坏问题就是问题本身的解对数据变化的比较敏感,反之属于好问题。

希望读者通过本实验对此有一个初步的体会。

数值分析的大部分研究课题中,如线性代数方程组、矩阵特征值问题、非线性方程及方程组等都存在病态的问题。

病态问题要通过研究和构造特殊的算法来解决,当然一般要付出一些代价(如耗用更多的机器时间、占用更多的存储空间等)。

问题提出:考虑一个高次的代数多项式∏=-=---=201)()20)...(2)(1()(k k x x x x x p (E1-1)显然该多项式的全部根为l ,2,…,20,共计20个,且每个根都是单重的(也称为简单的)。

现考虑该多项式方程的一个扰动0)(19=+xx p ε (E1-2)其中ε是一个非常小的数。

这相当于是对(E1-1)中19x 的系数作一个小的扰动。

我们希望比较(E1-1)和(E1-2)根的差别,从而分析方程(E1-1)的解对扰动的敏感性。

实验内容:为了实现方便,我们先介绍两个 Matlab 函数:“roots ”和“poly ”,输入函数u =roots (a )其中若变量a 存储1+n 维的向量,则该函数的输出u 为一个n 维的向量。

设a 的元素依次为121,...,,+n a a a ,则输出u 的各分量是多项式方程0...1121=++++-n n n n a x a x a x a的全部根,而函数b=poly(v)的输出b 是一个n +1维变量,它是以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)的解应当相差很小。

计算中你有什么出乎意料的发现?表明有些解关于如此的扰动敏感性如何?(2)将方程(E1-2)中的扰动项改成18x ε或其他形式,实验中又有怎样的现象出现?实验步骤:(1)程序function t_charpt1_1clcresult=inputdlg({'请输入扰动项:在[0 20]之间的整数:'},'charpt1_1',1,{'19'});Numb=str2num(char(result));if((Numb>20)|(Numb<0))errordlg('请输入正确的扰动项:[0 20]之间的整数!');return;endresult=inputdlg({'请输入(0 1)之间的扰动常数:'},'charpt1_1',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.56717i 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.8696+0i 17.2186+0i15.4988+0.0211828i 15.4988-0.0211828i 13.7707+0i 13.1598+0i11.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的图像如下:从实验的图形中可以看出,当ess 充分小时,方程E.1.1和方程E.1.2的解相差很小,当ess 逐渐增大时,方程的解就出现了病态解,这些解都呈现复共轭性质。

(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时没有出现复共轭,误差很小。

因此,扰动作用到x n 上时,n 越小,扰动引起的误差越小。

2、实验2。

多项式插值的振荡现象,即插值的龙格(Runge )现象问题提出:考虑在一个固定的区间上用插值逼近一个函数。

显然,拉格朗日插值中使用的节点越多,插值多项式的次数就越高、自然关心插值多项式的次数增加时,)(x L n 是否也更加靠近被逼近的函数。

龙格给出的一个例子是极著名并富有启发性的。

设区间]1,1[-上函数22511)(xx f +=实验内容:考虑区间]1,1[-的一个等距划分,分点为n i nix i ,...,2,1,0,21=+-= 则拉格朗日插值多项式为∑=+=ni iin x l x x L 02)(2511)( 其中的)(x l i ,n i ,...,2,1,0=是n 次拉格朗日插值基函数。

实验要求:(l )选择不断增大的分点数目,...3,2=n ,画出原函数)(x f 及插值多项式函数)(x L n 在]1,1[-上的图像,比较并分析实验结果。

(2)选择其他的函数,例如定义在区间[-5,5]上的函数x x g x xx h arctan )(,1)(4=+= 重复上述的实验看其结果如何。

(3)区间],[b a 上切比雪夫点的定义为1,...,2,1,)1(2()12(cos 22+=⎪⎪⎭⎫ ⎝⎛+--++=n k n k ab a b x k π 以121,...,,+n x x x 为插值节点构造上述各函数的拉格朗日插值多项式,比较其结果。

实验步骤:(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 k=1:n p=1.0; for j=1:n if (j ~= k)p = p*(z - x0(j))/(x0(k) - x0(j)); end ends = s + p*y0(k); end y(i) = s; endfunction t_charpt2promps = {'请选择实验函数,若选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' & Nb_f ~= 'h' & Nb_f ~= 'g')errordlg('实验函数选择错误!');return ;endresult = inputdlg({'请输入插值结点数N:'},'charpt_2',1,{'10'}); Nd = str2num(char(result));if (Nd <1)errordlg('结点输入错误!');return ;endswitch Nb_fcase'f'f=inline('1./(1+25*x.^2)'); a = -1;b = 1;case'h'f=inline('x./(1+x.^4)'); a = -5; b = 5;case'g'f=inline('atan(x)'); a = -5; b= 5;endx0 = linspace(a, b, Nd+1); y0 = feval(f, x0);x = 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=6n=7 n=8从图中可以看出,随着n的增大,拉格朗日插值函数在x=0附近较好地逼近了原来的函数f(x),但是却在两端x= -1和x=1处出现了很大的振荡现象。

通过分析图形,可以看出,当n为奇数时,虽然有振荡,但振荡的幅度不算太大,n为偶数时,其振荡幅度变得很大。

(2)将原来的f(x)换为其他函数如h(x)、g(x),结果如图所示。

其中h(x), g(x)均定义在[-5,5]区间上,h(x)=x/(1+x4),g(x)=arctan x。

h(x), n=7 h(x), n=8h(x), n=9 h(x), n=10g(x), n=8 g(x), n=9g(x), n=12 g(x), n=13分析两个函数的插值图形,可以看出:随着n的增大,拉格朗日插值函数在x=0附近较好地逼近了原来的函数f(x),但是却在两端x= -5和x=5处出现了很大的振荡现象。

相关主题