文档视界 最新最全的文档下载
当前位置:文档视界 › 数值计算实验报告MATLAB线性方程组的求解方法

数值计算实验报告MATLAB线性方程组的求解方法

评分

签名

日期

湖南商学院实验报告

课程名称数值分析

实验名称数值计算基本概念及算法基础

专业班级信科1301

指导老师胡桔州

组长姓名陈平学号 130320015 电话 150******** Q Q 1274839822 成员姓名周红学号 130320021 王双学号 130320022

贺嘉玲学号 130320026

2014—2015学年度第2学期

一、实验目的与要求

1、通过实验,进一步了解和掌握浮点数及其运算特点,了解阶段误差与舍入误差这两种主要误差及其对结果的影响;

2、通过实验,初步理解数值算法的复杂性、算法的稳定性、问题的病态性;

3、通过实验,比较数值解法与纯数学解法的区别,体会算法在工程与科学计算中的地位,并形成有意识比较算法、选择算法的良好习惯;从而进一步提高设计算法的能力。

二、实验环境

(包括硬件、软件配置)

硬件:屏幕尺寸:14英寸1366x768

CPU型号:Intel酷睿i34030U ... (2种)

CPU主频:1.9GHz ... (2种)

内存容量:4GB(4GB×1)DDR3

硬盘容量:500GB

显卡芯片:AMDRadeonR5M230

操作系统:预装Windows764bit

软件:开发环境:Window 7 旗舰版

开发工具:MATLABR2009a(MATLAB7.8)

三、项目内容

(包括相关准备知识)

浮点数的概念,浮点数的运算,截断误差和舍入误差的区别以及它们对结果的影响。以及用MATLAB比较两个不同算法的稳定性,通过研究和构造特定算法解决病态问题。

要做好实验,这些都是必须思考并掌握的技能。

(1)熟悉MATLAB的主要操作命令。

(2)学会简单的算法编程。

(3)掌握简单的程序调适。

(4)用MATLAB编程并利用函数解决问题。

四、实验步骤:

(要求详细对实验数据的处理、源程序代码、算法、实验结果、图表等进行描述,可以根据情况自己添加页)

一、浮点数及其运算特点

1、浮点数的分布

function f=fudianshu(beta,t,L,U)

m=1/beta:1/beta^(t):1/beta+(beta^(t-1)-1)/beta^(t);

n=length(m);

s=L:U;

for k=1:(U-L+1)

format rat

f(k,:)=beta^s(k)*m;

f(k,:)

plot(f,zeros(1,n), '.')

axis([0 4 -1 1])

pause

end

pause

subplot(2,1,2)

hist(f(:))

%其中beta为基数,t为精确度,L和U分别是指数部分的上界和下界。类似于二进制的情形。在命令窗口输入:

>> fudianshu(2,3,-2,2)

ans =

1/8 5/32 3/16 7/32

ans =

1/4 5/16 3/8 7/16

ans =

1/2 5/8 3/4 7/8

ans =

1 5/4 3/

2 7/4

ans =

2 5/2

3 7/2

ans =

1/8 5/32 3/16 7/32

1/4 5/16 3/8 7/16

1/2 5/8 3/4 7/8

1 5/4 3/

2 7/4

2 5/2

3 7/2

00.51 1.52 2.53 3.54

-1

-0.500.5100.51 1.52 2.53 3.5

2468

分析:从图中可以看出浮点数在接近其下界时比较稠密,而接近上界时则相对比较稀疏。我们要思考的是这种分布不均匀性会给数值计算带来问题吗?

假设x 是精确值,p 是对x 的近似,则a=|x-p|;称为p 近似x 的绝对误差。而r=a/|x|=|x-p|/|x|,称为p 近似x 的相对误差。可见用浮点数近似接近上界的实数时,由于这时浮点数的间隔大,所以近似的绝对误差则会比较大;然而只要被近似的数在在下界和上界之间,相对误差在各处是一样的。所以数值计算中更多的还是相对误差的问题,因此浮点数分布的不均匀并不会影响算法噢!

