1绪论-求根课件1 下载本文

第1章 绪 论

本章主要介绍科学计算的特点、数值分析基本知识和概念,它们对学习数值分析、了解科学计算原理,以及进行科学计算都是很有帮助的。

1.1 学习数值分析的重要性

例 1.1 将数列In?解:因为

?10xndx 写成递推公式的形式,并计算数列I1,I2,?的值。 x?5xn?5xn?1?5xn?1In??dx0x?5 n?111x1??xn?1dx?5?dx??5In?100x?5n1得到计算I n 的递推公式

In?由I0?1?5In?1nn?1,2,??1.1?

?1016dx?ln ,借助递推公式(1.1)可依次算出I1,I 2,??。 x?55 实际中,一般需要具体的数据,因此若取I0为准确到小数点后8位的近似值作为初始值,在字长为8的计算机上编程计算,可出现I12??0.32902110?102的结果,这显然是错误的!因为数列I n 的被积函数在积分区间上是非负的,故总应有I n≥0才对。(请读者可以

在自己的计算机上用递推公式(1.1)编程做一个数值实验,来检验当n较大时,I n的计算结果一定会出现负数的现象!)

现在很多科学研究和工程问题的解决都是借助计算机进行的。通常用计算机解决实际问题有四个步骤

1建立数学模型;○2选择数值方法;○3编写程序;○4上机计算。 ○

1.2计算机中的数系与运算特点

1.计算机的数系

数学理论告诉我们:实数集是稠密的无限集,其中任何一个非零实数可表示为

x??10c?0.a1a2a3???数学上可以方便的定义 ? 进制的浮点数

?1.2?

其中ai??0,1,2,3,4???,9?,c 为整数。式(1.2)表示的数x称为十进制浮点数。类似地,

x???c?0.a1a2a3??? ai??0,1,2,3,4???,??1?。

在计算机中,由于机器本身的限制,数学中的实数被表示为

x???c?0.a1a2a3???at ai??0,1,2,3,4???,??1? (1.3)

其中t是正整数,表示计算机的字长;c是整数,满足L≤c≤U,L和U为固定整数。对不

1

同的计算机,t、L和U是不同的,而?一般取为2,8,10和16。式(1.3)表示的数x称为t位?进制浮点数,其中c称为阶码,0.a1a2a3???at称为尾数,这样一些数的全体

F(?,t,L,U)???c?0.a1a2a3???atai??0,1,?????1?,L?c?U?

称为机器数系,它是计算机进行实数运算所使用的数系。

显然,机器数系是有限的离散集。该数集从总体看,在实数轴上分布不是均匀的,但从局部看,阶数相同的数又等距地分布在实数轴的某一段上。机器数系中有绝对值最大和最小的非零数M和m,例如在4位十进制浮点数系F(10,4,-99,99)中,绝对值最大的非零数M??10?0.9999,绝对值最小的非零数 m??10?0.0001。

若一个非零实数的绝对值大于M ,则计算机产生上溢错误,若其绝对值小于m,则计算机产生下溢错误。上溢时,计算机中断程序处理;下溢时,计算机将此数用零表示并继续执行程序。无论是上溢,还是下溢,都称为溢出错误。

通常,计算机把尾数为0且阶数最小的数表示数零。

2.计算机对数的接收与处理 实数集合是稠密无限的,而计算机使用的数系是有限的。为建立数学中的实数集和计算机中机器数系的对应,计算机采用下面方式进行这种对应。 设非零实数x是计算机接收的实数,则计算机对其的处理为 (1)若x?F(?,t,L,U) 则原样接收x ;

(2)若 x?F(?,t,L,U) 但 m?x?M,则用F(?,t,L,U) 中最接近x的数

99?99?fl(x)表示并记录x,以便后面处理。

计算机是怎样做数学运算的呢?首先要说明的是计算机本质上只能做加减乘除运算。两个数在计算机中参与运算的方式为:

(1) 加减法: 先对阶,后运算,再舍入; (2) 乘除法: 先运算,再舍入

计算机的运算一般是在计算机中央处理器的运算器中进行,而运算器中一般是多倍字长存储的,故计算机参加运算的数据允许有超出原机器数系字长的数出现。 例如,某计算机的数系F(10, 4,-99,99)的两个数x1=0.2337×10 -1和x2=0.3364×102 ,则运算过程如下

fl(x1?x2)?fl(0.2337?10?1?0.3364?102)对阶?fl(0.0002337?102?0.3364?102)fl(0.3366337?102)2运算?

舍入?0.3366?10运算fl(x1?x2)?fl(0.2337?10?1?0.3364?102)?fl(0.7861668?100)

舍入?0.7862?100?0.7862

在计算机的数系F(?,t,L,U) 中,把尾数的第一位a1≠0数称为规格化的浮点数。

2

1.3 误差

科学计算的实质是用近似的数据,通过近似有效地计算去获得可用的结果。其强调的是计算结果的可用性,而不是准确性.

1.误差的来源

误差可以来自很多场合,但其来源主要有4个方面 1).模型误差(也称描述误差);2).观测误差(也称数据误差); 3).截断误差(也称方法误差);4).舍入误差(也称计算误差)。

x2xnx2xnxe?1?x??????????近似公式 e?1?x??????

2!n!2!n!x舍入误差是计算机因数系有限,在接收和运算数据时引起的误差,它也是数值分析关注的

内容。本章例1.1的错误就是舍入误差造成的。

2.误差的定义 定义1.1 设x是准确值,x*是x的一个近似值,称差 x*- x为近似值x*的绝对误差,简称误差,记为e* 或e (x*) ,即e (x*)= x*- x。

由于准确值x通常是未知的,故误差e* 一般是计算不出来的,为此引入误差限来对误差进行估计。

*** 定义1.2 称满足 e?x?x??的正数? * 为近似值x*的误差限。

误差限一般是可以求出的,例如用具有毫米刻度的皮尺去测量某个物件的长度,测得的数据与物件的实际长度不会超过半个毫米。这半个毫米就是该物件长度的误差限。误差限 给出了准确值x 所在范围为x???x?x?? ,该范围常用x?x??表示。

误差限是误差绝对值的上界,显然误差限是不唯一的。由于合适的最小的误差限一般求不出来,实用中可以根据问题的特点选择一个合适的误差限即可。通常误差限一般取为获取该近似数仪器精度的半个单位。显然误差限?* 越小,近似值越精确。

绝对误差不能反映近似值x*的近似程度。例如某个量的准确值x1=9,其近似值x*1=10; 另一个量的准确值x2=999,其近似值x*2=1000,这两个量的绝对误差都是1,但显然x*2的近似程度比x*1好。为描述近似数的近似程度,需引入如下相对误差概念。

******e*x*?x 定义1.3 设x是准确值,x*是x的一个近似值,称? 为近似值x*的相对误

xxe*x*?x*?差,记为e*r或er (x*),即er?x??。 xx相对误差的绝对值越小,近似程度越高。例如,前面例子近似值x*1和x*2的相对误差

****分别为erx1?1/9和erx2?1/999 ,因为erx1?erx2 ,所以x*2比x*1逼近

????????程度好。

同样,由于准确数x通常是未知的,导致相对误差不可计算,为此引入估计相对误差的相对误差限概念。

x*?x 定义1.4 称满足 er???r*的正数?r* 为 近似数x* 的相对误差限。

x* 相对误差限不如绝对误差限容易得到,实用中常借助绝对误差限来估计。为方便估计

3

e*x*?xe*x*?x*?相对误差限,实用中常用*?代替er?x??进行误差限的计算,这样就*xxxx?**得到实用中容易计算的相对误差限?r?*。显然,该值越小,近似程度越高。

x相对误差在科学计算中的重要性比绝对误差大。

3.数值计算的误差

计算机中的数值运算本质上是加、减、乘、除四则运算,带有误差的数经过四则运算后误差的变化有如下估计关系

定理1.1 假设x* 和y* 分别是准确值x和y的一个近似值,则有: 1) 四则运算的绝对误差估计

e(x*?y*)?e?x*??e?y*?, e(x*?y*)?y*e?x*??x*e?y*?

e(x*y*e?x*??x*e?y*?y*)?? y*?22) 四则运算的相对误差估计

e?x*??y*er?y*?r(x*?y*)?x*erx, e*y*)?e***?y*r(x?r?x??er?y?

x*er(y*)?er?x*??er?y*?

证明 只证估计式e(x*?y*)?y*e?x*??x*e?y*?,其余类似可证之。

由定义有

e(x*y*)?x*y*?xy?x*y*?xy*?xy*?xy?y*?x*?x??x?y*?y??y*e?x*??xe?y*? ?y*e?x*??x*e?y*?(?x?x*)证毕。

e (x*)= x*- x = dx ,

e*x*?xr?x??x?dxx?dlnx 出近似数x*的绝对误差和相对误差与微分的关系

1)dx?e?x*? 2)dlnx?e*r?x?

由d(x?y)?dx?dy ,可得 e(x*?y*)?e?x*??e?y*?;

由 dln(xy)?dlnx?dlny,可得 e?x*?y*??e?x*??e?y*rrr?

例1.2 考查函数y = x n 的相对误差与自变量x的相对误差关系。 解 因为 lny?nlnx,所以dlny?ndlnx,得

4

er

??x???ne?x?*n*r

定理1.2设多元函数 u?f(x1,x2,???xn),向量自变量(x1,x2,???xn)的近似值为

