数值分析实验一：误差分析、误差传播及算法稳定性

毕节学院实验报告

实验名称： 误差分析、误差传播及算法稳定性 实验报告序号： 1 组 别

姓 名 朱海涛 同组实验者 周礼伟

实验项目

计算1

10

e e n x n I x dx -=?(0,1,)

n = 并估计误差

实验日期 2012年9月26日

实验类别

□√ 1、验证性实验或基础性实验； □ 2、综合性实验

□ 3、设计性实验； □ 4、创新性实验和研究性实验；

教师评语

实验成绩

指导教师（签名）

赖志柱

年 月 日

实验目的：

通过本实验对求解问题的算法进行好坏判断有一个初步了解，并加强对设计一个好算法的理解，体验数值计算稳定性，从而了解数值计算方法的必要性，体会数值计算的收敛性与收敛速度。

实验任务与要求：

计算11

e

e n x

n I x dx -=?

(0,1,)n = 并估计误差

（1）建立若干个（不少于两个）计算公式； （2）分析计算公式的理论误差；

（3）编写程序（推荐MATLAB ）实现（1）中的计算公式、输出结果并比较实际误差；

（4）任选正整数m n ，要求既从m I 计算n I ，又从n I 计算m I ，并分析您的结果。这里0m ≠且9n ≠。

小组分工合作说明：

实验过程及内容：

解：

由分部积分可得计算n I 的递推公式

1111

01,1,2,e 1.n n x I nI n I e dx e ---=-=??

?==-??? (1)

若计算出0I ，代入（1）式，可逐次求出 1

2

,,I I …

的值。

要算出0I 就要先算出1e -，若用泰勒多项式展开部分和

2

1

(1)(1)1(1),

2!

!

k

e

k ---≈+-+

++

…

并取k=19,用4位小数计算，则得

1

0.3679

e

-≈，截断误差

1

4

711|0.3679|10

8!

4

R e

--=-≤

<

?.计算过程中小数点后第5位的数字按四

舍五入原则舍入，由此产生的舍入误差这里先不讨论。当初值取为

00

0.6321I I ≈= 时，用（1）式递推的计算公式为 01

0.6321

A 1n n I I nI -?=?=-?

（），n=1,2，…。

计算结果见表1的n

I 列。用0

I 近似0I 产生的误差0

00

E I I =- 就是初值

误差，它对后面计算结果是有影响的. 从表1中看到18

I 出现负值，这与一切0n

I >相矛盾。实际上，由积

分估值得

1

1

1

1

1

01

01

1(im )(max )1

1

x

n n

n x x e

e m e x dx I e x dx n n ---≤≤≤≤=<<=

++?? （2）

因此，当n 较大时，用n

I 近似n I 显然是不正确的。这里计算公式与

每步计算都是正确的，那么是什么原因合计算结果出现错误呢？主要就是初值0

I 有误差0

00

E

I I =- ，由此引起以后各步计算的误差

n n n

E I I =- 满足关系 1,1,2,n n E nE n -=-=….

由此容易推得

0(1)!n

n E n E =-，

这说明0

I 有误差0E ，则n

I 就是0E 的n!倍误差。例如，n=19，若

4

01||10

2

E -=

?，则190||19!||2E E =?>。这就说明8

I 完全不能近似19I 了。

它表明计算公式（A ）是数值不稳定的。 我们现在换一种计算方案。由（2）式取n=19,得

1

19120

20

e

I -<<

，

我们粗略取1

*

99

1

1

(

)0.0684210

10

e

I I -≈

+

==，然后将公式（1）倒过来算，

即由*9

I 算出*8

I ，*7

I ，…，*0

I ，公式为

*

19**

10.0342()1(1),198n n I B I I n n -?=?

=?=-=??

，1， (1)

计算结果见表1的*

n

I 列。我们发现*0

I 与0I 的误差不超过410-。记

**

n n n

E I I =-，则**

1||||!

n E

E n =

，*

0E 比*n E 缩小了

n!倍，因此，尽管*9

E 较

大，但由于误差逐步缩小，故可用*n

I 近似n I 。反之，当用方案（A ）计算时，尽管初值0

I 相当准确，由于误差传播是逐步扩大的，因而

计算结果不可靠。此例说明，数值不稳定的算法是不能使用的。

程序如下：

function x11 = facto(n) %这个函数的功能是求 n 的阶乘； x11 = 1; if n == 0 x11 = 1; else

for i = 1:n x11 = x11*i; end end

function e_1 = telor(k)

%这个函数的功能是求e^(-1);用泰勒多项式展开式进行计算的， % k 是代表展开到第 k+1 项

e_1 = 1;

if k == 1

e_1 = 1;

else

for i=1:k

e_1 = e_1 +(-1)^i/facto(i);

end

end

function jifen(m)

I0=1-telor(19);

%第一种算法

I(1)=I0;

for i = 1:m

I1 = 1 - i*I0;

I(i+1)=I1;

I0=I1;

end

%第二种算法

Im=(1/2)*(1/(m+1)+telor(19)/(m+1));

B(1)=Im;

for i=1:m

In=(1/(m+1-i))*(1-Im);

B(i+1)=In;

Im=In;

end

disp(' n 第 1 种算法第 2 种算法 '); for i = 0:(length(B)-1)

fprintf('%4d %33.4f %12.4f |\n',i,I(i+1),B(m+1-i)); end

在Matlab命令窗口输入如下命令即可得到如表1的结果。

jifen(19)

表1 计算结果

n 第1种算法第2种算法

0 0.6321 0.6321

1 0.3679 0.3679

2 0.2642 0.2642

3 0.2073 0.2073

4 0.1709 0.1709

5 0.1455 0.1455

6 0.1268 0.1268

7 0.1124 0.1124

8 0.1009 0.1009

9 0.0916 0.0916

10 0.0839 0.0839

11 0.0774 0.0774

12 0.0718 0.0718

13 0.0669 0.0669

14 0.0627 0.0627

15 0.0590 0.0590

16 0.0555 0.0557

17 0.0572 0.0527

18 -0.0295 0.0508

19 1.5596 0.0342