线性方程组的求解

线性方程组的求解第1页，课件共76页，创作于2023年2月§1引言一般的线性方程组本章：即有唯一第2页，课件共76页，创作于2023年2月的求法求法直接求法（解析法）间接方法（数值方法、迭代方法）1.Grammer法则：2.Gauss消元法及其改进3.三角分解法（L-U方法）适用：当A为低阶稠密阵（n<150）.1.一般迭代法2.Jacobi法3.Seidel法4.SOR法适用：当A为高阶稀疏矩阵（n>150)第3页，课件共76页，创作于2023年2月二、预备知识1.几个概念：1°单位上(下)三角阵：主对角元素均为1的上（下）三角阵相关结论：其逆仍为单位上（下）三角阵

其积仍为单位上（下）三角阵2°置换阵（排列阵）：单位阵经过一次两行（两列）的互换形成矩阵为初等置换阵。初等置换阵的乘积即为置换阵初等置换阵一定为对称阵,3°三对角阵:对方阵当称A为三角阵.

4°可约与不可约矩阵：对n阶矩阵A，若存在置换阵P使其中，为r阶子块，为n-r阶子块，则称A为可约矩阵，否则不可约。第4页，课件共76页，创作于2023年2月5°对角占优矩阵：对n阶矩阵A，若或

成立，则称A为严格行或列对角占优矩阵。

若或且上述不等式

至少有一个严格成立，称A为弱对角行（列）占优矩阵。

相关结论：若A为严格行或列对角占优矩阵（或A为不可约行或列弱对角占优矩阵），则A可逆。论证：不妨设A为严格行对角占优矩阵反证：A不可逆，即AX=0有非零解设，对AX=0的第K个方程

第5页，课件共76页，创作于2023年2月6°矩阵的谱半径：

对n阶矩阵A,设为A的全部特征值,称为A的谱半径。

第6页，课件共76页，创作于2023年2月2.再论范数1°2°3°

主要讨论1°

相关结论：中的所有范数均等价,即任意取中的两个范数，则

一定存在两个正数C1以及C2使2°上的范数（乘积用的最多）

定义：中的范数常用第7页，课件共76页，创作于2023年2月例：但在中常用的是所谓的“算子范数”,称

为A的算子范数。易知：称为矩阵范数和向量范数相容。注意：向量范数与矩阵算子范数相容。常用的算子范数为：

行范数：

列范数：

第8页，课件共76页，创作于2023年2月相关结论：1)对任意n阶矩阵A，，（可以是算子范数，也可以是一般范数）

论证：设为A的特征值，为对应的特征向量

两边取算子范数：

第9页，课件共76页，创作于2023年2月若为一般矩阵范数，设

,λ对应的特征向量为ξ，ξ≠0Aξ=λξ,令矩阵,使B≠0,其中它一定存在。2)对任意n阶矩阵A，总存在一个算子范数3)对n阶矩阵A，若A的算子范数‖A‖<1,则E±A可逆,且论证：反证：若E±A不可逆，即（E±A）X=0有非零解ξ≠0即（E±A）ξ=0,即有‖A‖≥1,矛盾。

4）上的范数均等价,即对上的两个范数及,存在常数C1>0,C2>0.使第10页，课件共76页，创作于2023年2月3.向量列及矩阵列的收敛定义1：给定的向量列及向量

若

结论：

论证：提示取范数为定义2：

结论：‖·‖为任意范数

特别当n=1，A=a，

第11页，课件共76页，创作于2023年2月§2Gauss消元法一、基本思想回忆：方程组的初等变换：1)互换两个方程2)某个方程的非零倍数3)某个方程的非零倍数加到另一个方程上

对应线性方程组Ax=b，用矩阵的语言描述：即为增广矩阵（A,b）的初等行变换对列初等变换，除列交换之外，其余两种变换虽可作，但得到的新方程组与原来的方程组不等价

特别注意的是：列交换时须记录！

显然：方程组经过初等变换变成等价的方程组。第12页，课件共76页，创作于2023年2月下面讨论Gauss消元法：对给定的方程组：Ax=b。经过一系列的方程组的初等变换，将Ax=b转化为Ux=g（U为上三角矩阵，称Ux=g为上三角方程组）或Lx=f（L为下三角矩阵，称Lx=f为下三角方程组）然后回代即可求解，此即为Gauss消元法的基本思想。解释“回代”对Ux=g，即

从倒数第一个方程求解，依次求得：类似地Lx=f：

第13页，课件共76页，创作于2023年2月二、计算过程

消元过程一般化为上三角方程组计算过程

回代过程以上三角方程组为例，说明如下：对已知，即为：设doing第14页，课件共76页，创作于2023年2月设：

doing：继续以上过程，注意设得到上三角方程，，可以推出，无需假设然后回代即可求解。第15页，课件共76页，创作于2023年2月注Notes：①称为消元的主元素。如何保证

有两种途径：1）若A可逆，总可以通过行交换使得（k，k）位置上的元素非零2）Th.1:为A的第k阶顺序主子式.