***,则有多元函数f(x1,x2,???xn)的误差估计 (x1,x???,x)2n***?f?x1,x2,?,xn? 1)efx1,x2,?,xn??***????i?1n?xie?xi*?

2)?fx1,x2,?,xn??***????i?1n***?f?x1,x2,?,xn??xi***?f?x1,x2,?,xn???xi*?

*** 3)?rfx1,x2,?,xn??????i?1n??xi*?f?x,x,?,x*1*2*n?xi?*

证明 利用Taylor展式有

**e?u??du?f?x1,x2,?,xn*??f?x,1x,?2,xn?**?f?x1,x2,?,xn*?*?f?x1,*x,?2,xn?

??i?1n?xi?x*i?xi???i?1n?xie?x*i?

然后利用绝对误差、绝对误差限和相对误差限定义可得。

例 1.3 设有一长方体水池,测得其长、宽、深分别为50?0.01米,25?0.01米,20?0.01米,试按所给数据求出该水池的容积,并给出绝对误差限和相对误差限。

解 令L,W, H分别代表长方体水池的长、宽、深;V代表长方体水池的容积,有

V=V(L,W,H)=LWH

由题意有长方体水池的长、宽、深的近似值为:

L*=50米,W*=25米,H*=20米,?(L*)=?(W*)=?(H*)= 0.01米

按所给数据求出该水池的容积为:V*=V(L*,W*,H*)=L*W*H*=50?25?20=2500 (米3)

?V*?V*?V***?e?V??e?L??e?W??e?H*??L?W?H**?W*H*e?L*??L*H*e?W*??LWe?H*?*?V*?V*?V***???V????L????W????H*??L?W?H**?W*H*??L*??L*H*??W*??LW??H*?*?25?20?0.01?50?20?0.01?50?20?0.01?27.50?米3???r?V*????V*?V*?27.50?0.11% 2500故本题绝对误差限和相对误差限为27.50立方米和0.11%。

4. 计算机的舍入误差

5

设计算机的数系为F(?,t,L,U) ,m及M是其中绝对值最小及最大的正数,某数

x???c?0.a1a2???,a1?0满足m?x?M,x? F(?,t,L,U),则计算机经舍入处理?, 后以数fl(x) 接收,即fl(x)???c?a??0.a1a2???at?a???t??0.a1a2???at??1) efl?x??x?fl(x)?0.5??2)erfl?x??0?at?1??/2at?1??/2

因此计算机对x的舍入绝对误差和舍入相对误差有如下估计

??c?t

??x?fl(x)x0.5??c?t??0.5??1?t c0.1??由此可知,计算机对任何实数的舍入相对误差限与实数本身无关,只与计算机字长t有关,其值为0.5??1?t。因此通常定义数eps?0.5??1?t 为计算机的精度。给定一台计算机后,其舍入相对误差的精度也就给定了,任何运算要想得到比计算机精度更小的舍入相对误差,需要做特别处理才行。由于计算机的精度只与字长有关,计算机字长t越大,其精度越高。这也是在数值计算中,为了提高计算结果的精度,有些数值要用双字长处理的原因,双字长数据也称为双精度数。

1.4有效数字

定义1.5 若近似数x* 的误差限是其某一位上数字的半个单位,就说近似数x* 准确到该位;由该位自右向左数到x* 的第一个非零数字若有n位,就称近似数x*有n位有效数字。

由定义可知,若近似数x*有n位有效数字,只要最后一位有效数字不是零,则可以认为与准确数x相比, x*的前n-1位数字都是准确的,只是最后一位可能比x在该位数的值少1。因此,可以从有效数字推断出准确数的主要数据部分。

有效数字的数学描述是:设x*??0.a1a2???ak?10m,a1?0,al??0,1,2,???,9?,m为整数,k为不小于正整数n的整数。若有关系式

e?x*??x*?x?0.5?10m?n(1.4)

就称近似数x*有n位有效数字,此时x*有n位有效数字的值为?0.a1a2???an?10m。 式(1.4)常用来确定近似数有多少位有效数字。

可以证明,如果十进制准确数x经过四舍五入得到近似数x*,则x*的有效数字位为将x*写为规格化浮点数后的尾数的位数,例如x =0.00345, 四舍五入得x*=0.0035=0.35?10-2,可知x*有2位有效数字。

由定义可知,有效数字越多,绝对误差和相对误差就越小,因此近似数就越准确,这也是科学计算中要尽可能多保留有效数字的原因。

例1.4 求圆周率??3.1415926??? 的近似值x1?3.14和x2?3.141的有效数字。 解 x1?0.314?101 x2?0.314?1110m,?

1 6

由??x1?0.015926????10?2?0.15926????0.5?10?2,有m-n =-2 ,得n =3 ,即x1有3位有效数字;

由??x2?0.0005926????10?3?0.5926????0.5?10?2,有m-n =-2 ,得n=3 , 即x2有3位有效数字。

例 1.5 已知近似数x*有5位有效数字,试求其相对误差限。

解 因为x*有5位有效数字,可以设x*??0.a1a2???a5?10m,a1?1,于是有n=5和

x*?x?0.5?10m?5,考虑x*的相对误差:

x*?xx*0.5?10m?55?10?511?4?4????10??10 m0.a1a2?a5?10a12a12故有5位有效数字的x*相对误差限为0.5?10 -4 。

有效数字与相对误差有如下关系:

定理1.3设近似数x*??0.a1a2???ak?10m,a1?0,al??0,1,???,9?,m为整数,有: 1) 若x*有n位有效数字,则x*的相对误差

er?x2) 若x*的相对误差

*??x*?xx*?1?101?n2a1(1.5),

er?x*??则x*有n位有效数字。

x*?xx*?1?101?n2?a1?1?m?n(1.6)

*证明 1) 因为x*有n位有效数字,则有x?x?0.5?10,于是

er?x2) 由

*??x*?xx*0.5?10m?n0.51?n1?n???10??10 m0.a1a2?ak?100.a12a1x*?xx*?*1?101?n,有

2?a1?1?*0.a1a2?ak?10m11?nx?x?x??10??101?n2?a1?1?2?a1?1?a.a?ak?12?10m?n2?a1?1?证毕。

a1.a2?ak?a1?1?1?10m?n2

例 1.6 为保证某算式的计算精度,要求参与计算的323的近似值x*的相对误差小于0.1%,请确定x*至少要取几位有效数字才能达到要求。 解 先将323写成浮点数。 因为2?3323?,所以323?2.a2a3??0.a2??1,得到10a1=2。假设x*至少要2a3取n位有效数字才能保证相对误差小于0.1%,由定理1.3的(1.5)式,选择满足

111?101?n?0.1%得?101?n??101?n?0.1%的最小整数n即可。由

2?22a12?2

7

104?n?4,因而要n?4才行。满足此式的最小整数n为4,故x*至少要取4位有效数字

才能达到相对误差小于0.1%的要求。

例 1.7 若x =1.1062, y =0.947是经过四舍五入后得到的近似数,求x-y和xy的绝对误差限,指出x-y和xy各有几位有效数字并给出其用有效数字表示的值。

- 4 -3

解 由有效数字定义有:|dx|? 0.5?10 ,|dy|? 0.5?10,因为

|e(x-y)|= |d(x-y)|? |dx-dy|?|dx|+|dy|? 0.5?10-4+0.5?10-3<0.5?10-2

故x-y的绝对误差限为0.5?10,也说明x-y可以准确到小数点后两位。

直接计算有x-y=0.1592的值有2位有效数字, 利用四舍五入法,得出x-y具有2位有效数字的近似数为0.16;由

-2

|e(xy)|= |d(xy)|? |ydx+xdy|?|y||dx|+|x||dy|?0.947?0.5?10-4+1.1062?0.5?10-3<0.5?10-2

故xy的绝对误差限为0.5?10,说明xy可以准确到小数点后两位。直接计算有xy=1.04757的值有3位有效数字,利用四舍五入法,得出xy具有3位有效数字的近似数为1.05。

1.5数值分析研究的对象、内容及发展

数值分析是一门最古老的数学,因为它只涉及到数的四则运算,但它又是一门很年轻的数学,因为它是随计算机的诞生、发展和广泛应用而逐步发展形成的,其理论还涉及到现代数学的理论和内容。

数值分析属计算数学的范畴,有时也称它为计算方法、数值方法、计算数学等,其研究对象是各种数学问题的数值方法设计、分析,有关的数学理论和软件实现,它是一个数学分支。

数值分析的内容可分为两大方面:1.连续系统的离散化;2.离散型方程的数值求解。

-2

求定积分

?f(x)dx 是连续系统问题,数值分析常用公式 ?f(x)dx??Abbk?1bbnkf(xk)

微分方程的定解问题

??y'?f(x,y)(a?x?b) ?ya?y?a???令ym?y(xm),h?(b?a)/n,xm?a?mh(m?0,1,2,???,n)可转化为离散化方程

ym?1?ym?hf(xm,ym)

这样由y0?ya,可依次求出y1,y2,???,yn,从而求函数y?x?在xm的近似值

y?x1?,y?x2?,???,y?xn?,然后数据逼近方法求出y?x?的近似解。

串行计算方法或串行算法 ,

1.6 1. 数值问题

由一组已知数据(输入数据),求出一组结果数据(输出数据),使得这两组数据之间

8

数值分析中常用的一些概念

满足预先指定的某种关系的问题,称为数值问题。 2. 数值解

由近似公式计算出的解称为数值解。一般数值解是近似解。 3. 算法

