应用数值分析（三）

线性方程组的数值解法线性方程组Ax=b的一般数值解法：适用于低阶稠密方程组非零元素较多，零元素较少适用于大型稀疏方程组上万阶，零元素很多，非零元素很少评价求解Ax=b数值方法好坏的标准：§5.2线性方程组的性态及条件数

矩阵的基本运算(5)单位矩阵(1)矩阵加法(2)矩阵与标量的乘法(3)矩阵与矩阵乘法(4)转置矩阵(6)非奇异矩阵

称为非奇异矩阵.

如果均为非奇异矩阵，

设如果则称是的逆矩阵，记为且则(7)矩阵的行列式

设则的行列式可按任一行(或列)展开，其中为的代数余子式，即

的余子式.为元素行列式性质矩阵的特征值与谱半径设若存在数（实数或复数）和非零向量，使（1）则称为的特征值，为对应的特征向量，称为矩阵的谱半径.（2）有非零解，

由(1)知可使齐次线性方程组

的全体特征值记为的谱，记作，即故系数行列式，记

为矩阵的特征多项式，方程(3)称为矩阵的特征方程.（3）记行列式展开

因为次代数方程复数域中有个根故故矩阵个特征值是它的特征方程(3)的个根.且（4）记（5）称为的迹.

的特征值和特征向量的其他性质：(1)与有相同的特征值.(2)若非奇异，则的特征值为，特征向量为.(3)相似矩阵有相同的特征多项式.矩阵的特征方程为故特征值为2,2，-7

例

求的特征值及谱半径

的谱半径为解

设(1)对角矩阵

(2)三对角矩阵(3)上三角矩阵(4)海森伯格（Hessenberg）阵(5)对称矩阵

(6)埃尔米特矩阵(7)对称正定矩阵

特殊矩阵(12)设矩阵，若且至少有一个不等式严格成立，则称矩阵为弱对角占优阵，对所有不等式严格成立，则称矩阵

为严格对角占优阵。若(8)正交矩阵(9)酉矩阵(10)初等置换阵单位矩阵交换第行与第行(或交换第列与第列)(11)置换阵

（为交换第行与第行得到的矩阵）；（为交换第列与第列得到的矩阵）；由初等置换阵的乘积得到的矩阵.

定理1

设则下述命题等价：(1)对任何方程组有唯一解.(2)齐次方程组只有唯一解.(4)存在.(5)的秩(3)或的特征值定理2

设为对称矩阵.如果则为对称正定阵.定理3

设为对称正定阵，则(1)

为非奇异矩阵，且亦是对称正定阵.(2)

记为的顺序主子阵，则亦是对称正定矩阵，其中(3)

的特征值(4)

的顺序主子式都大于零，即

定理4(若尔当（Jordan）标准型)设为阶矩阵，则存在一个非奇异矩阵使得其中为若尔当块.(1)

当的若当标准型中所有若尔当块均为一阶时，此标准型变成对角矩阵.

(2)

如果的特征值各不相同，则其若尔当标准型必为对角阵线性方程组的数值解法§5.3Gauss消元法

用消元法解方程组

第2步.解

第1步.得等价的三角形方程组解为上述过程相当于例1

其中表示矩阵的第行.1)消元过程2023/10/27线性方程组的直接解法45Gauss消去法算法消元计算forforfor回代求解for2023/10/27线性方程组的直接解法46计算量/*AmountofComputation*/由于计算机中乘除/*multiplications/divisions*/

运算的时间远远超过加减/*additions/subtractions*/

运算的时间，故估计某种算法的运算量时，往往只估计乘除的次数，而且通常以乘除次数的最高次幂为运算量的数量级。

GaussianElimination:Stepk：设，计算因子且计算共进行n

1步(n

k)次(n

k)2

次(n

k)次(n

k)(n

k+2)次消元乘除次数：1次(n

i+1)次回代乘除次数：2023/10/27线性方程组的直接解法47计算量/*AmountofComputation*/由于计算机中乘除/*multiplications/divisions*/

运算的时间远远超过加减/*additions/subtractions*/

运算的时间，故估计某种算法的运算量时，往往只估计乘除的次数，而且通常以乘除次数的最高次幂为运算量的数量级。

GaussianElimination:Stepk：设，计算因子且计算共进行n

1步(n

k)(n

k+2)次消元乘除次数：回代乘除次数：GaussianElimination的总乘除次数为，运算量为级。n=20时,顺序Gauss消去法需3060次乘除法运算.

定理1

设其中

(1)

如果将约化为等价的三角形方程组则可通过高斯消去法

