数值分析讲义.ppt

上传人:来看看 文档编号:3853034 上传时间:2019-09-30 格式:PPT 页数:394 大小:6.23MB
返回 下载 相关 举报
数值分析讲义.ppt_第1页
第1页 / 共394页
数值分析讲义.ppt_第2页
第2页 / 共394页
数值分析讲义.ppt_第3页
第3页 / 共394页
数值分析讲义.ppt_第4页
第4页 / 共394页
数值分析讲义.ppt_第5页
第5页 / 共394页
亲,该文档总共394页,到这儿已超出免费预览范围,如果喜欢就下载吧!
资源描述

《数值分析讲义.ppt》由会员分享,可在线阅读,更多相关《数值分析讲义.ppt(394页珍藏版)》请在三一文库上搜索。

1、数 值 分 析,主讲 阎家斌,a,y,x,o,y=x,y=f(x),第1章 绪 论,1 数值分析研究的对象和内容,数值分析是研究科学计算中各种数学问题求解的数值计算方法。,用计算机进行科学计算解决实际问题的过程如下:,实际问题,数学模型,数值计算方法,程序设计,计算机计算求出结果,对数学模型建立数值计算方法,并对方法进行理论分析,直到编程上机计算出结果,以及对结果的分析,这就是数值分析研究的对象和任务。,如何评价不同算法的好坏呢?,一个好的算法应具有如下特点:,一个好的算法应具有如下特点:,(1)结构简单,易于计算机实现;,(2)理论上要保证方法的收敛性和数值稳定性;,(3)计算效率高:计算速

2、度快,节省存储量;,(4)经过数值实验检验,证明行之有效。,我们在学习的过程中,要注意掌握数值方法,随着计算机的飞速发展,数值分析方法已深入到计算,的基本原理和思想,要注意方法处理的技巧及其,与计算机的结合,要重视误差分析、收敛性和稳,定性的基本理论。,物理、计算力学、计算化学、计算生物学、计算经济学等,各个领域。本课仅限介绍最常用的数学模型的最基本的数,值分析方法。,2 误差的来源和分类,误差是描述数值计算之中近似值的精确程度,在数值计算中十分重要,误差按来源可分为模型误差、观测误差、截断误差和舍入误差四种。,1.模型误差 数学模型通常是由实际问题抽象得到的,一般带有误差,这种误差称为模型误

3、差。,2.观测误差 数学模型中包含的一些物理参数通常是通过观测和实验得到的,难免带有误差,这种误差称为观测误差。,3.截断误差 求解数学模型所用的数值方法通常是一种近似方法,这种因方法产生的误差称为截断误差或方法误差。,例如,利用ln(x+1)的Taylor公式:,实际计算时只能截取有限项代数和计算,如取前5项有:,这里产生误差(记作R5),4.舍入误差 由于计算机只能对有限位数进行运算,在运算中象,在数值分析中,我们总假定数学模型是准确的,因而不考虑模型误差和观测误差,主要研究截断误差和舍入误差对计算结果的影响。,等都要按舍入原则保留有限位,这时产生的误差称为舍入误差或计算误差。,3 绝对误

4、差、相对误差和有效数字,设x是精确值x*的一个近似值,记,e=x*-x,称e为近似值x的绝对误差,简称误差。,|e|,则称为近似值x的绝对误差限,简称误差限。,精确值x* 、近似值x和误差限之间满足:,x-x*x+,通常记为,x*=x,绝对误差有时并不能很好地反映近似程度的好坏,如,x*=10, x=1 ,y*=10000, y=5,虽然y是x的5倍,但在10000内差5显然比10内差1好。,如果满足,称er为近似值x的相对误差。,记,由于x*未知,实际使用时总是将x的相对误差取为,r =/|x|称为近似值x的相对误差限。|er|r.,设x=1.24是由精确值x*经过四舍五入得到的近似值,求x

5、的绝对误差限和相对误差限。,解 由已知可得: 1.235x*1.245,所以,=0.005, r= 0.0051.240.4%,例1,一般地,凡是由精确值经过四舍五入得到的近似值,其绝对误差限等于该近似值末位的半个单位。,定义1 设数x是数x*的近似值,如果x的绝对误差限,是它的某一数位的半个单位,并且从x左起第一个非零数,字到该数位共有n位,则称这n个数字为x的有效数字,也,称用x近似x*时具有n位有效数字。,例2,已知下列近似值的绝对误差限都是0.005,问它们,具有几位有效数字?,a=12.175,b=-0.10,c=0.1,d=0.0032,解 由于0.005是小数点后第2数位的半个单