由给定的已知量,经过有限次四则运算及规定的运算顺序,求出所关心未知量的数值解,这样所构成的整个计算步骤,称为算法。 数值分析本质上是研究和构造算法。

4.计算量

一个算法所需要的乘法和除法总次数称为计算量,常用N表示。

计算量的单位为flop,表示完成一次浮点数乘法或除法所需要的时间。算法的计算量可以衡量算法的优劣,因为它体现着算法的计算效率。

算法的计算量越小,则算法的计算效率越高,因而该算法也就越好。

注:在以往由于计算机做加减法通常要比乘除法快得多,故在一般情况下算法的计算量可以不考虑加减法的计算时间。但目前由于计算机技术的发展,计算机做加减法与做乘除法耗费的时间已相差不大,因此,一些数值分析的教材在计算量定义中加入了加减法次数。考虑到数值分析是一个科学计算的基础部分,为方便计,本书的计算量定义采用通常的做法。

例如假设A,B,C 分别为10×20,20×50,50×1的矩阵,计算D=ABC就有如下不同的算法和计算量:

算法1 D=(AB)C 计算量 N=10500 flop 算法2 D=A(BC) 计算量 N=1200 flop

显然算法2的计算量比算法1小,因而算法2比算法1要好。

5. 病态问题和良态问题 1) 病态问题

因初始数据的微小变化,导致计算结果的剧烈变化的问题称为病态问题。

病态问题也称坏条件问题,这类问题通常是问题本身固有的,其在函数计算、方程求根及方程组求解中都是存在的。例如,线性方程组

1111?x?x?x??122336?1113?1 (1.7) ?x1?x2?x3?23412?1147?1x?x?x??31425360?的准确解为 x1=x2=x3=1 ,把它的系数都舍入成两位有效数字做小的扰动后,原方程组变为

?x1?0.50x2?0.33x3?1.8??0.50x1?0.33x2?0.25x3?1.1 ?0.33x?0.25x?0.20x?0.78123?这个方程组的准确解为 x1=-6.222 ,x2=38.25,x 3= -33.65,此解与扰动前的解完全不同了。方

程组(1.7)的求解就是病态问题。

9

因为实际问题的数据都是近似的或经计算机做了舍入处理的,这都会引起原始数据的扰动。若所求解的问题正好是个病态问题,则采用通常算法计算就会出现很隐蔽的错误,导致不良的后果。因此,在做科学计算时要特别注意所求解的问题是否为病态问题。

病态问题的计算或求解应使用专门的方法或将其转化为非病态问题来解决。 2)良态问题

初始数据的微小变化只引起计算结果的微小变化的计算问题称为良态问题。

?2x1?x2?6?2x1?x2?6 例如对方程组?的常数项做微小扰动后变为?

x?2x??2x?2x??2.00522?1?1扰动前方程组的准确解为 x1=2 x 2 = -2,而扰动后方程组准确解为x1=1.999,x2=-2.002。

这两组解之间的差别是不大的。

数值分析中主要研究良态问题的数值解法。 6. 数值稳定算法

如果一个算法进行计算的初始数据有误差,而在计算过程中产生的舍入误差不增长,则称该算法为数值稳定算法,否则称为数值不稳定算法。 例1.8设计算机的数系为F(10,4, L,U),今有一批数据

x1?0.5055?10,x2?x3?????x11?0.4500 ,试求其和S??xi。

4i?111 解 算法1. 按 x i 角标由小到大的顺序计算

fl(x1?x2)?fl(0.5055?104?0.4500?100)?fl(0.5055?104?0.000045?104)

?fl(0.505545?104)?0.5055?104?x1同理,有fl((x1?x2)?x3)?fl(fl(x1?x2)?x3)?fl(x1?x2)?x1???

最后得S?0.5055?10

算法2. 按x i 的角标由大到小的顺序计算:

4fl(x11?x10)?fl(0.4500?100?0.4500?100)?0.9000?100

fl(x9?fl(x11?x10))?fl(0.4500?100?0.9000?100)?0.1350?101……

最后得S?0.5060?10

例1.9计算数列In? In?4?10xndx,n?1,2,?,100 值。先考察计算该数列的递推公式 x?51?5In?1的数值稳定性。如果该公式不是数值稳定的,试给出一个数值稳定的算法来n完成计算并证明其稳定性。

11**?5In?1,故其带有舍入误差的对应计算公式为In??5In?1,二式相nn**n**减得In?In?5In?1?In?1?5I0?I0,n?1,2,?,该计算公式由I0开始,依次计算出

解 因为 In?**,其在每次计算过程中,都将上次计算的舍入误差放大5倍,从式中可知计I1*,I2,?,I100算I n时的舍入误差为初始误差的5 n倍!因此该算法不是数值稳定算法。

10

下面采用逆序计算方式给出一个数值稳定的计算公式 将In?11In?5In?1,n?1,2,?,100转换为In?1??,n?100,99,?,1n5n5*,该计算

***公式由I100开始,依次计算出I99,I98,?,I1*,其舍入误差关系为In?1?In?1?1*In?In。5这说明,在每次计算过程中,都将上次计算的舍入误差缩小5倍,因此该算法是数值稳定算法。

*

该算法的关键是取计算的初值I100, 这里采用如下定积分估计的方法选取: 1001x11110011因为 , ?I100??dx?xdx???006?101x?5??5101???5?5?10111?11?11????0.001815,其初始误差为 ??2101?56?60601?11?*?3 I100?I100???0.00033??0.331?10??101?56??3**选择好后,用此做递推计算,依次计算出I99它们所有误差不会超过0.331?10。 ,I98,?,I1*,

*取均值I100?

1.7 科学计算中值得注意的地方

这4个值得注意的地方是

1. 避免两个相近的数相减; 2. 避免用接近零的数做除数; 3.控制舍入误差的积累和传播;4.简化计算过程。

?x?ye?x??xe?y?xxdx?xdy,所以有商的误差关系式 e???22yyy?y??x?当 y?0 时,分母y2就会更小,从而使e??可能变得很大。要避免以上两种情况出现,

y??因为d()?xy通常利用数学手段将相应的计算转化为等价的或近似的另一个公式来计算的方式。 例 1.10怎样计算下列算式更好? 1) S=tan2.01-tan2;

1,x??1;

x?1?x1?x?ex,x??1 3) 3)y?2x2) 2)y?解1)是两个相近数相减问题。为避免直接相减,利用数学公式有关系式

sin(2.01?2)sin0.01 S?tan2.01?tan2??cos2.01?cos2cos2.01?cos2使用

sin0.01来计算S,就避免了两相近数相减问题,由它得出计算结果会有更多

cos2.01?cos2的有效数字,因此比直接计算tan2.01-tan2要精确。

2)当x>>1表示x很大,此时2)式分母 x?1?x?0且是两个近似数相减的问题。因此不能直接计算。采用数学变形得到更好的计算公式:

11

y?1?x?1?xx?1?x?x?1?x??x?1?x??x?1?x 上式没有分母接近零的数作为除数和两相近数相减的计算发生。

3)当x<<1表示x接近于零,此时3)分母出现接近于零的数。因此不能直接计算。该式找不到数学变形公式,采用Taylor公式(取展开式的前3项)做之

1?x?ex1xx2e?31xx2y??????x????2x23!4!5!2624

求解的问题有时会涉及若干计算公式,若机械地直接按公式来顺序计算,计算量可能

会很大。要说明的是:某个算法的计算量大并不能保证计算结果的准确性,相反可能会出现舍入误差的积累,使计算结果失真。计算量大不但浪费机时,还费力不讨好,因此要尽量减少计算量。在这方面可用数学运算特点及公式推演将计算公式简化,并利用计算机数据运算及存贮特点达到减少计算量的目的,从而可用较少的时间得到更准确的计算结果。如计算x64在某点的值,数学上有如下算法 算法1: x13?x?x?x???x; 算法2 :x13?x2?x4?x8/x; 算法3: x13?x??x?x2?2

对算法1,其计算量N=12 flop;对算法2,由于x2k可用xk. xk计算,即由xk再用一次乘法就算出x2k,它的计算量N=6 flop;同理算法3的计算量为N=5flop。因此这3个算法中,算法3最好。

思考题

1. ?数值2.3569(2.3569?*,?<5;Or 2.3568?*,??5) 2.3590(2.3590?*,?<5;or2.3589?*,??5), 2.3500(2.3500?*,?<5;or2.3499?*,??5) 346782(34678?*,?<7;or34677?*,??7)

是一些具有5位有效数字的近似数,原来的准确数可能是多少?

2. ?设数列y0?28,按递推公式 yn?yn?1?783?10?2计算到y100,若783取5位有效数字参与计算,试问计算到 y100的绝对误差限是多少? 该计算过程稳定吗? 设y0?28,按递推公式 yn?yn?1???21783,100n?1,2,???

计算到y100,若取783?27.982(5位有效数字),试问计算到y100将有多大误差?

1??yn?yn?1?783解 ?100?y0?28?n?1,2,?? (1)

12

??yn?28?设 ??1?y?y??27.982,nn?1?100?n?1,2,??? (2)

?记en?yn?yn,将(1)和(2)相减,得

e0?0??1 ?e?e??(783?27.982),nn?1?100?递推可得 en??因而

n?1,2,???

n(783?27.982),100n?1,2,???

e100??(783?27.982)

e1001??10?32

3. 用99的5位有效近似数来计算10?99果各有几位有效数字:

??10,分析下面两种不同的算法得出的结

? 2)?10? 1)10?99??0.996247047?10

99???10?99??0.10013658?1010?1310?10?12

这个结果给你什么启示?

4. 由于计算机硬件的发展,目前计算机在做计算时,乘除法所用的时间与加减法相差不大,因此目前也有把算法的计算量定义为一个算法所需要的加减乘除法总次数称为计算量。若用此定义,两个矩阵的相乘的计算量是多少?

5请提出一个与本章内容有关的思考问题,并尽量尝试给出解答。

数值实验?

n1x1dx,n?1,2,? 值,观察随着n 用递推公式 In??5In?1编程计算数列In??0x?5n的增大总会出现负值情况。

习题一

1.

?

3?1.732050???的三个近似值分别为

x1?1.73,x1?1.7321,x1?1.7320,这些近似数各有几位有效数字。

2. ?下列各数是按四舍五入方法获得的近似数,指出它们各有几位有效数字?

10.006,-0.0105,0.3140×10 3,0.134×10 3

3?.下列各式如何计算更好?

13

11?x (当x??1时) ?1?2x1?x (2) y?sin??sin?(?与?很接近)

(1) y?11?x? (当x??1时) xx (4) y?107(1?cos2?) (5)y?x?tanxx??1

(3) y?x??b?b2?4ac 4. ?用一元二次方程求根公式x1,2?在字长为8位的计算机上求一元二次

2a299方程 x??10?4?x?4?10?0的根,将求出的计算解与方程的准确解做对比,对你的

计算结果给出解释。

5.若ac??1,能否用常用的求根公式x1,2出说明及更好的求根公式。

?b?b2?4ac来求根?若不能,请给?2ae*x*?xe*x*?x?6. ?说明把相对误差的计算公式er?x??用公式*?来代替的合xxxx**理性,并指出这种替代的条件

erer2? er?er?1?er1?er

证明 设x为精确值,x为其近似值,

?2x??xx??x, er(x)? e(x)?x?x, er(x)? xx??x??xx??xer?er??? xx

