解线性方程组的直接方法.ppt

《解线性方程组的直接方法.ppt》由会员分享，可在线阅读，更多相关《解线性方程组的直接方法.ppt（61页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、第第2 2章章 解线性方程组的直接法解线性方程组的直接法考虑一个简单的例。保持人体健康需要多种微量元素，主要通过食物摄取，并且不同的年龄，职业和身体条件要求的量也不同。每千克食品含量 菠菜 萝卜 苹果要求值 铁 a11 a12 a13 b1 锌 a21 a22 a23 b2 钙 a31 a32 a33 b3考虑 n元线性方程组(2.1)方程组(2.1)的矩阵形式为 A Ax=b=b其中 若矩阵A A非奇异,即det(A A)0,则方程组(2.1)有唯一解。所谓直接解法是指，若不考虑计算过程中的舍入误差，经过有限次算术运算就能求出线性方程组精确解的方法。但由于实际计算中舍入误差的存在，用直接解法

2、一般也只能求出方程组的近似解。在线性代数中我们学过用Cramer法则求解线性方程组。行列式,而每个n阶行列式按定义D=需要计算(n-1)n!次乘法,则Cramer法则至少需要(n2-1)n!次乘法,当n=20时,有(202-1)20!9.71020次乘法运算。如果用每秒钟计算1百万次乘除运算的计算机，约需要：6024603653000万年Cramer法则求n元线性方程组Ax=b的解,需要计算n+1个n阶9.71020106可见Cramer法则不是一种实用的直接法，下面介绍几种实用的直接法。1 Gauss1 Gauss消去法消去法 Gauss Gauss消元法是一种规则化的加减消元法，其基本思想

3、是通过逐次消元计算，把一般线性方程组的求解转化为等价的上三角形方程组的求解。对角上三角下三角1.1 1.1 顺序顺序GaussGauss消去法消去法为了清楚起见,先看一个简单的例子.考虑线性方程组单位上三角单位下三角上三角方程组 Ax=b:Ax=b:可按xn,xn-1,x1，顺序进行回代求解消去后两个方程中的x1得再消去最后一个方程的x2得消元结束,经过回代得解:上述求解的消元过程可用矩阵表示为:(A,bA,b)=这是GaussGauss消去法的矩阵计算形式,新的增广矩阵对应的线性方程组就是上三角形方程组,可进行回代求解。现在介绍求解线性方程组 Ax=b Ax=b 的顺序GaussGauss消