A可逆只能保证第n阶顺序主子式不为0，不能保证其他阶为0.②消元过程的矩阵分析：第16页，课件共76页，创作于2023年2月以此类推：

均为单位下三角矩阵归纳有：(单位下三角矩阵的逆为单位下三角矩阵，乘积也为单位下三角矩阵)由此引进如下更一般的三角分解的概念：定义：对n阶矩阵A，若存在下三角阵L及上三角阵U，使A=L•U，则称L•U为A的三角分解。

当L为单位下三角矩阵时，称L-U为Doolittel分解（D分解）

当U为单位上三角矩阵时，称L-U为Grout分解（G分解）第17页，课件共76页，创作于2023年2月Th.2：对n阶矩阵A，当A的顺序主子式≠0（k=1,2，…，n-1），则A一定

存在D分解且D分解唯一。于是有:Gauss消元顺利进行另，为使主元，只要A可逆，交换AX=b的两两个方程组，可使位置上的元素非零，此时Gauss消元法可描述为，存在置换阵P使P•A=L•U。第18页，课件共76页，创作于2023年2月③Gauss消元法的改进Gauss消元法的缺陷：Gauss消元法的改进：切入点：主元素

方法：选主元+消元+回代第19页，课件共76页，创作于2023年2月具体步骤：第20页，课件共76页，创作于2023年2月（b）列主元素消元法常用

第21页，课件共76页，创作于2023年2月消元④计算量

回代消元：步数12……n-1

除法n-1n-2……1

乘法n(n-1)(n-1)(n-2)……2×1

第22页，课件共76页，创作于2023年2月

§3三角求解求法一、基本思想第23页，课件共76页，创作于2023年2月第24页，课件共76页，创作于2023年2月LU解法计算量与Gauss消元法计算量相同第25页，课件共76页，创作于2023年2月二、特殊的LU解法对已知Ax=b，若A为某些特殊矩阵，LU解法亦有一些特殊变形，讨论以下两种情况。第26页，课件共76页，创作于2023年2月第27页，课件共76页，创作于2023年2月第28页，课件共76页，创作于2023年2月

§4、误差分析第29页，课件共76页，创作于2023年2月一、解对系数的敏感性关于Ax=b的求解，以上所述均为精确解，故误差分析只针对舍入误差，而舍入误差实质上是一种绝对误差。

对Ax=b，可能产生舍入误差的来源是：A或b的微小变化（A及b会产生舍入误差）。问题：A或b的扰动对Ax=b的解是否有影响？第30页，课件共76页，创作于2023年2月例题解得x1=2，x2=0。

扰动b

解得x1=1，x2=1。第31页，课件共76页，创作于2023年2月值得指出：

①A、b的扰动是一种客观存在；②A、b的扰动对Ax=b的解有时会产生破坏作用。

此即为解对系数的敏感性，显然，解对A、b的敏感性越大，对解的破坏性越大。反之亦然。

定义：对给定的线性方程组Ax=b，若A、b的微小扰动对Ax=b的解产生很大的变化，称为病态方程组，A为病态矩阵；反之，称Ax=b为良性方程组，A称为良态矩阵。第32页，课件共76页，创作于2023年2月二、病态性分析给定Ax=b（|A|≠0），则Ax*=b，x*唯一存在。分三种情形：1.b扰动，A不扰动2.A扰动，b不扰动3.A以及b同时扰动第33页，课件共76页，创作于2023年2月1.b扰动，A不扰动b→b+δb，则此时x*→x*+δxA（x*+δx）=b+δb∴

Aδx=δb

∴

δx=A-1δb∴‖δx‖≤‖A-1‖·‖δb‖

又‖b‖=‖Ax*‖≤

‖A‖·‖x*‖1/

‖x*‖≤

‖A‖/‖b‖∴≤

‖A-1‖‖A‖‖δb‖·1/‖b‖易知‖A-1‖‖A‖越大,对解得变化影响越大。反之亦然。第34页，课件共76页，创作于2023年2月2.A扰动，b不扰动A→A+δA，则此时x*→x*+δx

（A+δA）（x*+δx）=b

展开得（A+δA）·δx=-δAx*

又A+δA=A（E+A-1δA）另设||A-1δA||≤||A-1||·||δA||＜1有E+A-1δA可逆，则‖(E+A-1δA)-1‖≤1/(1-‖A-1δA‖)≤1/(1-‖A-1‖‖δA‖)

第35页，课件共76页，创作于2023年2月进一步有：

δx=-(A+δA)-1δA·x*=-(E+A-1δA)-1A-1δAx*

‖δx‖≤‖(E+A-1δA)-1‖‖A-1‖‖δA‖‖x*‖≤1/(1-‖A-1‖‖δA‖)·‖A-1‖‖δA‖‖x*‖∴‖δx‖/‖x*‖≤‖A-1‖‖δA‖/(1-‖A-1‖‖δA‖)

≤