11*?)(x?x)?xx *2(x?x)?xx*?(因而

2x??x2xx er?er?( )??e?rx*x*x?e(x)21er? ?er? e(x)1?er1?x2或

14

x??x2xx*2 er?er?( )?*?er?*x*xx?e(x)1e?r ?er?e(x)1?er1?*x22

7. ?在计算机数系F(10,4,-77,77)中,试求x1=0.14281×103和x2=-0.314159×101的机器浮点数fl(xi) ,(i =1,2)及相对误差。

8. 机器数系 F(10,8,L,U)中的三个数为

x?0.23371258?10-4,y?0.33678429?102,z??0.33677811?102 试按(x+y)+ z和x+(y+z)两种算法计算x+ y+ z的值,并将结果与精确结果比较。

?? 9?.表达式 ?3?8?与如下各式相等

?3?8???3????19601?69308?17?683,17?68?1???3,3?8???6,3?8?6,,

,19601?69308试说明在数系F(10,5,-50,50)中用上述表达式中哪一个获得的结果最好。

10000 10.选用更好的方法计算S??1。

2j?2j?1 11. ?设n次多项式

Pn(x)?anxn?an?1xn?1?????a1x?a0 试构造一个计算Pn(x)的算法,使其计算量尽可能小。 12?设 数列In?e?1?10xnendx,n?1,2,?, .证明 In?1?nIn?1,n?1,2,?并讨论该

算法的数值稳定性。请偿试给出一个计算此数列的数值稳定的算法并证明其稳定性。

第2章 非线性方程的求根方法

2.1引例

门的气压控制是指利用汽缸的活动来控制门的开合。设由汽缸控制的关闭状态的门如图2.1所示,图中z为门枢,a为门宽,x为门销,活赛直径为d,汽缸长为lo,汽缸内气体

15

压强为p0。当用力F推门时,则门会打开一个角度?,同时活塞下降距离c,且门销与z的水平距离始终不变化,汽缸内的压强增大为p,如图2.2所示。如果已知在绝热条件下,当绝热系数 ? =1.4,a =0.8m, b =0.25m, d =0.08m, lo=0.5m时,求施加了力F=25N后,门打开的角度?。

得到关于角度?的方程

?4Facos?(l0?btan?)???d2bp0l0?0

将已知数据代入,并整理,得关于?的函数方程

f(a)?(1?0.5tan?)1.4cos??0.32?0

因而所求问题变成求满足上述方程根的问题。不过这个根并不能用公式的方法求出。

上述案例归为求一元函数f(x)的零点或一元函数方程根的问题。怎样求函数方程的根,特别是求非线性方程的根,是本章要研究的内容。

本章探讨求函数方程根的常用数值方法的构造和原理,主要介绍非线性方程求根方法的有关知识和方法,重点论述二分法、简单迭代法、牛顿迭代法及其变形的原理、构造、收敛性等内容。

2.2问题的描述与基本概念

定义2.1 设f(x)为一元连续函数,称方程

f(x)=0 (2.1)

为函数方程;特别当f(x)不是x的线性函数时,称对应的函数方程为非线性方程。

非线性方程中,当函数f(x)是多项式函数时,称为代数方程,否则称为超越方程。称n次多项式构成的方程

anxn?an?1xn?1?????a1x?a0?0 为n次代数方程。

若存在一点x*使方程(2.1)成立,即f(x*)=0,则称x*是方程(2.1)的根,于是非线性方程求根问题,就是求使方程(2.1)成立的点x*,也即求函数f(x)的零点x*。例如

?b?b2?4ac一元二次方程P2(x)?ax?bx?c?0有求根公式x1,2?,那么x1, x2的

2a值就是一元二次方程的根,同时也是函数P2(x)的零点。

一般情况下,在非线性方程中,绝大部分是没有求根公式的,因此寻找有效的求方程(2.1)根的方法,哪怕是求近似根的方法是非常重要的。

2求根问题由根的存在性、根的范围和根的精确化三个子问题组成。

数值分析中所采用的求根方法一般分为区间法和迭代法两类。区间法,迭代法 构造收敛于根的数列?xk?,不同的是产生?xk?的方法不一样。

对同一个非线性方程可以用多种方法求根。衡量一种求根方法的优劣,可由其产生的数列?xk?收敛于根x*的快慢决定,为此引入如下收敛阶的概念。 定义2.2 设数列?xk?收敛于根x*,若存在正数 p 和C,满足

limk??xk?1?x*xk?x*p?C

16

则称?xk?的收敛阶为p或?xk?收敛于x*的速度是p阶的;同时称产生?xk?的求根方法具有p阶敛速。特别当p=1时,称方法线性收敛;当p=2时,称方法平方收敛;p>1时称方法超线性收敛。C称为渐近误差常数。

收敛阶实际上是微积分中的无穷小阶。令ek?xk?x*,则?k?epk,?k?ek?1在k趋于无穷大时是两个无穷小,对它们进行比较可知若?xk?的收敛阶为p ,则有? k是? k的p阶无穷小。由此可知,收敛阶p越大,则数列?xk?逼近x*的速度越快,从而该方法也就越好。因此收敛阶是评价求根方法好坏的重要标准。

本章所讨论的求根方法一次只能求一个根,对要求方程多个根的问题,可利用函数f(x)的特点先将根进行隔离,使每个区间只含有一个的根,再用这里的求根方法分别求之即可。

2.3 二分法

二分法也称对分区间法、对分法等,是最简单的求根方法,属于区间法求根类型。 基本思想

利用连续函数零点定理,将含根区间逐次减半缩小,构造出收敛点列?xk?来逼近根x*。 1. 构造原理

依据:若f(x)?C?a,b?,且满足f(a)f(b)<0,则在区间(a,b)至少有一点 ?,使

f(?)?0。二分法求根数列构?xk?的造过程为:

① 记含根区间I0?[a,b],取I0的中点x0?0.5(a?b),并计算f(x0); ② 判别f(x0)的值:

若f(x*0)?0,则x?x0,终止;

若f(x0)f(a)?0,则x*?[a,x0],取I1?[a,x0]; 若f(x0)f(a)?0,则x*?[x0,b],取I1?[x0,b]。 ③ 记I1?[a1,b1],取I1的中点x1?0.5(a1?b1),

判别?x1?x0???是否成立,?为给定的精度要求, 若成立,取x*≈ x1,终止;否则用I1代替I0,转①。 按上述步骤求根的方法称为二分法。

若记第k次二分区间处理得到的含根区间为Ik?[ak,bk],则有二分法对应的求根数列

算式为 xk?0.5a(k?bk,)k?0,1,?2,。 2.分析

因为根x*?Ik,k?0,1,?,且I0?I1?I2??,所以有xk的误差满足

x*?xk?12(b?a11kk)?22(bk?1?ak?1)???2k?1(b?a) (2.2) 17

于是当k??时,由是收敛于根x*的。

1(b?a)?0得到xk?x*,这说明由二分法产生的数列?xk?总k?121(b?a)??即可,解出k,有 k?12* 给定精度??0后,要x?xk??成立,取k满足

k?ln(b?a)?ln??1 (2.3)

ln2这样就可以保证进行k次二分计算得到的 xk 就是满足精度要求的近似根。

式(2.3)确定的k往往偏大,它主要用于理论估计。实用中一般用相邻的xk 与xk -1

差的绝对值是否小于?来决定二分区间的次数。这是因为ak或bk本身就是xk -1,由二分法的构造,有bk?ak?2xk?xk?1,故总成立有

x*?xk?xk?xk?1 (2.4)

*因而当xk?xk?1??时就有x?xk??,即此时xk是满足精度的近似根。

二分法的收敛速度不易用定义说明,注意到一般情况下xk+1比xk更靠近根x*,即有

xk?1?x*?xk?x*,故常说二分法具有线性敛速。

二分法的优点是不管含根区间[a,b]多大,总能求出满足精度要求的根,且对函数f(x)的要求不高,计算简单;缺点是不能求重根,且其收敛速度在数列xk越靠近根时越慢。

二分法一般用于为方程提供初始近似值。

例2.1用二分法求方程2x?5x?1?0在区间[1,2]内的根,绝对误差??10。 解 令f(x)?2x3?5x?1,则f'(x)?6x2?5,易验证f(1)?0,f(2)?0,且f'(x)在[1,2]内不变号。因此方程f(x)?0在[1,2]内有唯一根,由二分法算法有

3?2I0?[1,2],x0?0.5(1?2)?1.5

由f(x0)f(1)?0,得 I1?[x0,2]?[1.5,2],x1?0.5(1.5?2)?1.75 由f(x0)f(1)?0,得I2?[1.5,x1]?[1.5,1.75],x2?0.5(1.5?1.75)?1.625 类似计算,可得I4?[1.625,1.6875],x4?1.65625

8x7 I5?[1.6525,1.6,?]1.6718,75x5?x4?0.015625?10?2 55687?51.]6796, I6?[1.671875,1.,x688x6?x5?0.007813?0.7813?10 得近似解 x*?x6?1.679688

例2.2 写出二分法求非线性方程根的算法。

解 二分法算法可以由计算过程和误差控制给出,设非线性方程为f(x)?0,??0为给定的精度,因为f(x)是实数,为描述其为零的情况,引入充分小的数?1?0,用

?2?10?2

f(xk)??1表示f(xk)?0。此外,借助计算机的存贮特点,将二分法中的数列ak都存在

变量a中,bk都存在变量b中,用变量x记录xk,由此得二分法算法

输入 a,b,f(x),?,?1

输出 根x

步骤:① x?0.5(a?b)

② If ?f(x)???1, then 输出根x,stop ③ If f(a)f(x)?0,then b?x else a?x

④ If b?a?? then 输出根x,stop else 转①

18

2.4 简单迭代法

简单迭代法也称逐次代换法,是非线性方程求根方法中各类迭代法的基础,其中的概念、理论和处理技术可以很方便地应用到方程(组)求解问题中。

基本思想 利用对方程做等价变换根不发生变化的特点,将方程f(x)?0等价变形为x??(x),获得迭代计算公式xk?1??(xk),再由它算出逼近根x*的数列?xk?。

1. 构造原理

把方程f(x)?0等价变为方程x??(x),若x*为根,由f(x*)?0,有x*=?(x*)成立,若有?(x0)?x0,则x*?x0,根找到了。不过一般没有那么巧,通常有?(x0)?x0,此时,记x1??(x0),再用x1继续试探,如此重复,此过程可以描述为

① 将方程f(x)?0按某种方法改写成另一等价形式x??(x); ② 构造迭代公式xk?1??(xk);

③ 取定一个初值x0由迭代公式算出数列x1???x0?,x2???x1?,?。

按上述过程获得数列?xk?称为迭代数列,而称构造迭代点列公式?xk?的函数?(x)为迭代函数。用如上迭代数列来求根的方法称为简单迭代法。

x*点形象的称为?(x)的不动点,而称方程x??(x)为不动点方程。

2.简单迭代法的几何意义

方程x??(x)的根,在几何上就是直线y?x与曲线y??(x)交点的横坐标x,如图2-4所示。

*

(1) (2)

19

(3) (4)

图2-4简单迭代法的迭代图示

图2-4给出了简单迭代法几何解释,其中(1)、(2)是迭代收敛图示;(3)、(4)是迭代发散图示。

3.分析

假设迭代数列?xk?是收敛的,不妨设limxk?x,注意到迭代函数?(x)一般是连续的,

k??利用连续函数与极限的关系,有

x?limxk?1?lim?(xk)??(limxk)??(x)

k??k??k??这说明迭代数列?xk?的极限就是所求的根,即x?x,故用简单迭代法求根是可行的。

*由方程f(x)?0,可以得到很多的不动点方程,例如方程x?x?1?0,它的不动点方程有

5x?x5?1,x?x?4?x?3,x?51?x等等,

这说明对同一个方程f(x)?0可以构造很多的迭代格式的每个迭代格式产生的迭代数

列都是收敛的呢?答案是否定的。

方程x?10?2?0在x?1附近有一个根,若选择迭代方程x?lg(x?2)建立迭代格式xk?1?lg?xk?2?,取x0?1进行迭代计算得:

xx1?0.477121,x2?0.393947,?,x13?0.375812,x14?0.375812

显然此数列是收敛的,且得到1附近的近似根0.375812。但如果选择不动点方程x?10?2建立迭代格式xk?1?10k?2,也取x0?1进行计算,就有

99999998?2,? x1?8,x2?99999998,x3?10xx显然这个数列肯定不收敛。

上面的例子提出能否事先判断迭代数列收敛的问题,这里给出一个判别收敛的充分条件。

定理2.1 设迭代函数?(x)满足下列两个条件:

① 当x?[a,b]时,有?(x)?[a,b];

② 任取x1,x2?[a,b],存在与x1和x2无关的正常数L?1,满足

? (2.5) ?(x1)??(x2??)?Lx??x12则?(x)在[a,b]中有唯一的不动点x,且迭代公式xk?1??(xk)对任取x0?[a,b],产生的数列?xk?都收敛于x。

**证明 易证迭代函数 ?(x)?C[a,b]。作辅助函数

20

?(x)?x??(x)。

)b(?)。由中值定理,至少存在一个0显然?(x)?C[a,b]。由条件①知 ?(a?,即???(?),这说明?(x)在[a,??[a,b,使]?(?)?0b上有不动点]?。

如果?(x)在[a,b]上还有一个不动点?,???(?),利用条件②,有 ???(?)L????????????*矛盾,这就证明了满足定理条件的?(x)在[a,b]中有唯一的不动点,记为x。

*为证?xk?的收敛性,利用x是不动点、xk的迭代格式及条件②,有

????????(?)?xk?x*????(xk?1)??(x*)??L?xk?1?x*??L2?xk?2?x*????Lk?x0?x*?

k注意到0?L?1,在上式中令k??,可得L?0,于是,得出

?xk?x*??0 lim, 即有 limxk?x。

k??k??*定理得证。

例2.3证明迭代格式xk?1?3?0.5xk,x0??15产生的数列是收敛的。

证明 由迭代格式可知迭代函数为??x??3?0.5x,取其定义区间为实数R,显然有,任取x?R???x??R,另外任取x1,x2?R有

??(x1)??(x2)????????x1????????x2敛于其唯一根x,当然有本题结论成立。

*??0.5?x1??x2?0.5?x1?x2?

取L=0.5<1,则由定理2.1有任取x0?R,迭代格式xk?1?3?0.5xk产生的数列?xk?都收

定理2.1中条件②不易使用。注意到若?'(x)存在,且有??'(x)??L?1,x?[a,b],则由中值定理可得??'(x)??L?1,?(x1)??(x2)????'(?)??x1?x2??L?x1?x2?,这样由?x?[a,b],可以得到条件②,于是有如下:

推论2.1 设迭代函数?(x)满足

① 当x?[a,b]时,有?(x)?[a,b],② ??'(x)??L?1,x?[a,b]

则?(x)在[a,b]上有唯一的不动点x,且对任意初值x0?[a,b],迭代公式xk?1??(xk)产生的数列{xk}都收敛于x。

对于一个迭代,让迭代函数的导数在较大的范围满足??'(x)??L?1一般不容易做到,但找到较小的范围满足??'(x)??L?1还是较容易的。通常称在根的较小邻域内取初值x0才 能保证收敛的求根方法为局部收敛方法,

没有这个取初值限制的方法称为全局收敛方法。二分法是全局收敛的,简单迭代法一般是局部收敛的。

定理2.2 设x是迭代函数?(x)的不动点,且??(x)在点x处连续,则 ① 当???(x*)???时,迭代格式xk?1??(xk)局部收敛; ②当???(x*)???时,迭代格式xk?1??(xk)发散 证明 只给出①的证明。

由???(x*)???和??(x)在点x处连续性,存在一个正实数L<1和x的某个闭邻域

******D?:x?x*???,??0,使x?D时有??'(x)??L?1成立。

**注意到,当x?D时,由x??(x)及中值定理有

??(x)?x*????(x)??(x*?????'(?)??x?x*??L?x?x*???x?x*???

21

所以x?D时,有?(x)?D,由推论2.1可知迭代xk?1??(xk)产生的数列{xk}对任意

x0?D都收敛于不动点x*,故迭代格式xk?1??(xk)局部收敛。

利用定理2.2可以得到一个更容易使用的判别简单迭代法收敛的充分条件。

推论2.2设[a,b] 是包含迭代函数?(x)不动点x的区间,??(x)在其上连续且有

*???(x)??L???x?[a,b],则有对任意初值x0?[a,b],由迭代公式xk?1??(xk)产生的数列{xk}都收敛于x*。

4. 简单迭代法的误差估计和收敛速度

定理2.3 设定理2.1的条件成立,则有如下误差估计式

LLk*xk?xk?1 ② x?xk?x1?x0 ① x?xk?1?L1?L*证明 由迭代公式和定理2.1的条件,有

?xk?1?xk????xk??xk=??(xk)??(x*)??(x*)?xk?

**** =x?xk??(xk)??(x)?x?xk??(xk)??(x)

???x*?xL?x?k?k因为0?L?1,所以有 x?xk?另一方面

*** ???x??L??)?x?kx1xk?1?xk (2.6)

1?L?1k

?xk?1?x???(xk)??(x?)?L?xk??k1?k代入式(2.6),得结论①。

为证结论②,反复应用式(2.7)结果,有

?xk?1?x?Lx?k??k代入式(2.6),可以得到结论②。

?k1?x (2.7)

xx????k?L1?x?0定理2.3给出了收敛迭代数列{xk}的误差估计公式,利用它,在给定精度??0后,要

LLkxk?xk?1??或?x1?x0???即可。由定理2.3的②使x?xk??,只要计算到

1?L1?L*可以直接求出满足要求的根迭代次数k

k?ln???L??x1?x0/lnL (2.8)

由定理2.3的①式表示可以用刚算出的数列来估计误差。由这个误差限控制迭代次数,可用较少的迭代运算得到满足精度的近似根,特别当L?

1*时,有x?xk?xk?xk?1,此2

时可用更简单的不等式xk?xk?1??成立与否终止迭代。

由于xk?xk?1??具有简易处理且能较好描述收敛点列的特点,实用中一般不管是否

L?1成立,都用xk?xk?1??来作为终止条件,这样做通常也能求出满足精度要求的根。 2

32例2.4用简单迭代法求方程x?2x?10x?20?0在x?1附近的根,计算结果准确

22

到4位有效数字。

解 令f(x)?x3?2x2?10x?20

因为f(1)f(2)?0,在x?[1,2]时,f'(x)?3x2?4x?10?0,故原方程f(x)?0在

[1,2]内有唯一的根。将原方程改写为 x?20。

x2?2x?10为观察?(x)在[1,2]内的取值,注意到

20,得迭代函数 2x?2x?10?(x)? ?'(x)??20(2x?2) ]?0, x?[1,222(x?2x?10)20101020??(x)??2。 ,?(2)?,有1?139913故?(x)在[1,2]单调下降,由?(1)?于是有当x?[1,2]时,?(x)?[1,2]。此外,易得当x?[1,2]时

20(2x?2)20(2?2?2)120???L?1 2222(x?2x?10)(1?2?1?10)16920由推论2.1,可得迭代格式 xk?1?2 对任取x0?[1,2]产生的数列{xk}都收敛。

xk?2xk?10??'(x)??由x*?[1,2],得计算结果准确到4位有效数字时应有误差限????????。取x0?1进行迭代计算,有

??x1?1.538461538 ?x1?x0??0.538??0.5?10?2 x2?1.295019517 ?x2?x1??0.243442021?0.5?10?2

? ?

x10?1.36869397 ?x10?x9??0.365842?10?3?0.5?10?3

所以取 x*?x10?1.36869397。

注意到本题准确根为 1.3688081078214…,显然以上算出的近似根满足要求。 如果用式(2.8)估计求得迭代次数k,因为k要满足

???L??k?lnx1?x0ln?7/26000?49?13?10?3/lnL?ln/lnL??24.006?

169?7?2ln?120/169?结果给出要迭代25次才行,估计明显偏大。 例2.5建立一个迭代公式计算 g?敛性并求出g的值。

2?2?2??,要求分析所建迭代公式的收解 令xk?2???2?2?2,有k个开方号,则有limxk?g和迭代公式

k??xk?1?2?xk,由此得它的迭代函数为?(x)?2?x, (x?0)。因为当x?0时,有

11?'(x)????1,x?0。取区间[a,b]?[0,??),由推论2.1知,?(x)?0,且 ?22?x2迭代xk?1?2?xk对任意初值x0?0都产生收敛于g的数列,于是式xk?1?2?xk就是23

计算g的迭代公式。

此外,在xk?1?2?xk中令k??,得g2?2?g,它有两个解,g1?2,g2??1,

2?2?2??。

*由于xk?0,故g?0,得g?g1?2即2?*为说明简单迭代法的收敛速度,引入下面的定理。

定理2.4 设x是迭代函数?(x)的不动点,m为正整数,且?(m)?x?在x的邻域N(x*)内连续,并有关系:

??(x*)?0,???(x*)?0,?,?(m?1)(x*)?0,?(m)(x*)?0,则有

xk?1??(xk)产生的数列{xk}在N(x*)上是m阶收敛的,且有

*x*?xk?1?(m)(x) lim* (2.9) ?k??(x?x)mm!k***证明 因为?'(x)?0??'(x)?0?1,有由定理2.2可知:存在x的某个闭邻域

D?:x?x*???,??0,D?N(x*),迭代xk?1??(xk)产生的数列{xk}对任意x0?D都收敛

于不动点x。将?(x)在D内写成如下Taylor公式,并由假设有:

*?(x)??(x)??'(x)(x?x)??? ?x?特别有

?(xk)?x?*****?(m?1)(x*)(m?1)!(x?x)*m?1??(m)(?)m!

(x?x*)m

?(m)(?)m!(x?x*)m (?在x*与x之间)

?(m)(?k)m!由迭代关系xk?1??(xk)并整理上式,可得

(xk?x*)m (?k在x*与xk之间)

xk?1?x*?xk?x*?m??(m)(?k)m! (?k在x与xk之间) (2.10)

*上式令k??则得极限关系式(2.9)。

?xk?1?x*???(m)(?k)???(m)(x*)?由(2.10)式有 lim?lim??0,得满足定理条件的迭*mk??k??m!m!?xk?x?*代数列{xk}在N(x)内是m阶收敛的。

由于一般情况下?'(x*)?0,故简单迭代法通常是线性收敛的,但对某个具体的迭代函数可以达到超线性或更高阶收敛。

aa2例2.6 确定常数p,q,r使迭代函数f?x??px?q2?r5构造的迭代数列局部收

xx敛到??3a,a?0,并有尽可能高的收敛阶(要求给出这个阶数)。

解 本题要确定3个未知数,故应利用3个条件以获得3个方程才行。

aa2由题意有迭代格式为xk?1?pxk?q2?r5,它要局部收敛到?,则?就是其不动

xkxk点,即有关系??p??qa?2?ra2?5,将??3a代入并整理有:p?q?r?1。

此外要想有尽可能高的收敛阶,由定理2.4还应该有f?(?)?0,f??(?)?0,?。因为还差两个条件,故选取f?(?)?0,f??(?)?0即可,这样有

24

f?????p?2qf??????6q p?q?r?1 求解之,得p?q?a?3?5raa2?627?0

a?4?30r??0化简后得另外两个方程p?2q?5r?0,q?5r?0,这样得到关于p,q,r的方程组

p?2q?5r?0q?5r?0

51,r??,于是有具体迭代格式 995xk5aa2xk?1??2?5 (2.11)

99xk9xk因为f???(?)?10??2?0,故由(2.11)给出的迭代收敛阶最高为3。

5. 迭代法的加速

1)增量松弛法