(2)如果为非奇异矩阵，则可通过高斯消去法(及交换两行的初等变换)将方程组约化为？矩阵在什么条件下才能保证归纳法当时，成立.设时结论成立，即证明:充分性.且有可用高斯消去法将约化到，有即定理2

约化的主元素的充要条件是矩阵的顺序主子式即

由假设推出必要性.

例2

求解方程组用4位浮点数进行计算.精确解舍入到4位有效数字为

解法1

用高斯消去法计算解显然计算解是一个很坏的结果，不能作为方程组的近似解.原因：

在消元计算时用了小主元0.001，使得约化后的方程组元素数量级大大增长，经再舍入使得计算时发生了严重的相消情况，因此经消元后得到的三角形方程组就不准确了.解法2计算解为交换行，避免绝对值小的主元作除数.每一步选取系数矩阵(或消元后的低阶矩阵)中绝对值最大的元素作为主元素，以使高斯消去法具有较好的数值稳定性.线性方程组的数值解法§5.4

基于矩阵三角分解的方法5.4.1

矩阵的三角分解第1步消元其中第2步消元其中

一般第K步消元其中70条件:矩阵A经Gauss消元法后得到的上三角矩阵.

Gauss消元法的矩阵表示上三角阵单位下三角求矩阵的LU分解。取变换矩阵则有再取变换矩阵例解：

Doolittle分解中LU

元素的求解次序

A=LU三角分解的紧凑格式78对r=2,3,…,n,

矩阵三角分解的紧凑格式

例5

用直接三角分解法即Doolittle分解法解解从而求解得求解得直接LU分解元素的存储81解解84练2用直接三角分解法解解85

先求Ly=b得y

再求Ux=y

得x5.4.2

平方根法和改进的平方根法上三角阵单位下三角实对称正定矩阵的平方根法93对称

A为3阶对称正定阵,A=LLT,怎样求L?对称94例对下列矩阵进行Cholesky分解95

A为n阶对称正定阵,A=LLT

L中元素的求解次序

依次求L的第一列,第二列,…,第n列.元素仍然存放在矩阵的相应位置上97

求L的第一列

求L的第二列

对称正定阵A=LLT,求L的计算公式

设已求得L的前k－1列,现求L的第k列(k=3,…,n)101例9求的Cholesky分解.解102解

对称正定阵的LDLT分解中L,D的计算公式1(避免重复计算，便于编程)

对称正定阵的LDLT分解中L,D的计算公式2Step1Step2Step3Stepn储存比消去法节省一半，但比平方根法多用一个单元内存

对称正定阵的LDLT分解紧凑格式及存储

对称正定阵的LDLT分解求方程组（改进平方根法）例10用LDLT解方程组解例10用LDLT解方程组或紧凑格式112

对称正定阵的LDLT分解本质上是对A作Doolittle分解,即LU分解.LDLT分解中的D

即LU分解中的U的对角部分LDLT分解中的L

即LU分解中的L113

对称正定矩阵A的LU分解,计算量可以节省一半

求U的第1行

求L的第1列

对称正定阵的LDLT分解中L,D的计算

先对对称正定阵A作LU分解114

求U的第k行(k=2,3,…,n)

求L的第k列(k=2,3,…,n)

对称正定阵的LDLT分解中L,D的计算节省了计算量

求D5.4.3

三对角线性方程组的追赶法119例四阶三对角矩阵的三角分解（Gauss消去法）120例四阶三对角矩阵的三角分解单位下二对角阵上二对角阵122

三对角矩阵三角分解中LU的求解次序123对k=3,…,n

三对角矩阵三角分解中LU的求解

右端用矩阵乘法展开,比较两边元素得追赶法追的过程:n-1赶的过程:2n-1乘除运算量2n－2

系数矩阵为严格对角占优阵时,追赶法数值稳定

计算过程中撇开系数中大量的零,使计算公式大大简化,减少了计算量.追赶法总的乘除运算量5n－4.线性方程组的直接法小结三角分解的计算过程Step1Step2Step3Step4Step5Step6Step2n-1Step2(n-1)依次交替进行再计算的列先计算的行存储杜利特尔分解/*DoolittleFactorization*/：对方程组求解,只要得到了系数矩阵的三角分解形式,再利用前代算法和回代算法解两个三角方程组即得.线性方程组的数值解法对线性方程组要求寻找更经济、适用的数值解法.其中A为非奇异矩阵,当A为低阶稠密矩阵时,选主元消去法(直接法)是有效的.但对于大型稀疏矩阵方程组(A的阶数n很大

104，但零元素较多),方程组的阶数很高,则运算量将会很大并且大量占用计算机资源.利用迭代法求解是合适的.§5.5