2、双精度浮点数下的问题--舍入误差

>> format long >> x=4/3-1 x =

0.333333333333333 >> y=3*x y =

1.000000000000000 >> z=1-y z =

2.220446049250313e-016

按照理论计算,得到的z 值应当是0.而实际输出的非零。其实得到的z 就是机器精度。如果计算中产生一个小于realmin 的数则是“下溢”;如果产生了一个大于realmax 的数则是“上溢”。我们再来看一个多项式的计算: >> x=0.988:0.001:1.012;

>> y=x.^7-7*x.^6+21*x.^5-35*x.^4+35*x.^3-21*x.^2+7*x-1; >> plot(x,y)

0.985

0.990.9951 1.005 1.01 1.015

-4

-3

-2-1012345x 10

-14

从画出的图来看,这完全不像一个多项式!这个现象就是舍入误差在起作用。

3、浮点数不满足结合律和交换律

结合律

format long

a=rand(1,20000);

b=rand(1,20000);m=2^20*rand(1,20000); s1=sum(a)+sum(m.*b); s2=sum(a+m.*b); s1,s2 sub=s1-s2 输出结果为: s1 =

5.261188112279733e+009 s2 =

5.261188112279694e+009 sub =

3.910064697265625e-005

交换律

m=100000;x=1:m;y=m:-1:1;

s=m/(m+1);disp(['s=',num2str(s)])

s1=sum(1./(x.*(x+1)));s2=sum(1./(y.*(y+1))); fprintf('s1=%31.30f\ns2=%31.30f\n',s1,s2)

fprintf('s1-s=%31.30f\ns2-s=%31.30f\n',s1-s,s2-s) fprintf('s1-s2=%31.30f\n',s1-s2) 输出结果为: s=0.99999

s1=0.999990000100012160000000000000 s2=0.999990000099998940000000000000 s1-s=0.000000000000013100631690576847

分析:实验结果表明,浮点数加法不满足交换律和结合律。这是为什么呢?首先我们要了解计算机浮点数的加法是如何进行的。假设有两个浮点数:

d=+ -.d1*d2*...*dt*beta1,c=+-.c1*c2*...*ct*beta2。两个浮点数相加首先要比较它们的阶码。如果阶码相同则尾数相加,相加后尾数如果大于1,阶码进位;如果阶码不等,以相对大的阶码为标准,将阶码小的浮点数移位直到阶码一致,再按阶码相等时的规则进行相加。如果一个大数加一个非常小的数,小数的个数越多,计算的误差就越大。

二、算法的复杂性

1、测试多项式函数计算的不同算法在运行时间上的差别

①参考教材P35......

>> tic;

>> x=3;

>> f=2*x^5-5*x^4-4*x^3-6*x+7;

>> t=toc

t =

58.987017134798052

>> tic;

>> x=3;

>> x^2==x*x;

>> x^3==x^2*x;

>> x^4==x^3*x;

>> x^5==x^4*x;

>> f=2*x^5-5*x^4-4*x^3-6*x+7;

>> t=toc

t =

2.152971514854767e+002

>> tic;

>> x=3;

>> f=((((2*x-5)*x-4)*x+3)*x-6)*x+7;

>> t=toc

t =

10.009796315604856

>> x=3;

>> tic;

>> x=3;

>> f=((((2*x-5)*x-4)*x+3)*x-6)*x+7;

>> t=toc

t =

10.840845190320568

②>> tic;

>> A=rand(20,20);

>> t=toc

t =

31.961928037591385

>> tic;

>> A=rand(20,20);

>> b=det(A);

>> t=toc

t =

42.956819777617717

2、误差传播与算法的稳定性

算法一:

用递推公式E1=1/e。E n=1-n*E n-1。n=2,3, (20)

clear

ep(1)=exp(-1);

for n=2:9

ep(n)=1-n*ep(n-1);

ep(n)=round(10000*ep(n))/10000;

end

ep'

plot(ep,'b*');