迭代公式xk?1??(xk)产生的数列{xk}有时收敛的很慢或发散,利用定理2.4可以对这样的迭代做加速处理。实际上,如果迭代公式xk?1??(xk)产生的数列{xk}收敛的不好,可以考虑将增量?xk?xk?1?xk???xk??xk加在迭代值?(xk)上进行修正,以获得更快的收敛结果。进一步,引入带有变化因子C的增量C?xk代替?xk可以方便增量值的调控,此时经过修正的值为xk?1??(xk)?C(?(xk)?xk),它的迭代函数为

?(x)??(x)?C(?(x)?x)

*式中C??1,是待定常数。要得到收敛更快的迭代函数,当然选择C使?'(x)?0最理想。

?'(x*)*注意到?'(x)??'(x)?C(?'(x)?1)?0,可得 C?。但由于x未知,C不*1??'(x)*可计算。但由xk?x*,考虑用xk代替x,并令wk??'(xk),得到C的近似关系:

?'(xk)wk C? ?1??'x(k)?1wk将其代入新的不动点方程x??(x),可得加速的迭代格式:

*** xk?1??(xk)?wkw1?(x(k?)xk?)?xk(?)kxk (2.12)

1?wk1?wk1?wk实际计算表明这个迭代公式确实具有加速的效果,而且有时能将发散的迭代数列经过上述

加速处理后变为收敛的迭代数列。

数值计算方法中,常用近似值代替准确值来处理公式中的未知项,这种手法表面看起来很粗糙,但所得结果常常是很有用的,有时会得到很好的算法。

2) Aitken方法

Aitken方法在实用中经常使用,该方法不要求迭代函数??(x)存在。具体的构造过程为 假设求得xk后,在求xk?1前先用迭代格式xk?1???xk?校正2次:xk?1??(xk)和

?1?***?2??1??x。由微分中值定理有xk??(x)k?1?x???xk????x???????xk?x?和 ?1k?1?1? 25

?2??1??1?***?xk。 ?1?x??xk?1???x??????xk?1?x????若??(?)???(?)?L,由这二个方程,有

?1??2??1?****xk ?1?x?L?xk?x?和xk?1?x?Lxk?1?x??消去L并解出x*,有

x?xk?*??2(xk?1?xk)1xk?1?2xk?1?xk?2??1?