6、位,所以,a有4位有效数字1、2、1、7,b有2位有效数字1、0,c有1位有效数字1,d没有有效数字。,数x总可以写成如下形式,x=0.a1a2ak10m,x作为x*的近似值,具有n位(nk)有效数字当且仅当,由此可见,近似值的有效数字越多,其绝对误差越小。,为了使x*=,例3,问应取几位有效数字?,解 由于,x=0.a1a2ak10 ,a1=10,故取n=6,即取6位有效数字。此时x=1.41421。,注:精确值的有效数字可认为有无限多位。,的近似值的绝对误差小于10-5,,=1.4,则近似值x可写为,其中m是整数,ai是0到9中的一个数字,a10,4 数值计算中的若干原则,为了减少舍入误差

7、的影响,设计算法时应遵循如下的一些原则。,1.避免两个相近的数相减。,如果x* ,y* 的近似值分别为x,y,则z=x-y是z* =x* -y*,的近似值.此时,相对误差满足估计式,可见,当x与y很接近时,z的相对误差有可能很大。,在数值计算中,如果遇到两个相近的数相减运算,可,考虑改变一下算法以避免两数相减。,例如:,例4,求方程x2-64x+1=0的两个根,使它们至少具有四,位有效数字。( ),解 由求根公式有,对两个相近的数相减,若找不到适当方法代替,只能,在计算机上采用双倍字长计算,以提高精度。,2.防止大数“吃掉”小数,因为计算机上只能采用有限位数计算, 若参加运算的,数量级差很大,

8、在它们的加、减运算中,绝对值很小的数,往往被绝对值较大的数“吃掉”,造成计算结果失真。,例如,用八位十进制浮点计算A=26358713+0.8+0.2,按照加法浮点运算的对阶规则,应有,A=0.26358713108+0.000000008108+0.000000002108,由于采用八位数运算,于是有A=26358713,若改变计算顺序,计算0.2+0.8+26358713,则有,A=0.00000001108+0.26358713108=26358714,可见,在求和或差的过程中应采用由小到大的运算过程。,3.绝对值太小的数不宜作除数,由于除数很小,将导致商很大,有可能出现“溢出”现象.,

9、另外,设x* ,y* 的近似值分别为x,y,则z=xy是z*=x*y*,的近似值.此时,z的绝对误差满足估计式,可见,若除数太小,则可能导致商的绝对误差很大。,4.注意简化计算程序,减少计算次数,首先,若算法计算量太大,实际计算无法完成,例如用,Cramer法则求n元线性方程组Ax=b的解,需要计算n+1个n阶,行列式,而每个n阶行列式按定义D=,需要计算(n-1)n!次乘法,则Cramer法则至少需要(n2-1)n!,次乘法,当n=20时,有(202-1)20!9.71020次乘法运算。,如果用每秒钟计算1百万次乘除运算的计算机,约需要:,9.71020106,60,24,60,365,30

10、00万年,其次,即使是可行算法,则计算量越大积累的误差也越大,,因此,算法的计算量越小越好。例如计算n次多项式:,若直接逐项计算,大约需要乘法运算次数为,若将多项式改写为:,则只需n次乘法和n次加法运算。,5.选用数值稳定性好的算法,一种数值算法,如果其计算舍入误差积累是可控制的,,则称其为数值稳定的,反之称为数值不稳定的。,例如积分,利用分部积分法可得计算In的递推公式,In=1-nIn-1,n=1,2,,由于,n=0时,取I0具有四位有效数字的近似值I00.6321,递推可得:,对任何n都应有In0,但计算结果显示I80,可见,虽然I0的,近似误差不超过0.510-4,但随着计算步数的增加

11、,误差,明显增大。这说明这里的递推公式是数值不稳定的。,事实上,由于,In=1-nIn-1,和I*n=1-nI*n-1 ,n=1,2,,可得,In-I*n=-n(In-1-I*n-1)=(-1)nn!(I0-I*0),可见,随着计算步数的增加,误差迅速放大,使结果失真。,若将计算公式改写为,类似地可得,可见,近似误差Ik-I*k是可控制的,算法是数值稳定的。,例如,由于,取近似值,,递推可得:,可见,I0已精确到小数点后四位。,习题1(第10页) 1-1,1-2,1-3,1-4,第一章练习题,第2章 解线性方程组的直接法,本章讨论n元线性方程组,(2.1),的直接解法。方程组(2.1)的矩阵形

