第四章线性方程组的求解

第五章解线性方程组的方法,5.1引言5.2高斯消元法5.3三角分解法5.4误差分析5.5迭代求解法,1引言,一般的线性方程组,本章：,即,有唯一,的求法,求法,直接求法,（解析法）,间接方法,（数值方法、迭代方法）,1.Grammer法则：,2.Gauss消元法及其改进,3.三角分解法（L-U方法）适用：当A为低阶稠密阵（n150),二、预备知识,1.几个概念：1单位上(下)三角阵：主对角元素均为1的上（下）三角阵相关结论：其逆仍为单位上（下）三角阵其积仍为单位上（下）三角阵2置换阵（排列阵）：单位阵经过一次两行（两列）的互换形成矩阵为初等置换阵。初等置换阵的乘积即为置换阵初等置换阵一定为对称阵,3三对角阵:对方阵,当,称A为三角阵.,4可约与不可约矩阵：对n阶矩阵A，若存在置换阵P使,其中，,为r阶子块，,为n-r阶子块，则称A为可约矩阵，否则不可约。,5对角占优矩阵：对n阶矩阵A，若或成立，则称A为严格行或列对角占优矩阵。若或且上述不等式至少有一个严格成立，称A为弱对角行（列）占优矩阵。相关结论：若A为严格行或列对角占优矩阵（或A为不可约行或列弱对角占优矩阵），则A可逆。论证：不妨设A为严格行对角占优矩阵反证：A不可逆，即AX=0有非零解设，对AX=0的第K个方程,6矩阵的谱半径：对n阶矩阵A,设为A的全部特征值,称为A的谱半径。,2.再论范数123主要讨论1相关结论：中的所有范数均等价,即任意取中的两个范数，则一定存在两个正数C1以及C2使2上的范数（乘积用的最多）定义：,中的范数常用,例：但在中常用的是所谓的“算子范数”,称为A的算子范数。易知：称为矩阵范数和向量范数相容。注意：向量范数与矩阵算子范数相容。常用的算子范数为：行范数：列范数：,相关结论：1)对任意n阶矩阵A，（可以是算子范数，也可以是一般范数）论证：设为A的特征值，为对应的特征向量两边取算子范数：,若为一般矩阵范数，设,对应的特征向量为，0A=,令矩阵,使B0,其中它一定存在。2)对任意n阶矩阵A，总存在一个算子范数3)对n阶矩阵A，若A的算子范数A0,C20.使,3.向量列及矩阵列的收敛定义1：给定的向量列及向量若结论：论证：提示取范数为定义2：结论：为任意范数特别当n=1，A=a，,2Gauss消元法,一、基本思想回忆：方程组的初等变换：1)互换两个方程2)某个方程的非零倍数3)某个方程的非零倍数加到另一个方程上对应线性方程组Ax=b，用矩阵的语言描述：即为增广矩阵（A,b）的初等行变换对列初等变换，除列交换之外，其余两种变换虽可作，但得到的新方程组与原来的方程组不等价特别注意的是：列交换时须记录！显然：方程组经过初等变换变成等价的方程组。,下面讨论Gauss消元法：对给定的方程组：Ax=b。经过一系列的方程组的初等变换，将Ax=b转化为Ux=g（U为上三角矩阵，称Ux=g为上三角方程组）或Lx=f（L为下三角矩阵，称Lx=f为下三角方程组）然后回代即可求解，此即为Gauss消元法的基本思想。,解释“回代”对Ux=g，即从倒数第一个方程求解，依次求得：类似地Lx=f：,二、计算过程,消元过程一般化为上三角方程组计算过程回代过程以上三角方程组为例，说明如下：对已知，即为：设doing,设：doing：,继续以上过程，注意设得到上三角方程，可以推出，无需假设然后回代即可求解。,注Notes：称为消元的主元素。如何保证有两种途径：1）若A可逆，总可以通过行交换使得（k，k）位置上的元素非零2）Th.1:为A的第k阶顺序主子式.A可逆只能保证第n阶顺序主子式不为0，不能保证其他阶为0.消元过程的矩阵分析：,以此类推：均为单位下三角矩阵归纳有：(单位下三角矩阵的逆为单位下三角矩阵，乘积也为单位下三角矩阵)由此引进如下更一般的三角分解的概念：定义：对n阶矩阵A，若存在下三角阵L及上三角阵U，使A=LU，则称LU为A的三角分解。当L为单位下三角矩阵时，称L-U为Doolittel分解（D分解）当U为单位上三角矩阵时，称L-U为Grout分解（G分解）,Th.2：对n阶矩阵A，当A的顺序主子式0（k=1,2，n-1），则A一定存在D分解且D分解唯一。于是有:Gauss消元顺利进行另，为使主元，只要A可逆，交换AX=b的两两个方程组，可使位置上的元素非零，此时Gauss消元法可描述为，存在置换阵P使PA=LU。,Gauss消元法的改进Gauss消元法的缺陷：Gauss消元法的改进：切入点：主元素方法：选主元+消元+回代,具体步骤：,（b）列主元素消元法常用,消元计算量回代消元：步数12n-1除法n-1n-21乘法n(n-1)(n-1)(n-2)21,3三角求解求法一、基本思想,LU解法计算量与Gauss消元法计算量相同,二、特殊的LU解法对已知Ax=b，若A为某些特殊矩阵，LU解法亦有一些特殊变形，讨论以下两种情况。,4、误差分析,一、解对系数的敏感性,关于Ax=b的求解，以上所述均为精确解，故误差分析只针对舍入误差，而舍入误差实质上是一种绝对误差。对Ax=b，可能产生舍入误差的来源是：A或b的微小变化（A及b会产生舍入误差）。问题：A或b的扰动对Ax=b的解是否有影响？,例题,解得x1=2，x2=0。扰动b解得x1=1，x2=1。,值得指出：A、b的扰动是一种客观存在；A、b的扰动对Ax=b的解有时会产生破坏作用。此即为解对系数的敏感性，显然，解对A、b的敏感性越大，对解的破坏性越大。反之亦然。定义：对给定的线性方程组Ax=b，若A、b的微小扰动对Ax=b的解产生很大的变化，称为病态方程组，A为病态矩阵；反之，称Ax=b为良性方程组，A称为良态矩阵。,二、病态性分析,给定Ax=b（|A|0），则Ax*=b，x*唯一存在。分三种情形：1.b扰动，A不扰动2.A扰动，b不扰动3.A以及b同时扰动,1.b扰动，A不扰动,bb+b，则此时x*x*+xA（x*+x）=b+bAx=bx=A-1bxA-1b又b=Ax*Ax*1/x*A/bA-1Ab1/b易知A-1A越大,对解得变化影响越大。反之亦然。,2.A扰动，b不扰动,AA+A，则此时x*x*+x（A+A）（x*+x）=b展开得（A+A）x=-Ax*又A+A=A（E+A-1A）另设|A-1A|A-1|A|1有E+A-1A可逆，则(E+A-1A)-11/(1-A-1A)1/(1-A-1A),进一步有：x=-(A+A)-1Ax*=-(E+A-1A)-1A-1Ax*x(E+A-1A)-1A-1Ax*1/(1-A-1A)A-1Ax*x/x*A-1A/(1-A-1A)A-1AA/A/1-A-1AA/A)易知A-1A越大，对系数的扰动起的扩张作用越大，从而解的变化越大。反之亦然。,3.A以及b同时扰动,bb+b，AA+A，x*x*+x则（A+A）（x*+x）=b+b同2.知得：同样可看出A-1A对解的影响一样的。由此引出下面的条件数的概念,条件数,定义：对任意的n阶可逆矩阵A，称Cond(A)=A-1A为A的条件数。易知，Cond(A)越大，Ax=b病态程度越大。,常用的条件数Cond(A)=A-1ACond(A)2=A-12A2,条件数性质1.Cond(A)1Cond(A)=A-1AA-1A=12.对常数c0，Cond(cA)=Cond(A)3.若R为正交矩阵，则Cond(R)2=14.对任意的n阶可逆阵A及正交阵R，有Cond(A)2=Cond(RA)2=Cond(AR)25.对任意的n阶矩阵A、B有Cond(AB)Cond(A)Cond(B),三、病态性诊断及处理,已知对任意的n阶矩阵A,Cond(A)=A-1A但实际上的计算很困难例：书中P170，希尔伯特矩阵HnCond(H3)=748,诊断及处理方法,诊断1.出现小主元2.系数阵的行列式很小或某些行近似相关3.系数阵元素的大小分布差异性大处理方法1.高精度算法（双倍字节）2.预处理3.迭代处理,四、事后估计,设给定Ax=b，并设x为精确解x*的近似值，则有x*-x/x*Cond(A)r/b其中r为剩余向量r=Ax*-Ax=b-AxA(x*-x)=rx*-x=A-1r,5.5迭代求解法,一、基本思想,给定Ax=b，|A|0，Ax*=b等价变形x=Bx+f变形的可行性是显然的！a11x1+a12x2+a1nxn=b1x1=(1-a11)x1-a12x2-a1nxn+b1,任取x*的一个近似值X(0)=(x1(0),x2(0),xn(0)T视Bx+f为校正器！迭代x(1)=Bx(0)+fx(2)=Bx(1)+fx(k+1)=Bx(k)+f生成向量列x(k),若x(k)收敛，则其极限x*。,证明：x(k+1)=Bx(k)+flimx(k+1)=Blimx(k)+f,k=B+f故x=Bx+f的解又因为解唯一=x*,在x(k)收敛时，取充分大的自然数N，视x(N)x*,这种求解方法叫迭代法。命名：称x=Bx+f为迭代方程组称B为迭代矩阵称x(k+1)=Bx(k)+f(校正过程)为迭代格式称x(k)为迭代列,NOTES：,1.对给定的Ax=b，其等价迭代方程组不唯一!2.不同的初值x(0)或迭代矩阵B不同，得到不同的迭代列(有无穷多个迭代列)，但若收敛，其极限均为x*。3.如何使得x(k)收敛？,分析：因x(k)x*=Bx(k-1)+f-(Bx*+f)=B(x(k-1)x*)=B2(x(k-2)x*)=Bk(x(0)x*)故有x(k)收敛x(k)x*0（向量），kBk(x(0)x*)0（向量），kBk0（矩阵），k(B)1(谱半径小于1),二、几种常用的迭代格式（自始至终）,1。Jacobi格式（J-格式,三种形式,分量式,分量浓缩式,矩阵式）格式设计。已知Ax=b（|A|0）即a11x1+a12x2+a1nxn=b1a21x1+a22x2+a2nxn=b2an1x1+an2x2+annxn=bn,假设aii0,标准等价方程组x1=(-a12x2-a13x3-a1nxn+b1)/a11x2=(-a21x1-a23x3-a2nxn+b2)/a22xn=(-an1x1-an2x2-an,n-1xn-1+bn)/ann矩阵形式：X=BJX+fJ其中BJ=0-a12/a11-a13/a11-a1n/a11-a21/a220-a23/a22-a2n/a22-a31/a33-a32/a330-a3n/a33-an1/ann-an2/ann-an3/ann0,fJ=(b1/a11,b2/a22,bn/ann)TJ格式：X(k+1)=BJx(k)+fJx1(k+1)=(-a12x2(k)-a13x3(K)-a1nxn(k)+b1)/a11x2(k+1)=(-a21x1(k)-a23x3(K)-a2nxn(k)+b2)/a22xn(k+1)=(-an1x1(k)-an2x2(K)-an,n-1xn-1(k)+bn)/ann类似可写浓缩式：,J格式的矩阵分析：对Ax=b,2。Seidel格式（S-格式）定位：是对J格式的改进（保留Jacobi假设）目的：节约空间方法：在迭代计算中尽可能地利用最新的数据,具体格式设计：,S-格式的矩阵分析：起点：此为计算模式下的S-