令该式右端为xk?1,得到Aitken迭代格式

xk?1?xk?

Aitken算法

1.给定初值x0;

2.对k=0,1,2, ?,

??2(xk?1?xk)1xk?1?2xk?1?xk?2??1? (2.13)

1) 计算xk?1??(xk),xk?1??(xk?1) 2) 计算 xk?1?xk??2???2(xk?1?xk)1?1??2??1?xk?1?2xk?1?xk3?1?

上面两个方法都是加速收敛的方法,它们的加速效果也十分明显,有时不收敛的迭代

经过它们加速后可以变为收敛了。

例2.7 用加速方法求方程x?x?1?0在1.5附近的根。

3 解 将方程改写为x?x?1?0,获得迭代格式xk?1?xk?1,显然这个迭代是发散的。

3采用Aitken公式加速,其迭代格式为: x0?1.5;xk?1?x?1,xk?1?xk?1 xk?1?xk??2???2(xk?1?xk)1?1?3k?2?????13; ?1xk?1?2xk?1?xk计算得出x4?1.32480,x5?1.32472,x6?1.32472,故所求根为1.32472,加速效果很好。

简单迭代法的优点是理论丰富、算法简单、易于推广,缺点是不易找到收敛最快的迭代函数和只是局部收敛。简单迭代法主要用于迭代的理论分析上。

迭代格式xk?1??(xk)常称为单点迭代格式,因为它需要一个初值x0就能计算出迭代数列{xk}。将其推广成xk?1??(xk,xk?1,?,xk?l),l?1为整数,则得到多点迭代,这个迭代需要多于一个1个初值x1,x2,?,xl才能构造迭代数列{xk}。需要m个初值的多点迭代称为m点迭代,其在科学计算中也常出现。

?1?,k=0,1,2,…

2.5 Newton迭代法

Newton迭代法又称切线法。该方法是Newton在研究求根时问题时借助大胆近似技术得出的目前求根方法中收敛最快的方法。

基本思想

将函数f(x)做线性化处理〔即将f(x)在某点展开为Taylor级数后,取其线性部分

L(x)代替f(x)〕,把方程f(x)?0转化为对应的近似方程L(x)?0,再从L(x)?0中构

26

造迭代公式。

1. 构造原理

迭代点列的收敛是一个逐步逼近根的过程,于是可以尝试用近似手段构造迭代公式。Newton迭代法在这方面做得很成功,其具体过程为

① 先将f(x)在xk处做Taylor展开

f''x()x'k(x)?(xk?)kx?x(k2??)

2!② 取f(x)展式中关于x?xk的线性部分即L?x??f(xk)?f'(xk)(x?xk)代替(k?)f f(x)?fxf(x),将方程f(x)?0替换为近似方程f(x)?L?x??0,有

f(xx)(?xk)?f'(kkx?) 0 (2.10)

③ 从方程(2.10)中求出解,记为xk?1,得迭代关系式

f(xk) (2.11) f'(xk)迭代式(2.11)称为Newton迭代公式,用此公式求方程f(x)?0根的方法称为Newton

xk?1?xk?迭代法。

2. 分析

Newton迭代公式虽然是从近似方程(2.10)推导出的一种求根方法,但它本质上是一种简单迭代法。实际上,由式(2.11)可得出Newton迭代法的不动点方程

f(x)??(x) (2.12) f'(x)只要f'(x)?0,方程(2.12)可由f(x)?0等价变形得出。因此有关简单迭代法的结果可

x?x?以用于Newton迭代法中。因为?(x*)?0,可知Newton迭代法具有超线性收敛速度,但实际上它具有平方收敛速度。

***定理5 设f(x)?0,f'(x)?0,且f(x)在x的领域N(x)内有二阶连续导数,则Newton迭代格式

* xk?1?xk?至少是平方收敛的。

f(xk) f'(xk)*证明 由条件知f(x)在N(x)可以成立如下Taylor公式:

1f''(?k)(x?xk)2 (2.13) 2! xk?N(x*),?k在x与xk之间 f(x)?f(xk)?f'(xk)(x?xk)?*将x?x代入式(2.13),有

0?f(x*)?f(xk)?f'(xk)(x*?xk)?1f''(?)(x*?xk)2 2设f'(xk)?0,用f'(xk)同除上式两端,并利用Newton迭代公式,整理后可得 xx?1?x?所以有

27

*f''?() (xk?x*)22f'(xk)?xk?1?x*?f''(?)f''(x*) lim?lim?**2k??k??2f'(x)2f'(x)?xk?x?k说明Newton迭代法在f'(x*)?0时至少是平方收敛的。

Newton迭代法不仅给出了一种构造求非线性方程根迭代公式的方法,而且是采用近似 处理获得满意结果的典范,这种近似处理不仅没有背离原问题,而且获得了高阶的收敛方法。Newton迭代法使用的线性化技术是目前求非线性问题最常用的方法之一。

Newton迭代法具有很形象的几何意义。注意到直线y?f(xk)?f'(xk)(x?xk)是可知Newton迭代法是由过点(xk,f(xk))的切线与x轴的交点得f(x)在xk处的切线方程,

到下一个迭代点xk?1的,这也是称Newton迭代法为切线法的原因。

应该注意的是Newton迭代法的平方敛速在f'(x*)?0的条件下才有的,当f'(x*)?0时,Newton迭代法的收敛速度只能是线性的(留作习题)。

由于Newton迭代法是局部收敛的,故初值x0应充分靠近根x才能保证收敛,这在一般情况下不容易做到。实用中可先用二分法做求根预处理,二分若干次后得到较靠近根x的近似根x0,再用此根作为Newton迭代法的初值求根,这可以达到取长补短的作用。不过也有关于Newton迭代全局收敛的结论,这里不加证明地给出相应的定理。

定理6 设f(x)?C2[a,b],满足下列3个条件 ① f(a)f(b)?0; ② f'(x)?0,x?[a,b];

③ f''(x),x?[a,b]存在且不变号。

则在[a,b]内任取一点x0,只要f(x0)f''(x0)?0,由Newton迭代格式产生的数列{xk}一定收敛于[a,b]上的唯一根x。

例2.8 用Newton迭代法求例2.4中非线性方程在x?1附近的根,计算到

***?xk?1?xk??10?6。

解 令 f(x)?x?2x?10x?20

考虑区间[a,b]?[0,2],由f(0)f(2)?0及在[0,2]内有 f'(x) ?3x?4x?10?,0故f(x)在[0,2]内有唯一根,且可以用Newton迭代法计算。本题的Newton迭代计算格式为

32xk?2xk?10xk?20 xk?1?xk? 2xk?4xk?10取x0?0计算有

322x1?1.411764706 ?x1?x0??0.411764706?0.5?10?6 x2?1.369336471 ?x2?x1??0.42428235?10?1?0.5?10?6 x3?1.368808189 ?x3?x2??0.528282?10?3 x4?1.368808108 ?x4?x3??0.81?10?7?0.5?10?6

所以x*?x4?1.368808108即为所求。

把此题与例2.4相比可知Newton迭代法收敛速度是很快的。 例2.9试构造一个能求3的迭代公式。

解 令3?x,两边平方得3?x,进而得方程f(x)?x?3?0,此时3就是它的

28

22一个根。

取包含3的一个闭区间[0,2],则在此区间上有f(0)f(2)?0,f(x)?2x?0,因此3是方程f(x)?0的唯一根,为求3,建立如下的Newton迭代公式:

2xk?313 xk?1?xk??(xk?)

2xk2xk由于f''(x)?2?0,由定理6,取x0?2,此时有f(x0)f''(x0)?0,故所求迭代公式为

13 xk?1?(xk?), x0?2 (2.14)

2xk由公式(2.14)来求3,使开方运算变为可用加减乘除实现的运算,它是计算机上做开方运算的一个实际算法。

Newton迭代法的优点是收敛速度快,可以求重根,而且算法构造简单,缺点是要进行导数计算,对f(x)要求高和局部收敛。尽管如此,它仍是一种常用方法。

2.6 Newton迭代法的变形与推广

1. Newton迭代法的变形

1) 割线法

为克服Newton迭代法需要求导数f'(xk)的不便,采用f'(xk)的近似值

f(xk)?fx(k?1)

xk?xk?1代替f'(xk)的计算,则Newton迭代公式就变为如下割线法的迭代公式 xk?1?xk?xk?xk?1f(xk) (2.15)

f(xk)?f(xk?1) 这个计算公式是多点迭代公式,它需要两个初值才能产生迭代数列。可以证明割线法是超线性收敛的,收敛阶为1.618(参看参考文献【4】)。割线法的几何意义是用过两点

。 (xk,f(xk))和(xk?1,f(xk?1))的直线与x轴的交点作为下一个迭代点xk?1的(见图2.4)

2) Newton下山法

Newton迭代法太依赖于初值x0。若x0选择 不合适,会出现不收敛的情况,为解决此问题,

加入使函数绝对值单调减少的“下山”条件:

2, ?f(xk?1)???f(xk)? (k?0,1,?要改进Newton迭代公式,使它具有“下山”种性质,

可取迭代公式为

f(xk) (2.16) f'(xk)式中wk称为下山因子,它是一个可选择的参数,当wk?1时就是Newton迭代公式。公式

xk?1?xk?wk(2.16)称为Newton下山法迭代公式,计算中,下山因子的选择依赖于下山条件的成立, 进行每步迭代时都有一个试探过程。Newton下山法的具体计算过程可参考文献【2】。

Newton下山法可使迭代减少对初值x0的依赖,它的收敛速度变为线性收敛。

29

2.Newton迭代法的推广

把Newton迭代法的思想推广到求解非线性方程组中,可以得到求解非线性方程组的Newton迭代法。

n个变元的非线性方程组的向量形式为

??? F?x??0 (2.17)

?f1(x)??f1(x1,x2,?,xn)??x1?????????f(x)f(x,x,?,x)x2?2212n??????,x? 这里,F(x)?。式中fk(x1,x2,?,xn)是x1, ???????????????f(x)f(x,x,?,x)?n??n12n??xn?x2,?,xn的多元函数,且至少有一个不是x1,x2,?,xn的线性函数。

???*?非线性方程组的解是F?x??0的一组数x*?(x1*,x2*,?,xn*)T。

参照Newton迭代法构造思想,将方程组(2.17)中的每个函数fk(x1,x2,?,xn)做线性化处理,这用多元Taylor公式可以做到,如由

?(m)n??(m)?f(x)fk(x)?fk(x)??k(xj?xj(m))?? (k?1,2?,n, )?xj?1j可得近似的线性方程组

?(m)n??(m)?f(x)?f1(x)??1(xj?xj(m))?0?xjj?1???(m)n?(m)?f2(x)?f(x)?(xj?xj(m))?0?2? ? ?xjj?1????(m)??n(m)?f(x)??fn(x)(x?x(m))?0?njj??xj?1j?记yi?xi?xi(m),整理后,得到线性方程组

mm??f1x?m??f1x???f1x??m?y1?y2???yn??f1x????x1?x2?xn?m?f2x?m??f2x?m???f2x???m??y?y???y??fx12n2 ??x (2.18) ?x2?xn1?????m??m??m??fx?fx?fxnn?n?m?y?y???y??fx12nn??x?x2?xn1?????????????????????????它是关于向量y的线性方程组,若记

30

?f1x?m??m???fx?m??????(m)??fix????f,f,?,f(m)??212n??F'(x)??, F(x)?????xj???x1,x2,?,xn????n?n??fx?m??n????(m)?1??m?则式(2.18)的矩阵解为:y = -F'(x)Fx,注意到

?????????? ??????????(m)m()?y1??x1?x1??x1??x1??????(m)?(m)???y2??x2?x2x2??x2???(m)???? y??????x?x

???????????????????(m)(m)????yn??xn?xn??xn??xn??于是可得

??(m)????(m)?1???(m)x?x?[F'(x)]F(x)

??(m?1)记此解x为x,即得解非线性方程组(2.17)的Newton迭代公式

?(m?1)?m??m(?1??m(()??)? x (2.19) ?x?[F'(x)]F(x )这个公式形式上与非线性方程求根的Newton迭代法很像,只要给出初始向量组x?1??0?就

可以利用此式计算出迭代向量x,x,?取逼近非线性方程组的解。Newton迭代格式(2.9)是目前求解非线性方程组的常用方法之一。除此之外,求解非线性方程组的方法还有最速下降法、拟Newton法等等。在求解非线性方程组问题上,寻找更为有效的解法目前仍是数值分析研究的方向。

*

2.7 知识扩展阅读:不动点与压缩映射

不动点与压缩映射都与映射概念有关,所谓映射实际是函数的推广,它的定义为 定义8 设X,Y为两个非空集合,若存在一个对应规律T,对任何一个元素x?X,都有唯一的y?Y与之对应,则T称为X到Y映射,记为

T:X→Y 或 T:x→y (x?X) y称为x在T下的象,记为y=T(x)或y=Tx。

映射是用来研究集合之间的对应关系,当X=Y时,映射T称为集合X到自身的映射,此时也形象地称T为X的一个变换,高等数学中一元函数是实数集R上的一个变换。不动点是经映射作用后不变的点,即

*定义9 设D,X为两个非空集,且D?X,映射T:D?X,若X?D,满足

?2?Tx*?x*,则称x*是映射T的不动点。

不动点常与迭代点列联系在一起,迭代点列是映射T重复作用在一个点产生的象点列,即取x0?D,迭代点列为x1?Tx0,x2?Tx1,?,xn?Txn?1,?。

压缩映像还与距离空间有关,距离空间是一个特殊的集合,其上的元素有远近描述,它定义为

定义10 设S是一个非空集合,若任取x,y?S,都有一个实数d(x,y)与之对应,并满足

① d(x,y)?0,且d(x,y)?0?x?y

31

② 任取z?S,有d(x,y)?d(x,z)?d(y,z)

则称d(x,y)是x和y之间距离,而称有距离d的空间S为距离空间,记为(S,d)。

实数集R是用两数之差的绝对值表示两数之间的距离,即 d(x,y) (x,y?R) ??x? y?压缩映射的定义为:

定义11 设(S,d)是距离空间,T是S到自身的一个映射,如果存在实数??[0,1),满足任取x,y?S,都有

d(T(x),Ty(?)?)dx(y ,则称T为S上的一个压缩映射。

本章定理1中的迭代函数?(x)就是实数集R上的一个压缩映射。由压缩映射构造的迭代点列{x(k)}是一个不断逼近不动点的点列,这一点从本章定理1可看到。压缩映射和不动点都是近代数学中常用的概念,它在很多方面都有应用。

简 评 本章以怎样求非线性方程的根为主介绍了几个求根方法,这些方法的构造手法很有代表性,应该注意学习领会。对简单迭代法,定理1或定理1'只是一个充分条件,有些不满足定理1条件的迭代数列也可能收敛。此外用迭代法求根时,可以在编程上机计算中观察数列值的变化趋势来判别收敛性,这也是实用中常用的判别收敛的做法。

非线性方程组的求解方法仍在发展之中,本章只简单介绍了一个Newton迭代法,对非线性方程组求解感兴趣的同学可参看文献[10]。本章的最后介绍了几个与本章内容有关的现代数学的概念,有关现代数学的内容可参看文献[11]。

思考题

1.迭代技术可以应用到很多求解问题。迭代格式xk?1??(xk)求的是不动点(实数),把它应用在求不动函数上,可以形式的写为fk?1?x???fk?x?,你能用迭代技术求积分方程f?x??1????f?t?dt的未知函数f?x?吗?

0x 2.用Newton迭代法编程计算7的值,初值取2,结果准确到15位有效数字。请观察每次迭代值的有效数字的变化,你能发现什么规律?

3.?如果把定理2.1 的条件②“ 任取x1,x2?[a,b],存在与x1和x2无关的正常数L?1,满足 ??(x1)??(x2)??L?x1?x2? ”改为“任取x1,x2?[a,b],满足

??(x1)??(x2)???x1?x2? ”是否也有定理2.1的结论?

4.对同一个迭代法,如果选用相同的迭代函数不同的初值会出现什么情况? 5. 解非线性方程迭代法的整体收敛和局部收敛的主要区别是什么?

数值实验

1.写出Newton迭代法算法并编写通用程序程,求方程xe?1?0 在x0?0.5附近的根,要求误差不超过10-6。

2..?确定某一具体的迭代公式的收敛阶可以通过进行数值计算的方式来估计。给定非线性方程f?x??100x?x?0,已知它的一个根为x?0,假设求此根的迭代格式为

22 xk?1?0.99xk?xkx请设计一个数值方法来估计该迭代格式的收敛阶数。

32

习题二

.

1. .? 用二分法解方程2x3?5x?1?0,x??1,3?时,第一次二分后根所在区间是什么?进行两次二分后根所在区间又是什么?

2. .? 用二分法及简单迭代法求x?3x?1?0在x0=2附近的根,准确到3位有效数字。 3.

已知方程x?x?1?0在x=1.5附近有根,试用两种迭代格式来求此根,要求说

?3323明迭代的收敛性,计算结果绝对误差??10。

24. .?用数值方法求方程sinx?x?0的所有实根,要求根的绝对误差??0.005。

5..?证明迭代格式xk?1?4?sin2xk对任取初值x0?R产生的数列都是收敛的。 6. .?确定求方程3x?e?0正根的不动点迭代的收敛区间[a,b],并求出满足

2x13xk?1?xk?10?4的近似根,若要求近似根误差限??10?4,应该迭代几次?

f(zk)7. 对复变量z?x?iy的复值函数f(z)应用Newton迭代法zk?1?zk?来求复

f'(zk)根的近似值,试导出避免复数运算的对应迭代格式。

8. 使用某种方法求出本章引例中?的近似值,准确到4位有效数字。

9. 设a为已知数,使用Newton法导出求na的迭代公式,并求极限

nlimk??*a?xk?1a?xk?n?2

10. 设x为f(x)?0的单根,f''(x)连续,证明以下迭代格式至少三阶收敛

f(xk)?y?x?k?kf'(xk)???x?y?f(yk)k?1k?f'(xk)?11. 写出简单迭代法的求根算法。 12. .?证明迭代公式xk?k?0,1,2,?

xk?1?2xk?xk?3a?3x?a2k是计算a的三阶方法。

'13. 设函数f(x)的导数满足0?m?f(x)?M,x任意,证明:任取??(0,2),M迭代格式xk?1?xk??f(xk)对任意初值x0都收敛,并确定使迭代收敛最快的值?。

*14 设x是方程f(x)?0的m(?2)重根,试证明此时Newton迭代格式是线性收敛的,

而改进的Newton迭代格式xk?1?xk?m

f(xk)是平方收敛的。

f'(xk) 33