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

第5章解线性方程组的直接方法,5.1引言与预备知识5.2高斯消去法5.3矩阵三角分解法5.4向量和矩阵的范数5.5误差分析,这一章讨论线性方程组,5.1引言与预备知识,的数值解法.,在自然科学和工程技术中,很多问题归结为解线性方程组.例如电学中的网络问题，船体数学放样中建立三次样条函数问题，用最小二乘法求实验数据的曲线拟合问题，解非线性方程组问题，用差分法或者有限元方法解常微分方程、偏微分方程边值问题等都导致求解线性方程组，而这些方程组的系数矩阵大致分为两种，一种是低阶稠密矩阵(例如，阶数不超过150).另一种是大型稀疏矩阵(即矩阵阶数高且零元素较多).,5.1.1引言,有的问题的数学模型中虽不直接表现为含线性方程组,但它的数值解法中将问题“离散化”或“线性化”为线性方程组.因此线性方程组的求解是数值分析课程中最基本的内容之一.,关于线性方程组的解法一般有两大类：,但实际计算中由于舍入误差的存在和影响，这种方法也只能求得线性方程组的近似解.本章将阐述这类算法中最基本的和具有代表性的算法就是高斯消去法,以及它的某些变形和应用.这类方法是解低阶稠密矩阵方程组及某些大型稀疏矩阵方程组(例如，大型带状方程组)的有效方法.,1.直接法,经过有限次的算术运算,可以求得方程组的精确解(假定计算过程没有舍入误差).如线性代数课程中提到的克莱姆算法就是一种直接法.但该法对高阶方程组计算量太大,不是一种实用的算法.,就是用某种极限过程去逐步逼近方程组精确解的方法.迭代法具有计算机的存储单元较少、程序设计简单、原始系数矩阵在计算过程中始终不变等优点，但存在收敛条件和收敛速度问题.迭代法是解大型稀疏矩阵方程组(尤其是由微分方程离散后得到的大型方程组)的重要方法(见第6章).,为了讨论线性方程组数值解法，需复习一些基本的矩阵代数知识.,2.迭代法,5.1.2向量和矩阵,基本概念：,用Rmn表示全部mn实矩阵的向量空间；,用Cmn表示全部mn复矩阵的向量空间.,(由数排成的矩阵表，称为m行n列矩阵).,其中aj为A的第j列的m维列向量.同理,(m维列向量).,其中biT为A的第i行的n维行向量.,矩阵的基本运算：,(1)矩阵加法,(2)矩阵与标量的乘法,(3)矩阵与矩阵的乘法,(4)转置矩阵,(5)单位矩阵,其中,(6)非奇异矩阵,则称B是A的逆矩阵，记为A-1，且(A-1)T=(AT)-1.如果A-1存在，则A称为非奇异矩阵.如果A、B均为非奇异矩阵，则有(AB)-1=B-1A-1.,(7)矩阵的行列式,设ARnn，则A的行列式可按任一行(列)展开,其中Aij为aij的代数余子式，Aij=(-1)i+jMij，为Mij元素aij的余子式.,行列式性质：,5.1.3矩阵的特征值与谱半径,定义1设A=(aij)Rnn，若存在数(实数或复数)和非零向量x=(x0,x1,xn)TRn，使Ax=x,(1.1)则称为A的特征值,x为A对应的特征向量,A的全体特征值称为A的谱,记作(A),即(A)=1,2,n.记,称为A的谱半径.,由(1.1)式可知可使齐次线性方程组(I-A)x=0有非零解，故系数行列式det(I-A)=0，记,p()称为矩阵A的特征多项式，方程(1.3)称为A的特征方程.因为n次代数方程p()在复数域中有n个根1,2,n，故,由(1.3)式中的行列式展开可得,故A=(aij)Rnn的n个特征值1,2,n是它的特征方程(1.3)的n个根，并有,称trA为A的迹.,及,A的特征值和特征向量x还有以下性质：,AT与A有相同的特征值及特征向量x.,若A非奇异,则A-1的特征值为-1,特征向量为x.,相似矩阵B=S-1AS有相同的特征多项式.,例1求矩阵的特征值及谱半径,解矩阵A的特征方程为,得A的特征值为1=2=2,3=7，谱半径(A)=7.,5.1.4特殊矩阵,设A=(aij)Rnn.,(1)对角矩阵如果当ij时，aij=0.,(2)三对角矩阵如果当|i-j|1时，aij=0.,(3)上三角矩阵如果当ij时，aij=0.,(4)上海森伯格(Hessenberg)矩阵如果当ij+1时，aij=0.,(5)对称矩阵如果AT=A.,(6)埃尔米特(Hermit)矩阵设ACnn，如果当AH=A(AH是A的共轭转置矩阵).,(7)对称正定矩阵如果AT=A且对任意非零向量，xRn，(Ax,x)=xTAx0.,(8)正交矩阵如果A-1=AT.,(9)酉矩阵设ACnn，如果A-1=AH.,(10)初等置换阵由单位矩阵I交换第i行与第j行(或第i列与第j列)，得到的矩阵记为Iij，且.,IijA=B(为交换A第i行与第j行得到的矩阵)；,AIij=C(为交换A第i列与第j列得到的矩阵).,(11)置换阵由初等置换阵的乘积得到的矩阵.,定理1设ARnn，则下述命题等价：,(1)对任何bRn，方程组Ax=b有唯一解.,(2)齐次方程组Ax=0只有唯一解零解x=0.,(3)det(A)0.,(4)A-1存在.,(5)A的秩rank(A)=n.,定理2设ARnn为对称正定矩阵，则,(1)A为非奇异矩阵，且A-1亦是对称正定矩阵.,(2)设Ak为A的顺序主子阵，则Ak(k=1,2,n)亦是对称正定矩阵，其中,(3)A的特征值i0(i=1,2,n).,(4)A的顺序主子式都大于零，即Ak的行列式det(Ak)0(k=1,2,n).,定理3设ARnn为对称矩阵，如果det(Ak)0(k=1,2,n)，或A的特征值i0(i=1,2,n)，则A为对称正定矩阵.,有重特征值的矩阵不一定相似于对角矩阵，那么一般n阶矩阵A在相似变换下能简化到什么形状.,定理4(Jordan标准型)设A为n阶矩阵，则存在一个非奇异n阶矩阵P使得,为若当(Jordan)块.,其中,(1)当A的若当标准型中所有若当块Ji均为一阶时，此标准型变为对角矩阵.,(2)如果A的特征值各不相同，则其若当标准型必为对角矩阵diag(1,2,n).,5.2高斯消去法,本节介绍高斯消去法(逐次消去法)及消去法和矩阵三角分解之间的关系.虽然高斯消去法是一种古老的求解线性方程组的方法(早在公元前250年我国就掌握了解方程组的消去法)，但由它改进、变形得到的选主元素消去法、三角分解法仍然是目前计算机上常用的有效方法.我们在中学学过消去法,高斯消去法就是它的标准化的、适合在计算机上自动计算的一种方法.,5.2.1高斯消去法,设有线性方程组,或写成矩阵形式,简记为Ax=b.,解为,当m=n时，对三角形方程组的解非常容易的.,如:上三角矩阵所对应的线性方程组,下三角矩阵所对应的线性方程组,计算量（乘除法的主要部分）都为n2/2.,解为,因此，我们将一般的线性方程组化成等价的三角形方程组来求解.,首先举一个简单的例子来说明消去法的基本思想.,例2用消去法解线性方程组,解第1步，将方程(2.2)乘上-2加到方程(2.4)上去，消去(2.4)中的未知数x1，得到,第2步，将方程(2.3)加到方程(2.5)上去，消去(2.5)中的未知数x2，得到与原方程组等价的三角形方程组,显然，方程组是(2.6)是容易求解的，解为,上述过程相当于对方程的增广阵做初等行变换,其中ri用表示矩阵的第i行.,由此看出，用消去法解方程组的基本思想是用逐次消去未知数的方法把原方程组Ax=b化为与其等价的三角形方程组，而求解三角形方程组可用回代的方法求解.换句话说，上述过程就是用初等行变换将原方程组系数矩阵化为简单形式(上三角矩阵)，从而求解原方程组(2.1)的问题转化为求解简单方程组的问题.或者说，对系数矩阵A施行一些行变换(用一些简单矩阵左乘A)将其约化为上三角矩阵.这就是高斯消去法.,下面讨论求解一般线性方程组的高斯消去法.由,将(2.1)记为A(1)x=b(1)，其中,(1)第1步(k=1).,设a11(1)0，首先计算乘数mi1=ai1(1)/a11(1)(i=2,3,m).用-mi1乘(2.1)的第一个方程，加到第i个(i=2,3,m)方程上，消去(2.1)的从第二个方程到第m个方程中的未知数x1，得到与(2.1)等价的方程组,(2.7),简记为A(2)x=b(2)，,其中A(2),b(2)的元素计算公式为,(2)第k次消元(k=1,2,s,s=min(m-1,n).,设上述第1步，第k-1步消元过程计算已经完成，即已计算好与(2.1)等价的方程组，简记为A(k)x=b(k)，其中,设akk(k)0，计算乘数mik=aik(k)/akk(k)(i=k+1,m).用-mik乘(2.8)的第k个方程加到第i个(i=k+1,m)方程上，消去从第k+1个方程到第m个方程中的未知数xk，得到与(2.1)等价的方程组A(k+1)x=b(k+1)，,其中A(k+1),b(k+1)的元素计算公式为，,显然A(k+1)中从第1行到第k行与A(k)相同.,(3)继续上述过程，且设aii(i)0(i=1,2,s)，直到完成第s步消元计算.最后得到与原方程组等价的简单方程组A(s+1)x=b(s+1)，其中A(s+1)为上阶梯.,特别当m=n时，与原方程组等价的简单方程组为A(n)x=b(n)，即,由(2.1)约化为(2.10)的过程称为消元过程.,如果A是非奇异矩阵，且akk(k)0(k=1,2,n-1)，求解三角形方程组(2.10)，得到求解公式,(2.10)的求解过程(2.11)称为回代过程.,注意：设Ax=b，其中ARnn为非奇异矩阵，如果a11(1)=0，由于A为非奇异矩阵，所以A的第1列一定有元素不等于零，例如al10，于是可交换两行元素(即r1rl)，将al1调到(1,1)位置，然后进行消元计算，这时A(2)右下角矩阵为n-1阶非奇异矩阵，继续这过程，高斯消去法照样可进行计算.,总结上述讨论即有,定理5设Ax=b，其中ARnn.,(1)如果akk(k)0(k=1,2,n)，则可通过高斯消去法将Ax=b约化为等价的三角形方程组(2.10)，且计算公式为,消元计算(k=1,2,n-1),回代计算,(2)如果A为非奇异矩阵，则可通过高斯消去法(及交换两行的初等变换)将方程组Ax=b约化为等价的三角形方程组(2.10).,算法1(高斯算法)设ARmn(m1),s=min(m-1,n)，如果akk(k)0(k=1,2,s)，本算法用高斯方法将A约化为上三角形矩阵，且A(k)覆盖A，乘数mik覆盖aik.,对于k=1,2,s,(1)如果akk=0，则计算停止；,(2)对于i=k+1,m,(a)aikmik=aik/akk,(b)对于j=k+1,naijaij-mikakj.,显然，算法1第k步需要m-k次除法，(m-k)(n-k)次乘法，因此，本算法(从第1步到第s步消元计算总的计算量)大约需要s3/3-(m-k)s2/2+mns次乘法运算(对相当大的s).当m=n时,总共大约需要n3/3次乘法运算.,数akk(k)在高斯消去法中有着突出的作用，称为约化的主元素(简称主元).,算法2(回代算法)设Ux=b，其中U=(uij)Rnn为非奇异上三角矩阵，本算法计算Ux=b的解.,对于i=n,n-1,1,(1)xibi,(2)对于j=i+1,n,xixi-uijxj(in不算),(3)xixi/uii,这个算法需要n(n-1)/2次乘除法运算.,例子解方程组,解：消元,回代得,高斯消去法对于某些简单的矩阵可能会失败，例如,由此，需要对算法1进行修改，首先研究原来矩阵A在什么条件下才能保证akk(k)0(k=1,2,s).下面的定理给出了这个条件.,定理6约化的主元素aii(i)0(i=1,2,k)的充要条件是方阵A的顺序主子式Di0(i=1,2,k)，即,证明首先利用归纳法证明充分性.显然,当k=1时，结论成立，假设结论对k-1是成立的，对k由条件设Di0(i=1,2,k)，于是由归纳法假设我们有aii(i)0(i=1,2,k-1)，可用高斯消去法将A(1)约化到A(k)，即,利用行列式的性质，我们有,由条件有Di0(i=1,2,k)，利用(2.13)式，则有aii(i)0(i=1,2,k)，定理6的对k成立.,必要性，由条件aii(i)0(i=1,2,k)，利用(2.13)式亦可推出Di0(i=1,2,k).,定理6得证.,推论如果A的顺序主子式Di0(i=1,2,n-1),则,5.2.2矩阵的三角分解,下面我们借助矩阵理论进一步对消去法做些分析，从而建立高斯消去法与矩阵因式分解的关系.,设(2.1)的系数矩阵ARnn的各顺序主子式均不为零，对于A施行初等行变换相当于用初等矩阵左乘A，于是对(2.1)施行第1次消元后化为(2.7)，这时A(1)化为A(2)，b(1)化为b(2)，即,L1A(1)=A(2),L1b(1)=b(2)，,其中,一般第k步消元，A(k)化为A(k+1)，b(k)化为b(k+1)，相当于,LkA(k)=A(k+1),Lkb(k)=b(k+1)，,其中,重复以上过程，最后得到,(2.14),将上三角矩阵A(n)记为U，由(2.14)得到,其中,为单位下三角矩阵.,这就是说，高斯消去法实质上产生了一个将A分解为两个三角形矩阵相乘的因式分解，于是我们得到如下重要定理，它在解方程组的直接法中起着重要作用.,定理7(矩阵的LU分解)设A为n阶矩阵，如果A的顺序主子式Di0(i=1,2,n-1)，则A可分解为一个单位下三角矩阵L和一个上三角矩阵U的乘积，即A=LU,且这种分解是唯一的.,证明根据以上高斯消去法的矩阵分析，A=LU的存在性已经得到证明，现仅在A为非奇异矩阵的假定下来证明唯一性，当A为奇异矩阵的情况留作练习，设A有两种分解为,A=LU=L1U1，,其中L,L1为单位下三角矩阵，U,U1为上三角矩阵.,由于L-1,U1-1存在，故有L-1L1=UU1-1.,上式右边为上三角矩阵，左边为单位下三角矩阵，从而上式两边都必须等于单位矩阵I，故有L=L1，U=U1.证毕.,例3对于例2，系数矩阵,由高斯消去法m21=0，m31=2，m32=-1，故,从而得到求矩阵行列式的计算公式为,5.2.3列主元消去法,由高斯消去法知道,在消元过程中可能有akk(k)0的情况，这时消去法将无法进行；即使主元素akk(k)0但很小时，用其作除数，会导致其它元素数量级的严重增长和舍入误差的扩散，最后也使得计算解不可靠.,先看一个例子(同时参看书上p148-例4).,高斯消去法也称主元素消去法(条件detA0),即当akk(k)=0时，高斯消元法无法进行；或|akk(k)|0(i=1,2,n).由矩阵乘法及ljk=0(当j|c1|0;|bi|ai|+|ci|,ai,ci0,(i=2,3,n-1).|bn|an|0.,我们利用矩阵的直接三角分解法来推导解三对角线方程组(3.12)的计算公式.由系数阵A的特点，可以将A分解为两个三角矩阵的乘积，即A=LU,其中取L下三角阵,取U为单位上三角阵,这样求解方程组Ax=f的方法称为追赶法.设,其中为待定系数，比较(3.13)两边即得,下面我们用归纳法证明,(3.15)对i=1是成立的.现设(3.15)对i-1成立，求证对i亦成立.,由归纳法假设0|i-1|0,使得,考虑n元函数,记S=x|x|=1,xRn,则S是一个有界闭集.由于f(x)为S上的连续函数,所以f(x)于S上达到最大值和最小值,即存在x,xS使得,设xRn且x0,则x/|x|S,从而有,即,可以证明常用范数满足下述关系,以及,注意定理14不能推广到无穷维空间.由定理14可得结论:如果在一种范数意义下向量序列收敛时,则在任何一种范数意义下该向量序列均收敛.,定理15,其中|为向量的任一种范数.,证明显然,而对于Rn上任一种范数|,由定理14,存在常数c1,c20,使,于是又有,5.4.2矩阵范数,下面我们将向量范数概念推广到矩阵上去.可以将Rnn中的矩阵A=(aij)nn当作n2维向量,则由向量的2-范数可以得到Rnn中矩阵的一种范数,称为A的Frobenius范数.|A|F显然满足正定性、齐次性及三角不等式.,下面我们给出矩阵范数的一般定义.,定义5(矩阵的范数)如果矩阵ARnn的某个实值函数N(A)=|A|，满足条件:,(1)|A|0(|A|=0当且仅当A=0)(正定性)；,(2)|cA|=|c|A|,c为实数(齐次性)；,(3)对任意A,B有|A+B|A|+|B|(三角不等式);,(4)对任意A,B有|AB|A|B|(相容性)；,则称N(A)是Rnn上的一个矩阵范数(或模).,上面我们定义的F(A)=|A|F就是Rnn上的一个矩阵范数.,上述条件(即不等式)就称为矩阵范数与向量范数的相容性条件.,由于在大多数与估计有关的问题中,矩阵和向量会同时参与讨论,所以希望引进一种矩阵的范数,它是和向量范数相联系而且和向量范数相容的,即对任何向量xRn及矩阵ARnn都成立,|Ax|A|x|.,为此我们再引进一种矩阵的范数.,定义6(矩阵的算子范数)设xRn,ARnn,给出一种向量范数|x|v(如v=1,2或),相应地定义一个矩阵的非负函数,可验证|A|v满足定义5(见下面定理),所以|A|v是Rnn上矩阵的一个范数,称为A的算子范数.也称从属范数.,定理16设|x|v是Rn上的一个向量范数,则|A|v是Rnn上矩阵的范数,且满足相容条件,证明由(4.6)相容条件(4.7)是显然的.现只验证定义5中条件(4).,由(4.7),有,当x0,有,故,显然这种矩阵的范数|A|v依赖于向量范数|x|v的具体含义.也就是说,当给出一种具体的向量范数|x|v时,相应地就得到了一种矩阵范数|A|v.,定理17设xRn,A=(aij)Rnn,则有常用的算子范数,（称为A的行范数）,（称为A的列范数）,（称为A的2-范数）.,其中max(ATA)表示ATA的最大特征值.由于矩阵2-范数与ATA的特征值有关，所以又称为A的谱范数.证明见书p165.,由定理17看出,计算一个矩阵的|A|,|A|1还是比较容易的,而矩阵的2-范数在|A|2计算上不方便,但是矩阵的2-范数具有许多好的性质,它在理论上是非常有用的.,我们指出,对于复矩阵(即ACnn)定理17中(1),(2)显然也成立,对于(3)应改为,例7设,则,求相应的各种范数.,解,因为,令,得ATA的特征值,故,可见,参见书p166-例7.,补定义设ARnn的特征值为i(i=1,2,n),称,为A的谱半径.,例子设,求A的谱半径.,解,得A的特征值,故得A的谱半径为,定理18(特征值上界)设ARnn,|为任一种算子范数，则(A)|A|(对亦成立).(4.10)反之,对任意实数0,至少存在一种算子范数|A|,使|A|(A)+.(4.11),证明设是A的任一特征值,x为相应的特征向量,则Ax=x,由(4.7)得,注意到|x|0,即得,|A|.,定理后半部分证明可见文献2.,定理19如果ARnn为对称矩阵,则|A|2=(A).,(证明留作练习题或查有关参考书.),定理20如果|B|1,则IB为非奇异矩阵,且,其中|是指矩阵的算子范数.,证明用反证法.若det(I-B)=0,则(I-B)x=0有非零解,即存在x00使Bx0=x0,|Bx0|/|x0|=1,故|B|1,这与假设矛盾.又由(I-B)(I-B)-1=I,有,(I-B)-1=I+B(I-B)-1,从而,类似对加法有,(I+B)(I+B)-1=I,(I+B)-1=I-B(I+B)-1,得,求|A|1,|A|2,|A|,(A).,例子设,解:显然|A|1=4,|A|=4,这里,我们指出,对于实对称矩阵A,有(A)=|A|2.,显然的模最大，得,5.5误差分析,5.5.1矩阵的条件数,由于A(或b)元素是测量得到的,或者是计算的结果,在第一种情况A(或b)常带有某些观测误差,在后一种情况A(或b)又包含有舍入误差.因此我们处理的实际矩阵是A+A(或b+b),下面我们来研究数据A(或b)的微小误差对解的影响.即考虑估计x-y,其中y是(A+A)y=b的解.,考虑线性方程组Ax=b,其中设A为非奇异矩阵,x为方程组的精确解.,首先考察一个例子.,例8设有方程组,记为Ax=b,它的精确解为x=(2,0)T.,现在考虑常数项的微小变化对方程组解的影响,即考察方程组,也可表示为A(x+x)=b+b,其中b=(0,0.0001)T,y=x+x,x为(6.1)的解.显然方程组(5.2)的解为x+x=(1,1)T.,我们看到(5.1)的常数项b的第2个分量只有1/1000的微小变化,方程组的解却变化很大.这样的方程组称为病态方程组.,定义7如果矩阵A或常数项b的微小变化(小扰动),引起方程组Ax=b解的巨大变化,则称此方程组为“病态”方程组,其系数矩阵A称为“病态”矩阵(相对于方程组而言),否则称方程组为“良态”方程组,A称为“良态”矩阵.,应该注意,矩阵的“病态”性质是矩阵本身的特性,下面我们希望找出刻画矩阵“病态”性质的量.设有方程组,Ax=b,(5.3),其中A为非奇异矩阵,x为(5.3)的精确解.以下我们研究方程组的系数矩阵A(或b)的微小误差(小扰动)时对解的影响.,(1)现设A是精确的,x为Ax=b的精确解，当方程组右端有误差b,受扰解为x+x,则,A(x+x)=b+b,Ax=b,x=A-1b,|x|A-1|b|.(5.4),由(5.3)有,|b|A|x|.,于是由(5.4)及(5.5),得,定理21设A是非奇异矩阵,Ax=b0,且,A(x+x)=b+b,则,上式给出了解x的相对误差的上界,常数项b的相对误差在解中放大|A-1|A|倍.,(2)现设b是精确的,x为Ax=b的精确解，当A有微小误差(小扰动)A,受扰解为x+x,则,(A+A)(x+x)=b,(A+A)x=-(A)x.(5.6),如果A不受限制的话,可能A+A奇异,而,(A+A)=A(I+A-1A).,由定理20知,|A-1A|1时,(I+A-1A)-1存在.由(5.6)式,x=-(I+A-1A)-1A-1(A)x.,因此,设|A-1|A|1,即得,定理22设A是非奇异矩阵,Ax=b0,且,(A+A)(x+x)=b.,如果|A-1|A|1时,则(5.3)是“病态”的(即A是“病态”矩阵,或者说A是坏条件的,相对于方程组),当A的条件数相对的小,则(5.3)是“良态”的(或者说A是好条件的).注意,方程组病态性质是方程组本身的特性.A的条件数越大,方程组的病态程度越严重,也就越难用一般的计算方法求得比较准确的解.,例如对前面例8的矩阵作分析,由于条件数cond(A)1很大,可见矩阵A的病态程度十分严重,故由此方程组的解误差非常大.,计算其条件数,通常使用的条件数,有,(1)cond(A)=|A-1|A|;,(2)A的谱条件数;,当A为对称矩阵时,其中1,n为A的绝对值最大和绝对值最小的特征值.,条件数的性质:,(1).对任何非奇异矩阵,都有cond(A)v1.事实上,cond(cA)v=cond(A)v.,(2).设A为非奇异矩阵且c0(常数),则,(3).如果A为正交矩阵,则cond(A)2=1;如果A为非奇异矩阵,R为正交矩阵,则,cond(RA)2=cond(AR)2=cond(A)2.,由性质1知,1cond(A),即cond(A)总是大于等于1的数.条件数反映了方程组的“病态程度”.条件数越小,方程组的状态越好,条件数很大时,称方程组为病态方程组.但多大的条件数才算病态则要视具体问题而定,病态的说法只是相对而言.,条件数的计算是困难的,这首先在于要算A-1,而求A-1比解Ax=b的工作量还大,当A确实病态时,A-1也求不准确；其次要求范数,特别是求|A|2,|A-1|2又十分困难,因此实际工作中一般不先去判断方程组的病态.但是必须明白,在解决实际问题的全过程中,发现结果有问题,同时数学模型中有线性方程组出现,则方程组的病态可能是出问题的环节之一.,病态方程组无论选用什么方法去解,都不能根本解决原始误差的扩大,即使采用全主元消去法也不行.可以试用加大计算机字长,比如用双精度字长计算,或可使问题相对得到解决.如仍不行,则最好考虑修改数学模型,避开病态方程组.,如第3章中提到的拟合问题中出现的正规方程组,就往往呈现病态,此时解决问题的方法之一是避开正规方程组,采用正交多项式拟合的方法,尽管后者比前者在理论上和实际计算中都复杂得多.,例如,n阶希尔伯特矩阵,当n较大时,是有名的病态矩阵.在函数逼近时,有时得到方程组Hnx=b,当n稍大,用任何方法都难以解出理想的结果.考查它的条件数表可见Hn在n稍大时,即呈严重病态.下面我们通过例题来考察Hn的病态情况.,例9已知希尔伯特(Hilbert)矩阵,求条件数cond(H2)2,cond(H2)1和cond(H2).,解因为,下面再分析H3的情况,(2)考虑方程组,(1)计算H3的条件数cond(H3),|H3|=11/6,|H3-1|=408,所以cond(H3)=748.,同样可计算cond(H6)=2.9107,cond(H7)=9.85108.,由此看出当n越大时,Hn矩阵病态越严重.,H3x=(11/6,13/12,47/60)T=b,设H3及b都有微小误差(取3位有效数字)有,简记为(H3+H3)(x+x)=b+b.方程组H3x=b与(5.8)的精确解分别为x=(1,1,1)T,x+x=(1.089512538,0.48967062,1.491002798)T.,于是x=(0.0895,-0.5120,0.4910)T,这就是说H3与b相对误差不超过0.3%,而引起解的相对误差超过50%.,由上面的讨论,要判别一个矩阵是否病态需要计算条件数cond(A)=|A-1|A|,而计算A-1是比较费劲的,那么在实际计算中如何发现病态情况呢?,(1)如果在A的三角约化时(尤其是用主元素消去法解(5.3)时)出现小主元,对大多数矩阵来说,A是病态矩阵,例如用选主元的直接三角分解法解方程组(5.8)(结果舍入为3位浮点数),则有,(2)系数矩阵A的行列式值相对说很小,或系数矩阵某些行近似线性相关,这时A可能病态.,(3)系数矩阵A的元素间数量级相差很大,并且无一定规则,这时A可能病态.,用选主元素的消去法不能解决病态问题,对于病态方程组可采用高精度的算术运算(采用双倍字长进行运算)或者采用预处理方法,即将求解Ax=b转化为一等价方程组,选择非奇异矩阵P,Q使,cond(PAQ)cond(A).,一般选择P,Q为对角阵或者三角矩阵.,当矩阵A的元素大小不均时,对A的行(或列)引进适当的比例因子(使矩阵A的所有行或列按-范数大体上有相同的长度,使A的系数均衡),对A的条件数是有影响的.这种方法不能保证A的条件一定得到改善.,例10设,计算条件数cond(A).,解因为,所以得,现对A的第一行引进比例因子.如用,除第一个方程式,得Ax=b,即,而,于是得,当用列主元消去法解(5.9)时(计算到小数后3位数字),于是得到很坏的结果:x2=1,x1=0.,现用列主元消去法解(6.10),得到,从而得到较好的计算解:x1=1,x2=1.,设x为方程组Ax=b的近似解,于是可计算x的剩余向量r=b-Ax,当r很小时,x是否为Ax=b一个较好的近似解?下面定理给出了解答.,定理23(事后误差估计)设A为非奇异矩阵,x是方程组Ax=b0的