12、式为,Ax=b,其中,若矩阵A非奇异,即det(A)0,则方程组(2.1)有唯一解。,所谓直接解法是指,若不考虑计算过程中的舍入误差,,经过有限次算术运算就能求出线性方程组的精确解的方法。,但由于实际计算中舍入误差的存在,用直接解法一般也只,能求出方程组的近似解。,Cramer法则是一种不实用的直接法,下面介绍几种实,用的直接法。,1 Gauss消去法,Gauss消元法是一种规则化的加减消元法,其基本思,想是通过逐次消元计算,把一般线性方程组的求解转化为,等价的上三角形方程组的求解。,1.1 顺序Gauss消去法,为了清楚起见,先看一个简单的例子.,考虑线性方程组,消去后两个方程中的x1得,再

13、消去最后一个方程的x2得,消元结束,经过回代得解:,上述求解的消元过程可用矩阵表示为:,(A,b)=,这是Gauss消去法的计算形式,新的增广矩阵对应的线性,方程组就是上三角形方程组,可进行回代求解。,现在介绍求解线性方程组(2.1)的顺序Gauss消去法:,记,则,线性方程组(2.1)的增广矩阵为,第一步.设 ,依次用,乘矩阵的第1行加到第i行,得到矩阵:,其中,第二步.设 ,依次用,乘矩阵的第2行加到第i行,得到矩阵:,其中,如此继续消元下去,第n-1步结束后得到矩阵:,这就完成了消元过程。,对应的方程组变成:,对此方程组进行回代,就可求出方程组的解。,顺序Gauss消去法求解n元线性方程

14、组的乘除运算量是: n2-1,+(n-1)2-1,+22-1,+1+2+n,n=20时,顺序Gauss消去法只需3060次乘除法运算.,顺序Gauss消去法通常也简称为Gauss消去法.,顺序Gauss消去法中的,称为主元素.,主元素都不为零矩阵A的各阶顺序主子式都不为零.,1.2 主元Gauss消去法,(用十进制四位浮点计算):,(用Cramer法则可得精确解x1*=1.00010 ,x2*=0.99990),解 用顺序Gauss消去法, 消元得,回代得解:x2=1.00 ,x1=0.00,若将方程组改写成:,例1 解线性方程组,用顺序Gauss消去法, 消元得,回代得解:x2=1.00 ,