{‖A-1‖‖A‖‖δA‖/‖A‖}/{1-‖A-1‖‖A‖‖δA‖/‖A‖)}易知‖A-1‖‖A‖越大，对系数的扰动起的扩张作用越大，从而解的变化越大。反之亦然。第36页，课件共76页，创作于2023年2月3.A以及b同时扰动b→b+δb，A→A+δA，x*→x*+δx则（A+δA）（x*+δx）=b+δb同2.知得：同样可看出‖A-1‖‖A‖对解的影响一样的。由此引出下面的条件数的概念第37页，课件共76页，创作于2023年2月条件数定义：对任意的n阶可逆矩阵A，称Cond(A)=

‖A-1‖‖A‖

为A的条件数。易知，Cond(A)越大，Ax=b病态程度越大。①常用的条件数Cond(A)∞=

‖A-1‖∞‖A‖∞Cond(A)2=

‖A-1‖2‖A‖2第38页，课件共76页，创作于2023年2月②条件数性质1.Cond(A)≥1Cond(A)=

‖A-1‖‖A‖≥

‖A-1A‖=12.对常数c≠0，Cond(cA)=Cond(A)3.若R为正交矩阵，则Cond(R)2=1

4.对任意的n阶可逆阵A及正交阵R，有Cond(A)2=Cond(RA)2=Cond(AR)25.对任意的n阶矩阵A、B有Cond(AB)≤Cond(A)·Cond(B)第39页，课件共76页，创作于2023年2月三、病态性诊断及处理已知对任意的n阶矩阵A,Cond(A)=

‖A-1‖

‖A‖但实际上的计算很困难例：书中P170，希尔伯特矩阵HnCond(H3)∞=748第40页，课件共76页，创作于2023年2月诊断及处理方法诊断1.出现小主元2.系数阵的行列式很小或某些行近似相关3.系数阵元素的大小分布差异性大处理方法1.高精度算法（双倍字节）2.预处理3.迭代处理第41页，课件共76页，创作于2023年2月四、事后估计设给定Ax=b，并设x^为精确解x*的近似值，则有‖x*-x^‖/‖x*‖≤Cond(A)·‖r‖/‖b‖其中r为剩余向量r=Ax*-Ax^=b-Ax^A(x*-x^)=r∴x*-x^=A-1r

第42页，课件共76页，创作于2023年2月§5.5迭代求解法第43页，课件共76页，创作于2023年2月一、基本思想给定Ax=b，|A|≠0，Ax*=b等价变形x=Bx+f变形的可行性是显然的！a11x1+a12x2+…+a1nxn=b1x1=(1-a11)x1-a12x2-…-a1nxn+b1第44页，课件共76页，创作于2023年2月任取x*的一个近似值X(0)=(x1(0),x2(0),…,xn(0))T视Bx+f为校正器！迭代x(1)=Bx(0)+fx(2)=Bx(1)+f…x(k+1)=Bx(k)+f生成向量列﹛x(k)﹜,若﹛x(k)﹜收敛，则其极限x*。第45页，课件共76页，创作于2023年2月证明：x(k+1)=Bx(k)+f∴limx(k+1)=Blimx(k)+f,k→∞

ξ=Bξ+f

故x=Bx+f的解

又因为解唯一∴ξ=x*

第46页，课件共76页，创作于2023年2月

在﹛x(k)﹜收敛时，取充分大的自然数N，视x(N)≈x*,这种求解方法叫迭代法。命名：称x=Bx+f为迭代方程组

称B为迭代矩阵

称x(k+1)=Bx(k)+f(校正过程)为迭代格式

称﹛x(k)﹜为迭代列第47页，课件共76页，创作于2023年2月NOTES：1.对给定的Ax=b，其等价迭代方程组不唯一!2.不同的初值x(0)或迭代矩阵B不同，得到不同的迭代列(有无穷多个迭代列)，但若收敛，其极限均

为x*。3.如何使得﹛x(k)﹜收敛？第48页，课件共76页，创作于2023年2月

分析：因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（向量），k→∞Bk(x(0)–x*)→0（向量），k→∞Bk→0（矩阵），k→∞

ρ(B)<1(谱半径小于1)

第49页，课件共76页，创作于2023年2月二、几种常用的迭代格式（自始至终）1。Jacobi格式（J-格式,三种形式,分量式,分量浓缩式,矩阵式）①格式设计。已知Ax=b（|A|≠0）

即a11x1+a12x2+…+a1nxn=b1a21x1+a22x2+…+a2nxn=b2…an1x1+an2x2+…+annxn=bn

第50页，课件共76页，创作于2023年2月假设aii≠0,标准等价方程组x1=(-a12x2-a13x3-…-a1nxn+b1)/a11x2=(-a21x1-a23x3-…-a2nxn+b2)/a22

…xn=(-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/ann…0第51页，课件共76页，创作于2023年2月fJ=(b1/a11,b2/a22,…,bn/ann)TJ格式：X(k+1)=BJx(k)+fJ

x1(k+1)=(-a12x2(k)-a13x3(K)-…-a1nxn(k)+b1)/a11x2(k+1)=(-a21x1(k)-a23x3(K)-…a2nxn(k