xlabel('算法一');

%axis([0 100 -10 10])

1234

56789

-1

012345678算法一

假设E n =I n -Ln,则有E n =-n*E n-1=(-1)^n!E n 。因此误差是逐渐放大的。

算法二、递推公式E 20=0,E n-1=(1-E n-1)/n,n=19,18, (1)

% 误差传播与算法稳定性 % 稳定的算法 clear ep(500)=0; for n=500:-1:2

ep(n-1)=(1-ep(n))/n;

ep(n-1)=round(10000*ep(n-1))/10000; end ep

plot(ep,'b-'); axis([0 500 -0.10 1]) ep(1:9)'

xlabel('算法二');

%axis([0 100 -10 10])

50

100

150

200

250300

350

400

450

500

-0.1

00.10.20.30.40.50.60.70.80.9算法二

假设E n =I n -Pn ,则有|En-1|=1/n*|En|,E k =(-1)^(n-k)*(k!/n!)*E n 。因此误差是逐渐缩小的。

分析:由此我们可以判断算法二比算法一计算稳定。

3、病态问题的初步讨论

% 病态问题实验

% 代数多项式的根及其扰动后多项式的根 p=poly(1:20); ep=zeros(1,21); ep(3)=1.0e-5; re=roots(p+ep) plot(re,'b+'); hold on

plot(1:20,0,'r*'); hold off

0510152025

-4

-3-2-101234

分析:从图中我们可以看出,虽然扰动项的系数非常小,从个人感觉上说对多项式的解是不会有太大的影响。但通过计算我们发现,有很多解却差异极大。即有很多解对这种扰动行为极其敏感。

这个实验告诉我们,不但算法有好坏之分,被求解的问题也有好坏之分,甚至在同一问题中有好的部分和坏的部分。数值分析的大部分领域中,如线性代数方程组、矩阵特征值问题等等,都有关于病态问题的深入研究。病态问题在付出一些代价(如耗用更多机器时间,占用更多存储空间等)后,能通过研究和构造特殊的算法来解决。

三、思考题

1、函数sinx 有幂级数展开,利用幂级数计算sinx 的matlab 程序为: function s=powersin(x) s=0; t=x; n=1;

while s+t~=s; s=s+t;

t=-x^2/((n+1)*(n+2))*t;

n=n+2;

End

(a)解释上述程序的终止准则

(b)对于x=pi/2,11pi/2,21pi/2,计算的精度是多少?分别需要计算多少项?答:

(a)终止准则为终止条件为t=0,但是t=-x^2/((n+1)*(n+2))*t始终不为0。不过在计算中,当t小于计算精度时,计算机识别为0,此时计算才会终止。

(b)x=pi/2; n =23;t= -1.2539e-018; s=1.0000 ;

x=11*pi/2; n =75;t= -2.6232e-017; s= -1.0000 ;

x=21*pi/2; n =121; t=6.4693e-018; s=0.9999 ;

2、已知数列{Xn}满足递推公式 Xn=10X(n-1)-1,n=1,2,...

若取X0=1.41(3位有效数字),问按上述递推公式,从X0计算到X10误差有多大?这个计算过程稳定吗?若不稳定请设计一种稳定的递推的算法。

在M文件窗口输入

function xk=suanfayi(x);

x(1)=1.41;

for n=2:10

x(n)=10*x(n-1)-1;

disp(x(n));

end

输出结果为

13.1000

130

1299

12989

129889

1298889

12988889

129888889

1.2989e+009

function xk=suanfaer(x);

x(10)=129888889;

for n=9:1:1

x(n)=(x(n+1)+1)/10;

disp(x(n));

end

t=0.1

n=1:10

e=n/10-n*t

运行结果:

t =

1/10

n =

Columns 1 through 5

1 2 3 4 5

Columns 6 through 10

6 7 8 9 10

e =

Columns 1 through 5

0 0 * 0 0

Columns 6 through 10

* * 0 0 0

五、实验小结:

陈平: 通过这次实验,让我又对MATLAB有了新的认识。数值计算中误差是不可避免的,截断误差、舍入误差等等。浮点数不仅是有限的,而且它们在实数轴上的分布也是不均匀的,更重要的事,这种分布不均匀性并不影响算法的讨论。但是浮点数是有限精度的,这会导致许多问题。比如“上溢”、“下溢”,由于有舍入误差,所以眼见未必为实噢!浮点数的加法在计算过程中由于有对阶,因此不满足结合律和交换律,通过实验我们得出这个结论。但是我们的思考绝不止于加法,我们进一步思考乘法呢?除法呢?也会像浮点数的加法一样不满足这两大运算规律吗?我们在编写算法时,要记住误差衰减的算法才是我们追求的。问题和算法一样,有“优”、“劣”之分。病态问题需要通过研究和构造特殊的算法来解决。思考题由于我们组的同学都十分热情,所以两个思考题都坚持做完。

这次实验让我们组每个同学都收益良多,团队的力量永远大于个体。

我们组的周红同学,细心研究浮点数,为我们这次实验增色不少。王双同学积极查阅资料,主动思考,为我们这组的实验奉献很多。贺嘉玲同学在实验中既耐心又有激情,总是积极思考,实验中提出很多十分棒的想法。作为这次实验的组长,我十分开心,我们这组的成员每一个人都很主动,大家一起讨论问题,找到解决问题的方法,并得出结论。对于每个问题,我们都应该积极去思考并解决,这样才能不断的学习到新的知识技能。思考+动手实践=真正的学习。

周红: 通过这次试验我知道了数值计算中应该注意的几个问题:

1.避免两个相近的数相减,容易造成有效数字的严重丢失

(解决方案:1.改变计算公式:更有效;2.多保留几位有效数字)

2.舍入误差对计算结果影响小,舍入误差不增长

3.避免大数加小数,这样小数容易被忽略,误差较大(此点在浮点数第二题有所体现)如:s=1./(n*(n+1));

s1=1/2+1/6+1/12……+1./(n*(n+1));

s2=1./((n+1)*n)+1./(n*(n-1))+……+1/12+1/6+1/2;

用s2方法更好.

4.先化简再计算,减少运算次数,避免误差积累

另外,浮点数的计算不满足交换律和结合律,计算多项式时,使用秦九韶算法可以降低复杂度,减少运算次数

王双:这次的实验题给我一种很陌生的感觉,之前没有做过也没有听说过,很多概念也看不太懂,就自己慢慢着磨通过跟同学讨论帮助了我了解。在做题的过程中发现问题解决问题,慢慢就有收获了。在这次实验中,我知道了关于浮点数的一些运算规则,打破了原有的错误的认识。对误差传播与算法稳定性有了一点了解,还有病态问题实验,运用了之前的基础可又要比之前的复杂些。

贺嘉玲:通过这次的实验,我对算法有了全新的认识。(1)对于浮点数而言,首先它的分布是不均匀的但是不会影响算法,其次它的的加法并不满足交换律和结合律,这就导致了相同的两个数据用大数加小数的误差会比用小数加大数的误差大,其根本原因是浮点数在

计算式首先做的事比较其阶码。(2)对于算法的复杂性而言,在计算多项式能够通过计算时间明显的感觉到算法优劣性对计算时间的影响。优秀的算法(如秦九韶算法)和普通的直接计算的算法相比,它能够讲效率提高一倍。(3)对于病态问题而言,算法的优劣能够直接影响到一些问题的运算结果。以这一次的求解多项式的根为例,“坏”的算法会让解严重偏离实际解,这时只能通过特殊的算法才能解决问题。(4)对于稳定性而言,是针对算法比较其误差是逐步扩大还是逐步缩小的,只有稳定的算法才能够算得上是优秀的算法。这一次的实验让我明白了每一个普通的运算其算法都凝聚着编程人员的心血,值得我们钦佩。同样,我们的每一次实验都凝聚着所有成员的心血。

相关文档