15、x1=1.00,为了提高计算的数值稳定性,在消元过程中采用选择主元的方法.常采用的是列主元消去法和全主元消去法.,给定线性方程组Ax=b, 记A(1)=A,b(1)=b,列主元Gauss消去法的具体过程如下:,首先在增广矩阵B(1)=(A(1),b(1)的第一列元素中,取,然后进行第一步消元得增广矩阵B(2)=(A(2),b(2).,再在矩阵B(2)=(A(2),b(2)的第二列元素中,取,然后进行第二步消元得增广矩阵B(3)=(A(3),b(3).,按此方法继续进行下去, 经过n-1步选主元和消元运算,得到增广矩阵B(n)=(A(n),b(n).则方程组A(n)x=b(n)是与原方程组等价的

16、上三角形方程组,可进行回代求解.,易证,只要|A|0,列主元Gauss消去法就可顺利进行.,采用十进制四位浮点计算,分别用顺序Gauss消去法和列主元Gauss消去法求解线性方程组:,例2,方程组具有四位有效数字的精确解为 x1*=17.46, x2*=-45.76, x3*=5.546,解 1. 用顺序Gauss消去法求解,消元过程为,回代得: x3=5.546, x2=100.0, x1=-104.0,2. 用列主元Gauss消去法求解,消元过程为,回代得: x3=5.545, x2=-45.77, x1=17.46,可见,列主元Gauss消去法是在每一步消元前,在主元所在的一列选取绝对值

17、最大的元素作为主元素.而全主元Gauss消去法是在每一步消元前,在所有元素中选取绝对值最大的元素作为主元素.但由于运算量大增,实际应用中并不经常使用.,2 直接三角分解法,2.1 Gauss消去法的矩阵表示,对矩阵,若,令,记,则有,A(2)=,记,令,若,则有,A(3)=,如此进行下去, 第n-1步得到:,A(n)=,其中,也就是:,A(n)=Ln-1A(n-1),其中,=Ln-1Ln-2A(n-2),= Ln-1Ln-2L2L1A(1),所以有:,A=A(1)= L1-1L2-1Ln-1-1A(n)=LU,而且有,其中L=L1-1L2-1Ln-1-1, U=A(n) .,L称为单位下三角矩

18、阵; U是上三角矩阵.,式 A=LU称为矩阵A的三角分解.,2.2 Doolittle分解法,设n阶方阵A的各阶顺序主子式不为零,则存在唯一单位下三角矩阵L和上三角矩阵U使A=LU .,证明 只证唯一性, 设有两种分解 A=LU,则有,=E,所以得,于是 Ax=b LUx=b,令 Ux=y 得,定理2.1,下面介绍矩阵三角分解的Doolittle分解方法,则得,对k=2,3,n,计算,设,akj=lk1u1j+lk2u2j+lkk-1uk-1j+ukj,aik=li1u1k+li2u2k+likukk,对k=2,3,n,计算,由,可得,这就是求解方程组Ax=b的Doolittle三角分解方法。

19、,利用三角分解方法解线性方程组,解 因为,所以,例3,先解,得,再解,得,解线性方程组Ax=b的Doolittle三角分解法的计算量约为1/3n3,与Gauss消去法基本相同.其优点在于求一系列同系数的线性方程组Ax=bk ,(k=1,2,m)时,可大大节省运算量.,例如,求上例中矩阵A的逆矩阵.可分别取常向量 b1=(1,0,0)T, b2=(0,1,0)T, b3=(0,0,1)T,由,所以,为了提高数值稳定性, 可考虑列主元三角分解法, 设已完成A=LU的k-1步分解计算, 矩阵分解成,设,令 rkri,相当于取,为第k步分解的主元素.,但要注意方程组的常数项也要相应变换.,例如,用列主

20、元三角分解解例3中方程组.则有,设A为对称正定矩阵, 则有唯一分解A=LU, 且ukk0.,则有 A=LDM,又因为 (LDM)T=MTDLT=LDM 所以 M=LT,=LDLT,则有,2.3 平 方 根 法,分解A=GGT称为对称正定矩阵的Cholesky分解.,Ax=b 转换为 Gy=b , GTx=y,若记G=(gij), 则有: 对k=1,2,n,实际计算时,可采用紧凑格式,_平方根法.,解三角方程 Gy=b , GTx=y 可得,解,例4 解线性方程组,平方根法是求对称正定系数线性方程组的三角分解法,对称正定矩阵的Cholesky分解的计算量和存贮量均约为一般矩阵的LU分解的一半.

21、且Cholesky分解具有数值稳定性.,追赶法是求三对角线性方程组的三角分解法.即方程,三对角矩阵A的各阶顺序主子式都不为零的一个充分条件是:,|a1|c1|0 ; |an|dn|0 ; |ai|ci|+|di| , cidi 0 ,i=2,3,n-1.,在此条件下, A=LDM=TM , 称之为矩阵A的Crout分解.,对三对角矩阵A进行Crout分解,有,2.4 追 赶 法,其中,解三角方程 Ty=b , Mx=y 可得,称之为解三对角方程组的追赶法.,解,例5 解线性方程组,当满足条件 |a1|c1|0 ; |an|dn|0 ; |ai|ci|+|di| , cidi 0 ,i=2,3,

22、n-1. 时, 追赶法是数值稳定的, 追赶法具有计算程序简单, 存贮 量少,计算量小的优点.,3 向量和矩阵的范数,3.1 向量的范数,定义2.1 设是向量空间Rn上的实值函数, 且满足条件:,(1)非负性: 对任何向量xRn ,x0 ,且x=0当 且仅当x=0,(2)齐次性: 对任何向量x Rn 和实数 , x=| |x,(3)三角不等式: 对任何向量x ,yRn x+yx+y,则称为空间Rn上的范数,x为向量x的范数.,记x=(x1,x2,xn)T, 常用的向量范数有:,向量的1-范数:x1=|x1|+|x2|+|xn|,向量的2-范数:x2=,向量的-范数:x=,例6 设向量x=(2,-

23、4,3,1)T, 求向量范数xp ,p=1,2, .,解 由定义x1=10 , x2=,x=4 .,虽然不同范数的值可能不同,但它们间存在等价关系.,定理2.2 (范数的等价性),对于 Rn 上的任何两种范数 和 ,存在正常数m,M,使得 m x x M x , xRn,常用的三种向量范数满足如下等价关系 x x1 nx , xRn,定义2.2 设向量序列,k=1, 2,向量,如果,则称向量序列x(k)收敛于向量x*, 记作,易见,3.2 矩阵的范数,定义2.3 设是以n阶方阵为变量的实值函数,且满足条件:,(1)非负性: A0 ,且A=0当且仅当A=0,(2)齐次性: A=| |A, R,(

24、3)三角不等式:A+BA+B,(4)三角不等式:ABAB,则称A为矩阵A的范数.,记A=(aij) , 常用的矩阵范数有:,矩阵的1-范数:A1,也称矩阵的列范数.,矩阵的2-范数:A2,也称为谱范数.,矩阵的-范数:A,也称为行范数.,矩阵的F-范数:AF,例7 设矩阵,求矩阵A的范数Ap ,p=1,2, ,F.,解 A1=4 , A=5 , AF,设是一种向量范数, 则定义,称之为由向量范数派生的矩阵算子范数.,矩阵的算子范数满足,AxAx, xRn,把满足上式的矩阵范数称为与向量范数相容的矩阵范数.,对于p=1,2,矩阵范数Ap是由向量范数xp派生的矩阵算子范数,所以Ap是与xp相容的矩

25、阵范数.但AF不是一种算子范数,却与x2是相容的.,设是一种算子范数, 则,矩阵的范数与矩阵的特征值之间也有密切的联系.,设是矩阵A的特征值,x是对应的特征向量,则有,Ax= x,利用向量和矩阵范数的相容性, 则得,|x=x=AxAx,于是 |A,设n阶矩阵A的n个特征值为1, 2, , n, 则称,为矩阵A的谱半径.,对矩阵的任何一种相容范数都有, (A)A,另外, 0, 一种相容范数, 使 A (A)+,任何两种矩阵范数也具有等价性 m A A M A , ARnn,矩阵序列的收敛性也定义为,4 线性方程组固有性态与误差分析,4.1 线性方程组的固有性态,考虑线性方程组,其精确解为 x*=

26、(-9800b1+9900b2 , 9900b1-10000b2)T,若把线性方程组变为,解为x=(-9800b1+9900b2-19700 , 9900b1-10000b2+19900)T,可见 x-x*=(-19700 , 19900)T,解的误差可能放大到常数项的误差的近2万倍。,若把线性方程组变为,则线性方程组可能无解.,这种由于原始数据微小变化而导致解严重失真的线性方程组称为病态方程组, 相应的系数矩阵称为病态矩阵.,设线性方程组 Ax=b,系数矩阵是精确的, 常数项有误差b, 此时记解为x+x ,则 A(x+x) =b+b,于是 Ax =b,所以 x = A-1b A-1 b,又由

27、于 b= Ax A x,因此 x b AA-1bx,即,再设b是精确的, A有误差A, 此时记解为x+x ,则 (A+A)(x+x) =b,则有 Ax +A(x+x) =0,所以 x =-A-1A(x+x) =0,于是 x A-1Ax+x,也就是,记 Cond(A)=AA-1, 称为方程组Ax=b或矩阵A的条件数.,经常使用的条件数有 Condp(A)=ApA-1p p=1,2,。,当A为对称矩阵时,有 Cond2(A)=|1|/|n| 其中1,n分别是A的按绝对值最大和最小的特征值。,例如,对前面方程组的系数矩阵A有,Cond1(A)= Cond(A)= 39601, Cond2(A)392

28、06,由于计算条件数运算量较大, 实际计算中若遇到下述情况之一,方程组就有可能是病态的:,(1) 矩阵元素间数量级差很大 ,且无一定规律;,(2) 矩阵的行列式值相对来说很小 ;,(3) 列主元消去法求解过程中出现量级很小的主元素;,(4) 数值求解过程中 , 计算解x的剩余向量r=b-Ax已经很小, 但x仍不符合要求.,4.2 预条件和迭代改善,1. 线性方程组的预条件处理,对病态方程组Ax=b ,考虑线性方程组,其中,称之为预条件方程组, 显然与原方程组等价.可逆矩阵C称为预条件矩阵.矩阵C应满足条件,(1)条件数,(2)方程组Cz=d容易求解。,比Cond(A)明显小.,对于一般的矩阵A

29、没有十分有效的方法去选择预条件矩阵. 当A是对称正定矩阵时,可取C .,2. 线性方程组的迭代改善,设已求得方程组Ax=b的近似解x(1), 计算剩余向量,r(1)=b-Ax(1),再求解余量方程组Ax=r(1), 得到解,则x(1)的迭代改善解为:,练习题,第50页 习题2,2-2(1), 2-6(1)(用三位浮点计算),2-3(1), 2-4, 2-5, 2-8 ,2-10,2-11, 2-12, 2-14, 2-16, 2-17,第3章 解线性方程组的迭代法,迭代法的基本思想是,把n元线性方程组,(3.1),或,Ax=b,改写成等价的方程组,,或x=Mx+g,迭代法是从某一取定的初始向量

30、x(0)出发,按照一个适当的迭代公式 ,逐次计算出向量x(1), x(2),使得向量序列x(k)收敛于方程组的精确解.迭代法是一类逐次近似的方法.其优点是,算法简便,程序易于实现.,由此建立方程组的迭代公式 x(k+1)=Mx(k)+g , k=0,1,2, (3.2) 其中M称为迭代矩阵。,对任意取定的初始向量x(0),由(3.2)式可逐次算出迭代向量x(k),k=1,2,如果向量序列x(k) 收敛于x*,由(3.2)式可得,x*=Mx*+g,从而x*是方程组x=Mx+g的解,也就是方程组Ax=b的解.,这种求解线性方程组的方法称为迭代法 ,若迭代序列x(k)收敛,则称迭代法收敛,否则称迭代

31、法发散.,1 Jacobi迭代法和Gauss-Seidel迭代法,Jacobi方法是由方程组(3.1)中第k个方程解出x(k),得到等价方程组:,从而得迭代公式,式(3.3)称为Jacobi迭代法,简称为J迭代法.,则J迭代法可写成 x(k+1)=Bx(k)+g k=0,1,2,可见 ,J迭代法的迭代矩阵为,若记,J法也记为,G-S迭代法也可记为,式(3.4)称为Gauss-Seidel迭代法,简称为G-S迭代法.,若在J迭代法中,充分利用新值, 则可以得到如下的迭代公式,方程组的精确解为x*=(1,1,1)T.,解 J迭代法计算公式为,例1 用J法和G-S法求解线性方程组,取初始向量x(0)

32、=(0,0,0)T,迭代可得,计算结果列表如下:,可见,迭代序列逐次收敛于方程组的解, 而切迭代7次得到精确到小数点后两位的近似解.,G-S迭代法的计算公式为:,同样取初始向量x(0)=(0,0,0)T, 计算结果为,由计算结果可见,G-S迭代法收敛较快.取精确到小数点后两位的近似解,G-S迭代法只需迭代3次,而J迭代法需要迭代7次.,为了进一步研究,从矩阵角度来讨论上述迭代法.,对线性方程组Ax=b,记,D=diag(a11,a22,ann),则有 A=D-L-U,于是线性方程组 Ax=b 可写成 (D-L-U)x=b,等价于 Dx=(L+U)x+b 或 x=D-1(L+U)x+D-1b,由

33、此建立J迭代法迭代公式 x(k+1)=D-1(L+U)x(k)+D-1b k=0,1,2,或写成 x(k+1)=Bx(k)+g k=0,1,2,其中,G-S迭代法迭代公式可写成 x(k+1)=D-1Lx(k+1)+D-1Ux(k)+D-1b,讨论迭代法 x(k+1)=Mx(k)+g k=0,1,2,Dx(k+1)=Lx(k+1)+Ux(k)+b,(D-L)x(k+1)=Ux(k)+b,x(k+1)=(D-L)-1Ux(k)+(D-L)-1b,所以G-S迭代法可以写成 x(k+1)=Gx(k)+g k=0,1,2,其中,G=(D-L)-1U , g=(D-L)-1b,2 迭代法的收敛性,的收敛性

34、.,记误差向量e(k)=x(k)-x*,则迭代法收敛就是e(k)0.,由于,x(k+1)=Mx(k)+g k=0,1,2,x*=Mx*+g k=0,1,2,所以,e(k+1)=Me(k) , k=0,1,2,递推可得,e(k)=Mke(0) , k=0,1,2,可见,当k时, e(k)0 Mk O.,对任意初始向量x(0),迭代法收敛(M)1.,定理3.1,证 若Mk0, 则k(M)=(Mk)Mk0, 所以(M)1.,若(M)0,使得(M)+1.则 MkMk (M)+)k 0.,若M1,则对任意x(0),迭代法收敛,而且,定理3.2,证 由于 x(k+1)=Mx(k)+g x(k)=Mx(k-

35、1)+g x*=Mx*+g,所以 x(k+1) -x(k)=M(x (k) -x(k-1) ) , x(k+1) x*=M(x (k) x* ),于是有 x(k+1) -x(k) Mx (k) -x(k-1) x(k+1) x*Mx (k) x*,x(k+1) -x(k) =(x (k+1) x* )-(x(k) x* ) x (k) x*-x(k+1) x*,x(k+1) -x(k) =(x (k+1) x* )-(x(k) x* ) x (k) x*-x(k+1) x*,(1-M)x(k) x*,所以,定理3.2只是收敛的充分条件,并不必要,如,则M1=1.2,M=1.3,M2=1.09,

36、MF=1.17,但(M)=0.81,所以迭代法是收敛的.,由(3.5)式可见,M越小收敛越快,且当x (k) -x(k-1) 很小时,x(k) x*就很小,实际上用x (k) -x(k-1) 作为,迭代终止的条件.例如,x (6) -x(5) =0.011339,x(7) x(6)=0.0056695,由(3.6)式可得:,若使x(k) x* ,只需,可以事先估计达到某一精度需要迭代多少步., 即,用J迭代法求例1中方程组的解,取x(0)=(0,0,0)T,若使误差x(k)-x*10-5,问需要迭代多少次?,解 由例1知,x(1)=(1.4,0.5,1.4)T,于是有,x(1)-x(0)=1.

37、4 ,B=0.5 .,例2,k应满足,故取k=19, 即需要迭代19次.,3 J迭代法和G-S迭代法的收敛性,定理3.3 J迭代法收敛(B)1;若B1J迭代法收敛; G-S迭代法收敛(G)1;若G1 G-S迭代法收敛;,定义3.1 若n阶矩阵A=(aij)满足:,则称矩阵A是严格对角占优矩阵.,引理 若A是严格对角占优矩阵, 则det(A)0.,证 A=D-L-U=D(E-D-1(L+U)=D(E-B),因此, (B)B1,故=1不是B的特征值,det(E-B)0.,定理3.4 设A是严格对角占优矩阵,则解线性方程组Ax=b的J迭代法和G-S迭代法均收敛.,因为A是严格对角占优矩阵, 所以de

38、t(D)0, 而且,所以, det(A)0.,证 由于B1, 所以J迭代法收敛.,设是G的任一特征值, 则满足特征方程,det(E-G)= det(E-(D-L)-1U) = det(D-L)-1)det(D-L)-U)=0,所以有 det(D-L)-U)=0,若|1, 则矩阵(D-L)-U 是严格对角占优矩阵, 这与 det(D-L)-U)=0矛盾, 所以|1,于是(G)1.,定理3.5,设A 是对称正定矩阵, 则解方程组Ax=b的 (1)J迭代法收敛2D-A也正定; (2)G-S迭代法必收敛.,试建立一个收敛的迭代格式,并说明收敛性.,解 按如下方法建立迭代格式,例3 已知解线性方程组,由

39、于迭代矩阵的行范数小于1,故此迭代法收敛.,改写成,将Jacobi迭代法,4 逐次超松弛迭代法-SOR方法,写成向量形式就是 x(k+1)=x(k)+D-1(b-Ax(k) , k=0,1,2,Gauss-Seidel迭代法也可写成,或写成向量形式 x(k+1)=x(k)+D-1(b+Lx(k+1)+(U-D)x(k) , k=0,1,2,构造迭代公式,此迭代法称为SOR方法,其中参数称为松弛因子,当1时称为超松弛迭代, 当1时称为欠松弛迭代.,其矩阵形式 x(k+1)=x(k)+D-1(b+Lx(k+1)+(U-D)x(k) , k=0,1,2,于是有 Dx(k+1)=Dx(k)+(b+Lx

40、(k+1)+(U-D)x(k),所以 x(k+1)=(D-L)-1 (1-)D+Ux(k)+(D-L)-1b ,k=0,1,2,因此,SOR方法的迭代矩阵为 =(D-L)-1 (1-)D+U,SOR方法收敛()1;若1,则SOR方法收敛.,定理3.7 若SOR方法收敛, 则02.,定理3.6,证 设SOR方法收敛, 则()1,所以 |det()| =|12 n|1,而 det() =det(D-L)-1 (1-)D+U),=det(E-D-1L)-1 det(1-)E+D-1U),=(1-)n,于是 |1-|1, 或 02,定理3.8,设A是严格对角占优矩阵,则解方程组Ax=b的SOR方法,当

41、01时收敛.,定理3.9 设A是对称正定矩阵, 则解方程组Ax=b的SOR方法,当02时收敛.,证 设是的任一特征值, y是对应的特征向量, 则,(1-)D+Uy= (D-L)y,于是 (1-)(Dy,y)+(Uy,y)=(Dy,y)-(Ly,y),由于A=D-L-U是对称正定的, 所以D是正定矩阵, 且L=UT. 若记(Ly,y)=+i, 则有,(Dy,y)=0,(Uy,y)=(y,Ly)=(Ly,y),=-i,0(Ay,y)=(Dy,y)-(Ly,y)-(Uy,y) =-2,所以,当02时,有,(-+)2-(-)2= (2-)(2-) = (2-)(2-)0,所以|21, 因此()1,即S

42、0R方法收敛.,可得 =2/,设是B的任一特征值, y是对应的特征向量, 则,(L+U)y=Dy,于是 (Ly,y)+(Uy,y)=(Dy,y),当A对称正定时,即2-0,而 (2D-A)y,y)=(Dy,y)+(Ly,y)+(Uy,y) =+2,即,当A对称正定时,Jacobi迭代法收敛2D-A正定.,SOR方法收敛的快慢与松弛因子的选择有密切关系.但是如何选取最佳松弛因子,即选取=*,使()达到最小,是一个尚未很好解决的问题.实际上可采用试算的方法来确定较好的松弛因子.经验上可取1.41.6.,用SOR方法解线性方程组,解 SOR方法迭代公式为,方程组的精确解是x*=(2,1,-1)T.,

43、例4,取x(0)=(0,0,0)T,=1.46,计算结果如下:,从结果可见 ,迭代20次时已获得精确到小数点后五位的近似解.如果取=1.25,则需要迭代56次才能得到具有同样精度的近似解;如果取=1,则需迭代110次以上.,练习题,第75页 习题3 3-2, 3-3, 3-4, 3-7, 3-8, 3-9,第4章 解非线性方程的迭代法,本章讨论求非线性方程 (x)=0 (4.1) 的根的问题.,其中(x)是高次多项式函数或超越函数.如 (x)=3x5-2x4+8x2-7x+1 (x)=e2x+1-xln(sinx)-2 等等.,1 二 分 法,设(x)在区间a,b上连续且(a)(b)0,根据连

44、续函数的介值定理,区间a,b上必有方程(x)=0的根,称a,b为方程(x)=0的有根区间.,得到新的有根区间a1,b1,设(x)在区间a,b上连续且(a)(b)0 .,0,a,b,y,x,y=(x),记a0=a,b0=b,计算,若|(x0)|0,取a1=x0,b1=b0,而且有根区间a1,b1长度是有根区间a0,b0长度的一半,x0,再对有根区间a1,b1重复上面运算, 即: 计算,若|(x1)|0, 取a2=x1 ,b2=b1,得到新的有根区间a2,b2.,x1,而且有根区间a2,b2长度是有根区间a1,b1长度的一半.一直进行下去,直到求出有根区间ak,bk.,此时,再计算,或者有|(xk

45、)| ,或者有,可见,k趋向无穷大时,xk收敛于 .,而且,若要|xk-| ,只要,此时可取近似根xk .,在计算过程中,若出现|(xk)|1,或bk-ak2 .则可取xk作为方程(x)=0的近似根,终止运算.,例1,用二分法求x3+4x-10=0在区间1,2内根的近似值,并估计误差.,解 这里(x)=x3+4x-7, (1)(2)=-180,所以(x)=0在1,2区间有唯一根.,取x0=1.5,由于(x0)=2.375,得新有根区间1,1.5,x1=1.25,由于(x1)=-0.0468,得新有根区间1.25,1.5,x2=1.375,由于(x2)=1.0996,得新有根区间1.25,1.3

46、75,x3=1.3125,由于(x3)=0.511,得新有根区间1.25,1.3125,.,x9=1.254882813,得有根区间1.254882813,1.255859375,x10=1.255371094, (x10)=-0.000105285,取x10=1.255371094作为方程根的近似值,且有,只需k5ln210-115.61.即需取x16.,如果取精度=10-5,则要使,二分法要求函数在区间a,b上连续,且在区间两端点函数值符号相反,二分法运算简便、可靠、易于在计算机上实现。但是,若方程(x)=0在区间a,b上根多于1个时,也只能求出其中的一个根。另外,若方程(x)=0在区间a,b有重根时,也未必满足(a)(b)0. 而且由于二分法收敛的速度不是很快,一般不单独使用,而多用于为其他方法提供一个比较好的初始近似值.,2.1 简单

展开阅读全文
相关资源
猜你喜欢
相关搜索

当前位置:首页 > 其他


经营许可证编号:宁ICP备18001539号-1