雅可比迭代法和高斯-赛德尔迭代法

将介绍迭代法的一些基本理论及雅可比迭代法，高斯-赛德尔迭代法，超松弛迭代法，而超松弛迭代法应用很广泛。

则(1)式和(2)式同解,称(1)与(2)等价.

例1求解方程组记为Ax=b，其中此方程组的精确解是x*=(3,2,1)T.现将(1.2)改写为迭代法的基本概念

下面举简例，以便了解迭代法的思想.或写为x=B0x+f，其中取x(0)=(0,0,0)T.将这些值代入(1.3)式右边(若(1.3)式为等式即求得方程组的解，但一般不满足)，得到新的值x(1)=(x1(1),x2(1),x3(1))T

=(3.5,3,3)T

，再将x(1)分量代入(1.3)式右边得到x(2)，反复利用这个计算程序，得到一向量序列简写为x(k+1)=B0x(k)+f，其中k表示迭代次数(k=0,1,2,

).

迭代到第10次有一般的计算公式(迭代公式)

从此例看出，由迭代法产生的向量序列x(k)逐步逼近方程组的精确解是x*=(3,2,1)T.即有

对于任何一个方程组x=Bx+f(由Ax=b变形得到的等价方程组)，由迭代法产生的向量序列x(k)是否一定逐步逼近方程组的解x*呢？回答是不一定.请考虑用迭代法解下述方程组但x(k)并不是所有的都收敛到解x*！x*=(-3,-4)T

这种方式就称为迭代法.以上过程称为迭代过程.迭代法产生一个序列则称迭代法收敛,否则称为发散.(3)

如果对于任意其极限存在,即

对线性方程组(2),采用以下迭代步骤:

迭代法的基本思想是通过构造迭代格式产生迭代序列，由迭代序列来逼近原方程组的解，因此，要解决的基本问题有：1.如何构造迭代格式

2.迭代序列是否收敛设线性方程组(1)的一般形式为Jacobi迭代G-S迭代SOR迭代法(4)5.5.1

Jacobi迭代法由(4)建立Jacobi迭代公式为(5)将(5)式转化为矩阵形式矩阵形式为称(5)、(6)式为解线性方程组(1)的Jacobi迭代法(J法)(6)(5)用Jacobi迭代法求解方程组,误差不超过1e-4解:例2依此类推,得方程组满足精度的解为x12迭代次数为12次x4=3.02411.94780.9205d=0.1573x5=3.00031.98401.0010d=0.0914x6=2.99382.00001.0038d=0.0175x7=2.99902.00261.0031d=0.0059x8=3.00022.00060.9998d=0.0040x9=3.00031.99990.9997d=7.3612e-004x10=3.00001.99990.9999d=2.8918e-004x11=3.00002.00001.0000d=1.7669e-004x12=3.00002.00001.0000d=3.0647e-0055.5.2Gauss-Seidel迭代法分析Jacobi迭代法(5)的迭代过程称(7)、(8)式为解线性方程组(1)的Gauss-Seidel迭代法(G-S法)

雅可比迭代法不使用变量的最新信息计算xi(k+1)，而由高斯-塞德尔迭代公式可知，计算x(k+1)的第i个分量xi(k+1)时，利用了已经计算出的最新分量xj(k+1)(j=1,2,

,i-1).

可看作雅可比迭代法的一种改进.高斯—塞德尔迭代公式每迭代一次只需计算一次矩阵与向量的乘法.G-S迭代法有一个明显的优点是：在计算时，只需一组工作单元，以便存放近似解.而J迭代法需要两组工作单元，分别存放

和

.例3解:用Gauss-Seidel迭代法求解方程组,误差不超过1e-4x1=2.50002.09091.2273d=3.4825x2=2.97732.02891.0041d=0.5305x3=3.00981.99680.9959d=0.0465x4=2.99981.99971.0002d=0.0112x5=2.99982.00011.0001d=3.9735e-004x6=3.00002.00001.0000d=1.9555e-004x7=3.00002.00001.0000d=1.1576e-005通过迭代,至第7步得到满足精度的解x7Jacobi迭代法和Gauss-Seidel迭代法均为简单迭代法后面我们会看到，甚至有这样的方程组，Jacobi迭代收敛，而Gauss-Seidel迭代却是发散的.不能说G-S法比J法更好.从例2和例3看出,Gauss-Seidel迭代法的收敛速度比Jacobi迭代法要快（达到同样的精度所需迭代次数较少).但这个结论，在一定条件下才是对的.

练习

用Jacobi和Gauss-Seidel迭代法解方程组

Jacobi迭代

Gauss-Seidel迭代解：