4、去法:记则,线性方程组 Ax=bAx=b 的增广矩阵为 第一步.设 ,依次用乘矩阵的第1行加到第i行(i=2,3,n),得到矩阵:其中 第二步.设 ,依次用乘矩阵(A(A(2)(2),b,b(2)(2)的第2行加到第i行,得到矩阵:其中如此继续消元下去,第n-1步结束后得到矩阵:这就完成了消元过程。对应的方程组变成：对此方程组进行回代，就可求出方程组的解。顺序Gauss消去法求解n元线性方程组的乘除运算量是：n=20时,顺序Gauss消去法只需3060次乘除法运算.顺序Gauss消去法通常也简称为GaussGauss消去法消去法.顺序Gauss消去法中的 称为主元素主元素.主元素都不为零矩阵A

5、 A的各阶顺序主子式都不为零.例例 用高斯消去法求解Ax=b,其中解回代求解：x3=3,x2=2,x1=2 1.2 1.2 主元主元GaussGauss消去法消去法 (用十进制四位浮点计算):此问题的精确解为x1*=1.0001,x2*=0.9999。解解 用顺序Gauss消去法,消元得 回代得解:x2=1.00,x1=0.00 若将方程组改写成:例例1 解线性方程组 再用Gauss消去法,消元得 回代得解:x2=1.00,x1=1.00 为了提高计算的数值稳定性,在消元过程中采用选择主元的方法.常采用的是列主元消去法和全主元消去法列主元消去法和全主元消去法.给定线性方程组AxAx=b b,记

6、A A(1)=A A,b b(1)=b b,列主元列主元GaussGauss消去法消去法的具体过程如下:首先在增广矩阵B B(1)=(A A(1),b b(1)的第一列元素中,取 然后进行第一步消元得增广矩阵B B(2)=(A A(2),b b(2).再在矩阵B B(2)=(A A(2),b b(2)的第二列元素中,取 然后进行第二步消元得增广矩阵B(3)=(A(3),b(3).按此方法继续进行下去,经过n-1步选主元和消元运算,得到增广矩阵B(n)=(A(n),b(n).则方程组A(n)x=b(n)是与原方程组等价的上三角形方程组,可进行回代求解.易证,只要|A|0,列主元Gauss消去法就

7、可顺利进行.采用十进制四位浮点计算,分别用顺序Gauss消去法和列主元Gauss消去法求解线性方程组:例例2 方程组具有四位有效数字的精确解为 x1*=17.46,x2*=-45.76,x3*=5.546 解解 1.1.用顺序Gauss消去法求解,消元过程为回代得:x3=5.546,x2=100.0,x1=-104.0 2.2.用列主元Gauss消去法求解,消元过程为回代得:x3=5.545,x2=-45.77,x1=17.46 可见，列主元Gauss消去法是在每一步消元前,在主元所在的一列选取绝对值最大的元素作为主元素.全主元Gauss消去法是在每一步消元前,在所有元素中选取绝对值最大的元素

8、作为主元素.但由于运算量大增,实际应用中并不经常使用.2 2 直接三角分解法直接三角分解法 2.1 Gauss2.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).A=LU,LA=LU,L为单位下三角矩阵;U U是上三角矩阵.分解式 A=LU称为矩阵A的直接三角分解直接三

9、角分解.2.2 Doolittle2.2 Doolittle分解法分解法 设n阶方阵A的各阶顺序主子式不为零,则存在唯一单位下三角矩阵L和上三角矩阵U使A=LU.证明证明 只证唯一性,设有两种分解 A=LU则有 =E所以得设 A A有分解 A=LUA=LU,于是Ax=b LUx=b,令 Ux=y,则得 定理定理2.1这是与Ax=b 等价的线性方程组 下面介绍矩阵三角分解的Doolittle分解方法,利用矩阵乘法规则依次求 U U的第k行，L L的第k列(k=1,2),对k=2,3,n,计算 设 akj=lk1u1j+lk2u2j+lkk-1uk-1j+ukj aik=li1u1k+li2u2k

10、+likukk对k=2,3,n,计算要按顺序先求U U的第1行，L L的第1列，U U的第2行，L L的第列，,U U的第k行，L L的第k列，。可导出计算公式：在完成分解 A=LUA=LU后，依次求解方程组 Ly=bLy=b和 Ux=yUx=y，既可得这就是求解方程组AxAx=b b的Doolittle三角分解方法。利用三角分解方法解线性方程组 解解 令利用矩阵乘法规则，先求 U U的第1行元素，可得：例例32=u112=u12-3=u13,再求L L的第2列。利用矩阵乘法规则建立方程得到：2=2l21,3=2l31,则有再求U U的第2行。利用矩阵乘法规则建立方程得到：-1=2+u22,3

11、=-3+u23,则有:u22=-3,u23=6,然后，再求U U的第3行，2=-9/2+10+u33,最终得 ,得先解再解 ,得 解线性方程组Ax=b的Doolittle三角分解法的计算量约为1/3n3,与Gauss消去法基本相同.其优点在于求一系列同系数的线性方程组Ax=bk,(k=1,2,m)时,可大大节省运算量.求解：Ly=bLy=b和Ux=yUx=y例如,求如下 3阶矩阵A A的逆矩阵：这等价于求X X使AX=EAX=E。设。设X=(xX=(x1 1,x,x2 2,x,x3 3)，可分别取3维单位向量：e e1=(1,0,0)T,e e2=(0,1,0)T,e e3=(0,0,1)T由

12、矩阵分块乘法：AxAx1 1=e=e1 1,Ax,Ax2 2=e=e2 2,Ax,Ax3 3=e=e3 3,利用三角分解法得:依次求解 LyLyi i=e=ei i,Ux,Uxi i=y=yi i,i=1,2,3.则 X=(xX=(x1 1,x,x2 2,x,x3 3)为所求.所以 为了提高数值稳定性,可考虑列主元三角分解法,设已完成A=LU的k-1步分解计算,矩阵分解成设 令 rkri 相当于取 为第k步分解的主元素.但要注意方程组的常数项也要相应变换.例如,用列主元三角分解解例3中方程组.则有 设A为对称正定矩阵,则有唯一分解A=LU,且ukk0.则有 A=LDM又因为 (LDM)T=MT

13、DLT=LDM 所以 M=LT =LDLT则有 2.3 平平 方方 根根 法法 分解A=GGT称为对称正定矩阵的Cholesky分解分解.若记G G=(gij),则有:对k=1,2,n 按列求矩阵 G G的元素。实际计算时,可采用紧凑格式 Ax=bGGTx=bGy=b,GTx=y,称为平方根法 Ax=bGGTx=bGy=b,GTx=y,解三角方程可得 解解 令 例例4 解线性方程组 平方根法是求对称正定系数线性方程组的三角分解法,对称正定矩阵的Cholesky分解的计算量和存贮量均约为一般矩阵的LU分解的一半.且Cholesky分解具有数值稳定性.利用矩阵乘法规则可计算得 追赶法是求三对角线性

14、方程组的三角分解法.即方程组 三对角矩阵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 2.4 追追 赶赶 法法其中Ax=bTMx=b Ty=b,Mx=y ,解三角方程可得称之为解三对角方程组的追赶法追赶法.计算量为 0(n).例例5 解线性方程组解解.令利用矩阵乘法规则依次计算得到 当满足条件|a1|c1|0;|an|dn|0;|ai|ci|+|di|,cidi 0,i=2,3,n-1.时,追赶

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

16、范数:x2=向量的向量的-范数范数:x=例例6 6 设向量x=(2,-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*,记作 易见,按|x|1范数：3.2 3.2 矩阵的范数矩阵的范数 定义定义2.32.3 设是以n阶方阵为自变量的实值函数

17、,且满足条件:(1)非负性非负性:A0,且A=0当且仅当A=0 0 (2)齐次性齐次性:A=|A,R (3)三角不等式三角不等式:A+BA+B (4)乘积不等式乘积不等式:ABAB则称A为矩阵A的范数.记A=(aij),常用的矩阵范数有:矩阵的矩阵的1-1-范数范数:A1 ,也称矩阵的列范数列范数.矩阵的矩阵的2-2-范数范数:A2 ,也称为谱范数谱范数.矩阵的矩阵的-范数范数:A ,也称为行范数行范数.矩阵的矩阵的F F-范数范数:AF 例例7 设矩阵求矩阵A的范数Ap,p=1,2,F.解解 A1=4,A=5,AF 设是一种向量范数,则定义称之为由向量范数派生的矩阵算子范数矩阵算子范数.矩阵

18、的算子范数满足 AxAx,xRn把满足上式的矩阵范数称为与与向量范数相容的矩阵范数向量范数相容的矩阵范数.对于p=1,2,矩阵范数Ap是由向量范数xp派生的矩阵算子范数,所以Ap是与xp相容的矩阵范数.但AF不是一种算子范数,却与x2是相容的.设是一种算子范数,则 矩阵的范数与矩阵的特征值之间也有密切的联系.设是矩阵A的特征值,x是对应的特征向量,则有 Ax=x利用向量和矩阵范数的相容性,则得|x=x=AxAx 于是|A 设n阶矩阵A的n个特征值为1,2,n,则称 为矩阵A的谱半径谱半径.对矩阵的任何一种相容范数都有 (A)A 另外,0,一种相容范数,使 A(A)+任何两种矩阵范数也具有等价性

19、 m A A M A,ARnn 矩阵序列的收敛性定义为4 4 线性方程组固有性态与误差分析线性方程组固有性态与误差分析4.1 4.1 线性方程组的固有性态线性方程组的固有性态考虑线性方程组其精确解为 x*=(-9800b1+9900b2 ,9900b1-10000b2)T若把线性方程组变为解为x=(-9800b1+9900b2-19700 ,9900b1-10000b2+19900)T可见 x-x*=(-19700 ,19900)T 解的误差可能放大到右端项误差的近2万倍。若把线性方程组变为因系数行列式为0，则线性方程组可能无解.这种由于原始数据微小变化而导致解严重失真的线性方程组称为病态方程

20、组病态方程组,相应的系数矩阵称为病态矩阵病态矩阵.设线性方程组 Ax=b系数矩阵是精确的,右端项有误差b,此时记解为x+x,则 A(x+x)=b+b于是 Ax+A x=b+b,Ax=b,x=A-1b所以 x=A-1b A-1 b又由于 b=Ax A x因此 x b AA-1bx即可见误差是否可控制取决于A A-1的大小。再设b是精确的,A有误差A,此时记解为x+x,则 (A+A)(x+x)=b则有 Ax+Ax+A(x+x)=b,Ax+A(x+x)=0所以 x=A-1A(x+x)于是 x A-1Ax+x也就是记 Cond(A)=AA-1,称为方程组Ax=b或矩阵A的条件数。如果 Cond(A)1

21、，则称A为病态矩阵。经常使用的条件数有 Condp(A)=ApA-1p p=1，2，。当A为对称矩阵时，有 Cond2(A)=A2A-12=|1|/|n|其中1，n分别是A的按绝对值最大和最小的特征值。例如，对前面方程组的系数矩阵A有 Cond1(A)=Cond(A)=39601，Cond2(A)39206 由于计算条件数运算量较大,实际计算中若遇到下述情况之一,方程组就有可能是病态的:(1)矩阵元素间数量级差很大,且无一定规律;(2)矩阵的行列式值相对来说很小;(3)列主元消去法求解过程中出现量级很小的主元素;(4)数值求解过程中,计算解x的剩余向量r=b-Ax已经很小,但x仍不符合要求.4

22、.2 4.2 预条件和迭代改善预条件和迭代改善1.1.线性方程组的预条件处理线性方程组的预条件处理 对病态方程组Ax=b,考虑线性方程组其中 称之为预条件方程组预条件方程组,显然与原方程组等价.(1)条件数 (2)方程组Cz=d容易求解。比Cond(A)明显小.对于一般的矩阵A没有十分有效的方法去选择预条件矩阵.当A是对称正定矩阵时,可取C .2.2.线性方程组的迭代改善线性方程组的迭代改善 设已求得方程组Ax=b的近似解x(1),计算剩余向量 r(1)=b-Ax(1)再求解余量方程组Ax=r(1),得到解 则x(1)的迭代改善解为:可逆矩阵C称为预条件矩阵预条件矩阵.矩阵C应满足条件练习题练习题第第50页页 习题习题22-2(1),2-6(1)（用二位浮点计算）（用二位浮点计算）练习题练习题第第50页页 习题习题22-3(1),2-4,2-5,2-8练习题练习题第第50页页 习题习题22-10,2-11,2-12,2-14,2-16,2-17课课 堂堂 练练 习习 1.对矩阵A=进行LU分解.2.设矩阵A=,求(A)和Cond(A).解 1.由于 ,故 A=2.由于|A-E|=2-3-10=(+2)(-5)所以:1=-2,2=5,于是(A)=5.A-1=,Cond(A)=|A|A-1|=76/10=4.2 课间休息休息