kx1(k)x2(k)x3(k)10.720.830.8420.9711.071.15…………111.0999931.1999931.299991121.0999981.1999981.299997取x(0)=(0,0,0)T计算结果如下：kx1(k)

x2(k)x3(k)10.720.9021.1644…………81.0999981.1999991.3

Jacobi迭代Gauss-Seidel迭代5.5.3迭代法的收敛性设解线性方程组的迭代格式(10)(11)将(10)与(11)相减,得因此迭代法收敛的充要条件可转变为则定理(迭代法收敛的充要条件)

迭代格式收敛的充要条件为(12)

设Ax=b，其中A=D+L+U为非奇异矩阵,且对角矩阵D也奇异，则

(1)解线性方程组的雅可比迭代法收敛的充要条件是

(J)<1，其中J=-D-1(L+U).

(2)解线性方程组的高斯-塞德尔迭代法收敛的充要条件是

(G)<1，其中G=-(D+L)-1U.

(3)雅可比迭代法收敛的充分条件是||J||<1.(4)高斯-塞德尔迭代法收敛的充分条件是||G||<1.

由定理5.5.1和定理5.5.2可知解：方程组的迭代格式为取，问雅可比迭代法是否收敛？若收敛，需要迭代多少次，才能保证各分量的误差绝对值小于10-6？例4：用雅可比方法解下列方程组迭代阵所以要保证各分量的误差绝对值小于10-6，需要迭代14次

求迭代矩阵的谱半径，需要求迭代矩阵的特征值，对高阶矩阵是比较麻烦的。利用定理5.5.1去判别比较方便一些。但在科学及工程计算中，要求解的线性方程组Ax=b,其矩阵A常常具有某些特性.例如，A具有对角占优性质或A是对称正定矩阵等，根据这些矩阵的性质，可以比较方便地判别这些迭代法的收敛性。两种特殊情况1.A为严格(按行)对角占优阵2.A是对称正定矩阵

定义5.5.3(对角占优阵)设A=(aij)n×n.如果A的元素满足称A为严格(按行)对角占优阵.矩阵A的每一行对角元素的绝对值都严格大于非对角元素绝对值之和

定理5.5.3(对角占优定理)

如果A=(aij)n×n为严格对角占优矩阵，则A为非奇异矩阵.

定理5.5.4设方程组Ax=b，如果A为严格对角占优阵，则解Ax=b的Jacobi迭代法,Gauss-Seidel迭代法均收敛.证:(1)对于Jacobi迭代法,其迭代矩阵为根据定理5.5.1Jacobi迭代法收敛(2)对于G—S迭代法,其迭代矩阵为不能使用定理5.5.1,而用定理5.5.2.即从而因此由于可得矛盾由定理5.5.2知G—S迭代法收敛。定理5.5.5设矩阵A对称,且所有对角元素aii>0,则解线性方程组Ax=b的Jacobi迭代法收敛的充要条件是A及2D-A均为正定矩阵,其中D=(a11,…,ann).(2)解线性方程组Ax=b的Gauss-Seidel迭代法收敛的充分条件是A正定.

如果线性方程组系数矩阵A对称正定，则有以下的收敛定理.

例5已知方程组判断雅可比迭代法和高斯—塞德尔法的敛散性？

解

因为系数矩阵A的顺序主子式所以A是正定矩阵.也是正定矩阵，故Jacobi迭代法收敛.

再由A是对称正定矩阵，故

Gauss—Seidel迭代法也收敛.又由于

例6

在线性方程组Ax=b中，证明：当-1/2<a<1时高斯-塞德尔法收敛，而雅可比法只在-1/2<a<1/2时才收敛.即A正定，所以当-1/2<a<1时高斯-赛德尔法收敛.

证明

因为当-1/2<a<1时，得A的顺序主子式

对雅可比法迭代矩阵故只在

(J)=|2a|<1，即|a|<1/2时雅可比法收敛.当a=0.8时,G-S法收敛，

(J)=1.6>1，J法不收敛.Ax=b，例8判别下列方程组用J法和G-S法求解是否收敛因此不能用定理5.5.1，只能用定理5.5.2判断(1)Jacobi法的迭代矩阵解:即Jacobi迭代法收敛(2)Gauss-Seidel法的迭代矩阵同样用定理5.5.2判断即Gauss-Seidel迭代法发散该例说明G-S法发散时而J法却收敛！例9

讨论用Jacobi迭代法和Gauss-Seidel迭代法求解方程组Ax=b时的收敛性，如果收敛，并比较哪种方法收敛较快，其中解:(1)对Jacobi方法，迭代矩阵

(2)对Gaus