Transcript 第五章
第5章
解线性方程组的直接方法
1
5.1
引言与预备知识
5.2
高斯消去法
5.3
矩阵三角分解法
5.4
向量和矩阵的范数
5.5
误差分析
5.1
5.1.1
引言与预备知识
引言
线性方程组的数值解法一般有两类：
1.
直接法
经过有限步算术运算，可求得方程组精确解的方法(若
计算过程中没有舍入误差).
但实际计算中由于舍入误差的存在和影响，这种方法
也只能求得线性方程组的近似解.
2
2.
迭代法
是用某种极限过程去逐步逼近线性方程组精确解的方法.
3
5.1.2
向量和矩阵
用 R mn 表示全部 m n 实矩阵的向量空间， C mn 表
示全部 m n 复矩阵的向量空间.
A R mn
a11 a12 a1n
a a a
2n
A (aij ) 21 22
a
a
a
mn
m1 m 2
这种实数排成的矩形表，称为 m 行 n 列矩阵.
x1
x
xRn x 2
xn
称为 n 维列向量.
4
写成列向量的形式
A a1 a2 an ,
其中 ai 为 A的第 i 列.
也可写成行向量的形式
b1T
T
b2
A ,
bT
m
其中 biT 为 A的第 i 行.
5
矩阵的基本运算：
(1) 矩阵加法 C A B,
cij aij bij ( A, B,CR mn )
(2) 矩阵与标量的乘法 C A,
cij aij .
(3) 矩阵与矩阵乘法 C AB,
n
cij aik bkj ( AR mn , BR n p ,C R m p ).
k 1
(4) 转置矩阵 AR mn , C AT , cij aij .
(5) 单位矩阵 I e1 e2 en R nn ,
6
其中
ek 0,,0,1,0,,0 , k 1,2,,n.
T
(6) 非奇异矩阵
设 A R nn , B R nn .如果 AB BA I , 则称 B 是
1 T
T 1
1 且
的逆矩阵，记为
(
A
)
(
A
) .
A
,
A
如果 A1 存在， 则称 A为非奇异矩阵.
如果 A, B R nn 均为非奇异矩阵， 则 ( AB) 1 B 1 A1.
(7) 矩阵的行列式
设 A R nn , 则 A的行列式可按任一行(或列)展开，
7
即
n
det ( A) aij Aij
(i 1,2,, n),
j 1
其中 Aij为 aij 的代数余子式，Aij (1)i j M ij , M 为元素
ij
aij的余子式.
行列式性质：
(a) det ( AB) det ( A)det ( B), A, BR nn .
(b) det ( AT ) det ( A), AR nn .
(c) det (cA) c n det ( A), cR, AR nn .
(d) det ( A) 0 A是非奇异矩阵.
8
5.1.3
定义1
矩阵的特征值与谱半径
设 A aij R nn，若存在数 （实数或复数）
和非零向量 x ( x1 , x2 , , xn )T R n，使
Ax x,
（1.1）
则称 为 A的特征值，x为 A对应 的特征向量， A 的全
体特征值记为 A的谱，记作 ( A) ，即 ( A) {1 , 2 ,, n }.
记
( A) max i .
1i n
称为矩阵 A的谱半径.
（1.2）
9
由(1.1)可知 可使齐次线性方程组
(I A) x 0
有非零解，故系数行列式 det( I A) 0 ，记
a11
a21
p ( ) det( I A)
an1
a12
a22
an 2
a1n
a2 n
（1.3）
ann
n c1n 1 cn 1 cn 0.
p ( )称为矩阵 A的特征多项式，方程(1.3)称为矩阵 A 的特
征方程.
10
因为 n 次代数方程 p ( ) 在复数域中有 n个根
1 , 2 ,, n ,
故
p( ) ( 1 )( 2 ) ( n ).
由(1.3)中的行列式展开可得
n
c1 1 2 n aii ,
i 1
cn ( 1) n 12 n ( 1) n det A.
故矩阵 A a R nn 的 n 个特征值 , , , , 是它
ij
1
2
n
的特征方程(1.3)的 n个根.
11
并且
det( A) 12 n ,
（1.4）
及
n
n
tr A aii i .
i 1
（1.5）
i 1
称 trA为 A的迹.
A 的特征值 和特征向量 x 还有一下性质：
(1) AT 与 A有相同的特征值 及特征向量 .
(2)若 A 非奇异，则 A1 的特征值为 1，特征向量为 x .
(3)相似矩阵 B S 1 AS 有相同的特征多项式.
12
例1
求 A 的特征值及谱半径
1
A 2
2
解
2
2
4 .
2
2
4
矩阵 A的特征方程为
2
2
1
det( I A) 2
2 4
2
4 2
3 32 24 28
( 2) 2 ( 7) 0,
故 A特征值为 1 2 2,
3 7, A 的谱半径为 ( A) 13 7.
5.1.4
特殊矩阵
设 A (aij ) R nn .
(1) 对角矩阵 如果当i j
时，aij 0.
(2) 三对角矩阵 如果当i j 1，aij 0.
(3) 上三角矩阵 如果当i j
时，aij 0.
(4) 上海森伯格（Hessenberg）阵
如果当i j 1 时，aij 0.
(5) 对称矩阵 如果AT A.
14
(6) 埃尔米特矩阵 设A Cnn , 如果AH A.
(7) 对称正定矩阵
如果 (a) AT A,
(b) 对任意非零向量 x R n , ( Ax, x) xT Ax 0.
(8) 正交矩阵 如果A1 AT .
(9) 酉矩阵 设A Cnn ,
如果A1 AH .
(10) 初等置换阵
由单位矩阵 I 交换第 i 行与第 j 行(或交换第 i 列与
第 j 列)，得到的矩阵记为 I ij ,且
15
~（为交换 第 行与第 j 行得到的矩阵）；
A i
I ij A A
AI ij B（为交换 A第 i 列与第 j列得到的矩阵）；
(11) 置换阵 由初等置换阵的乘积得到的矩阵.
定理1
设 A R nn 则下述命题等价：
(1) 对任何 b R n , 方程组 Ax b有唯一解.
(2) 齐次方程组 Ax 0只有唯一解 x 0.
(3) det( A) 0.
(4) A1存在.
(5) A的秩 rank ( A) n.
16
定理2
设 A R nn 为对称正定阵，则
(1) A为非奇异矩阵，且 A1亦是对称正定阵.
(2) 记 Ak为 A的顺序主子阵，则
Ak (k 1,2,, n) 亦是对称正定矩阵，其中
a11
Ak
a
k1
a1k
akk
(k 1,2,, n).
(3) A 的特征值 i 0(i 1,2,, n).
(4) A 的顺序主子式都大于零，即
det( Ak ) 0 (k 1,2,,n).
17
定理3
设 A R nn 为对称矩阵.如果 det ( Ak ) 0
(k 1,2,,n), 或 A 的特征值 i 0(i 1,2,, n), 则 A为
对称正定阵.
定理4(若尔当（Jordan）标准型) 设 A为 n 阶矩阵，
则存在一个非奇异矩阵 P 使得
J1 (1 )
J 2 (2 )
1
P AP
,
J r (r )
18
其中
i
J i (i )
1
i
i
1
i
ni ni
r
ni 1(i 1,2,, r ), 且 ni n.
i 1
为若尔当块.
(1) 当 A的若当标准型中所有若尔当块 J i 均为一阶时，
此标准型变成对角矩阵.
19
(2) 如果 A 的特征值各不相同，则其若尔当标准型必为
对角阵 diag (1 , 2 ,, n ).
20
5.2
高斯消去法
21
5.2.1
高斯消去法
设有线性方程组
a11 x1 a12 x2 a1n xn b1 ,
a x a x a x b ,
21 1
22 2
2n n
2
am1 x1 am 2 x2 amn xn bm .
（2.1）
或写为矩阵形式
a11 a12 a1n x1 b1
a21 a22 a2 n x2 b2
,
a
m1 am 2 amn xn bm
简记为 Ax b.
22
例2
用消去法解方程组
x1 x2 x3 6,
4 x2 x3 5,
2 x 2 x x 1.
2
3
1
解
( 2.2)
(2.3)
( 2.4)
第1步.将方程(2.2)乘上 2加到方程(2.4)上去，
消去(2.4)中的未知数 x1 , 得到
4 x2 x3 11.
(2.5)
第2步. 将方程(2.3)加到方程(2.5)上去，消去方程
(2.5)中的未知数 x2 ,
23
得到与原方程组等价的三角形方程组
x1 x2 x3 6,
4 x2 x3 5,
2 x3 6.
(2.6)
显然，方程组(2.6)是容易求解的，解为
x (1,2,3)T .
上述过程相当于
1
A b 0
2
1
1
4
1
2
1
1
6
5 0
0
1
1
1
4
1
4
1
6
5
11
24
1
0
0
1
1
4
1
0
2
1
5
6
(2) r1 r3 r3
r2 r3 r3
其中用 ri 表示矩阵的第 i行.
由此看出，用消去法解方程组的基本思想是用逐次消
去未知数的方法把原方程组 Ax b 化为与其等价的三角
形方程组，而求解三角形方程组可用回代的方法.
上述过程就是用行的初等变换将原方程组系数矩阵化
为简单形式(上三角矩阵)，从而将求解原方程组(2.1)的
问题转化为求解简单方程组的问题.
25
或者说，对系数矩阵 A施行一些左变换(用一些简单矩
阵)将其约化为上三角矩阵.
下面讨论求解一般线性方程组的高斯消去法.
将(2.1)记为 A(1) x b (1) , 其中
A(1) (aij(1) ) (aij ),
b(1) b.
(1)x b ,
x1 a(12
).设
a
a11
k x2 1
(1)
第1步
a111
0,首先计算乘数
n n
1
a x a x a x b ,
21 1
22 (1)2
2n n
2
(1)
m
a
/
a
(
i
2
,
3
,
, m) .
i1
i1
11
（2.1）
m
i
用
乘(2.1)的第一个方程，加到第
个方程上，
aim1 1 x1 am 2 x2 amn xn bm .
26
(i 2,3, , m)消去(2.1)的从第二个方程到第 m个方程中
的未知数 xi , 得到与(2.1)等价的方程组
(1)
a11
0
0
(1)
a12
( 2)
a22
( 2)
am
2
a1(1n) x1 b1(1)
( 2)
( 2)
a2 n x2 b2
.
( 2 )
( 2)
amn
xn bm
（2.7）
简记为 A( 2) x b ( 2) ,
其中 A( 2) , b ( 2)的元素计算公式为
( 2)
(1)
(1)
a
a
m
a
ij
ij
i1 1 j
( 2)
(1)
(1)
bi bi mi1b1
i, j 2,3, , n,
i 2,3, , n.
27
(2) 第 k 次消元
(k 1,2, , n).
设上述第1步，…，第 k 1步消元过程计算已经完成，
即已计算好与(2.1)等价的方程组
(1)
a11
(1)
a12
( 2)
a22
简记为 A( k ) x b ( k ) .
a1(1k)
a2( 2k)
(k )
akk
(k )
amk
a1(1n) x1 b1(1)
( 2)
( 2)
a2 n x2 b2
,
(k )
(k )
akn xk bk
( k )
b(n)
x
amn
n
n
（2.8）
28
(k )
设 akk
0,计算乘数
(k )
(k )
mik aik
/ akk
(i k 1, , n) .
用 mik 乘(2.8)的第 k个方程，加到第 i个方程(i k 1, , n),
消去从第 k 1个方程到第 n个方程中的未知数 xk , 得到与
(2.1)等价的方程组 A( k 1) x b ( k 1) .
A( k 1) ,b ( k 1) 元素的计算公式为
( k 1)
(k )
(k )
a
a
m
a
ij
ij
ik kj
( k 1)
(k )
(k )
b
b
m
b
i
ik k
i
i, j k 1, , n,
i k 1, , n,
显然 A( k 1) 中从第1行到第 k 行与 A(k )相同.
（2.9）
29
(k )
(3) 继续上述过程，且设
akk
0(k 1,2, , n 1),
直到完成第 n 1步消元计算.
最后得到与原方程组等价的简单方程组 A( n ) x b ( n ) ,
(1)
a11
(1)
a12
( 2)
a22
a1(1n) x1 b1(1)
( 2)
( 2)
a2 n x2 b2
.
( n )
(n)
ann
xn bn
（2.10）
由(2.1)约化为(2.10)的过程称为消元过程.
30
如果 AR nn 是非奇异矩阵，且 akk( k ) 0(k 1,2,,n 1),
求解三角形方程组(2.10)，得到求解公式
(n)
xn b ( n ) / ann
,
(k )
x
(
b
k
k
n
(k )
(k )
a
x
)
/
a
kj j kk (k n 1,n 2,,1).
j k 1
（2.11）
(2.10)的求解过程(2.11)称为回代.
如果 a11 0, 由于 A 为非奇异矩阵，所以 A 的第一列一
定有元素不等于零.
31
例如 ai 1 0, 于是交换两行元素(即 r1 ri ），将 ai 1
1
1
1
调到(1,1)位置，然后进行消元计算，这时 A( 2) 右下角矩阵
为 n 1 阶非奇异矩阵.
继续这过程，高斯消去法照样可进行计算.
32
定理5
设 Ax b, 其中 AR nn .
(1) 如果 akk( k ) 0(k 1,2,,n), 则可通过高斯消去法
将 Ax b约化为等价的三角形方程组(2.10).
且计算公式为：
(a) 消元计算 (k 1,2,,n 1)
(1)
a11
(1)
(1)
b1(1)
a
a
(k )
( k1)n x1
12
mik aik / akk
(i k 1
, n),
,
( 2)
( 2)
( 2)
a22
a2 n x2 b2
（2.10）
.
a ( k 1)
( k
)
ikakj
1,,n),
aij( k ) m
ij
(i, j k
( n )
b(n)
x
ann
n
n
b ( k 1) b ( k ) m b ( k ) (i k 1,, n).
i
ik k
i
33
(b) 回代计算
(n)
xn b ( n ) / ann
,
n
(i )
(i )
(i )
x
(
b
a
x
)
/
a
i
i
ij
j
ii
j i1
(i n 1, n 2,,1).
(2) 如果 A 为非奇异矩阵，则可通过高斯消去法(及交
换两行的初等变换)将方程组 Ax b约化为(2.10).
3
3
n
n
n
以上消元和回代过程的总乘除法次数为 n 2 ,
3
3
3
3
2
3
(1)
(
1
)
(
1
)
(
1
)
n an x b
a11 na12 n
加减法次数为
1n . 1 1
2)
3
a ( 2 ) 2 3 a (3
x b ( 2 )
22
2
.
( n )
(n)
ann
xn bn
2n
2
（2.10）
34
高斯消去法对于某些简单的矩阵可能会失败，例如
0
A
1
1
.
0
由此，需要对前述算法进行修改，首先需要研究原来的矩
(k )
0(k 1,2,).
阵 A在什么条件下才能保证 akk
35
定理6
约化的主元素 aii(i ) 0(i 1,2,,k ) 的充要条件
是矩阵 A的顺序主子式 Di 0(i 1,2,,k ).即
D1 a11 0,
a11
Di
ai1
证明
a1i
0
(i 1,2,, k ).
（2.12）
aii
首先利用归纳法证明定理6的充分性.
显然，当 k 1时，定理6成立.
现设定理6充分性对 k 1是成立的，求证定理6充分性
对 k 亦成立.
36
设 Di 0(i 1,2,,k ), 于是由归纳法假设有 aii( i ) 0
(i 1,2,, k 1), 可用高斯消去法将 A(1)约化到 A(k )，即
A(1) A( k )
(1)
a11
(1)
a12
( 2)
a22
a1(1k)
a2( 2k)
(k )
akk
(k )
ank
a1(1n)
( 2)
a2 n
,
(k )
akn
(k )
ann
且有
(1)
a11
D2
0
(1)
a12
(1) ( 2 )
a
11 a22 ,
( 2)
a22
37
(1)
a11
Dk
a1(1k)
(1) ( 2 )
(k )
a11
a22 akk
.
（2.13）
(k )
akk
(k )
由设 Di 0(i 1,2,,k ), 利用(2.13)式， 则有 akk
0,
定理6充分性对 k 亦成立.
显然，由假设 aii(i ) 0(i 1,2,,k ), 利用(2.13)式亦可
推出 Di 0(i 1,2,,k ).
38
推论
如果 A 的顺序主子式 Dk 0(k 1,2,,n 1), 则
(1)
a11
D1 ,
(k )
akk
Dk / Dk 1
(k 2,3,,n).
39
5.2.2
矩阵的三角分解
下面借助矩阵理论进一步对消去法作些分析，从而建
立高斯消去法与矩阵因式分解的关系.
设(2.1)的系数矩阵 AR nn 的各顺序主子式均不为零.
由于对 A施行行的初等变换相当于用初等矩阵左乘 A，于
是对(2.1)施行第一次消元后化为(2.7)，这时 A(1) 化为化
为 A( 2)，
, b (1)化为 b ( 2 ) ,即
L1 A(1) A( 2) ,
L1b (1) b ( 2) ,
40
其中
1
m21
L1 m31
m
n1
1
1
.
1
)
一般第 k 步消元， A(k )化为 A( k 1), b ( k )化为 b ( k 1，
相当于
Lk A( k ) A( k 1) ,
Lk b ( k ) b ( k 1) ,
其中
41
1
Lk
1
mk 1,k
mn ,k
1
.
1
重复这过程，最后得到
Ln1L2 L1 A(1) A( n ) ;
Ln1L2 L1b (1) b ( n ) .
（2.14）
记上三角矩阵 A(n )为 U，由(2.14)得到
1
1
1
2
1
n1
A L L L U LU ,
42
其中
L L11 L21Ln11
1
m21
m31
m
n1
1
m32
1
mn 2
mn 3
1
为单位下三角矩阵.
这就是说，高斯消去法实质上产生了一个将 A分解为
两个三角形矩阵相乘的因式分解，于是得到如下重要定理，
它在解方程组的直接法中起着重要作用.
43
定理7（矩阵的LU分解）
设 A为 n阶矩阵，如果 A
的顺序主子式 Di 0(i 1,2,,n 1),则 A可分解为一个单位
下三角矩阵 L和一个上三角矩阵 U的乘积，且这种分解是
唯一的.
证明
现在在
根据以上高斯消去法的矩阵分析，存在性已得证，
为非奇异矩阵的假定下证明唯一性，
A
设
A LU L1U1 ,
其中 L,L1 为单位下三角矩阵， U ,U1为上三角矩阵.
44
由于 U 11存在，故
L1 L1 UU11.
上式右边为上三角矩阵，左边为单位下三角矩阵，从而上
式两边都必须等于单位矩阵，故 L L1 ,U U1 ,
唯一性得证.
45
例3
对于例2，系数矩阵
1
A 0
2
1
1,
1
1
4
2
由高斯消去法，m21 0,m31 2,m32 1, 故
1
A 0
2
0
1
1
0 1
0 0
0
1
1
4
0
1
1 LU .
2
46
5.2.3
高斯列主元素消去法
(k )
由高斯消去法知道，在消元过程中可能出现 akk
0
这时消去法将无法进行；
(k )
即使主元素 akk
0, 但很小时，用其作除数，会导致
其他元素数量级的严重增长和舍入误差的扩散，最后也使
得计算解不可靠.
例3
求解方程组
0.001
1.000
2.000
2.000
3.712
1.072
3.000 x1 1.000
4.623 x2 2.000 .
5.643
x3 3.000
47
用4位浮点数进行计算. 精确解舍入到4位有效数字为
x* (0.4904, 0.05104, 0.3675)T ,
解法1
用高斯消去法
0.001
A b 1.000
2.000
2.000
3.000
3.712
4.623
1.072
5.643
0.001
0
0
2.000
3.000
2004
3005
4001
6006
1.000 m21 1.000 / 0.001
1000
2.000
m 2.000 / 0.001
3.000
31
2000
1.000
m 4001/ 2004
1002 23
1.997
2003
48
0.001
0
0
2.000
3.000
2004
3005
0
5.000
1.000
1002 ,
2.000
计算解为
x (0.400, 0.09980, 0.4000)T .
显然计算解是一个很坏的结果，不能作为方程组的近似解.
其原因是我们在消元计算时用了小主元 0.001，使得
约化后的方程组元素数量级大大增长，经再舍入使得在计
算(3,3)元素时发生了严重的相消情况，因此经消元后得
到的三角形方程组就不准确了.
49
解法2
交换行，避免绝对值小的主元作除数.
2.000
A b 1.000
r1 r3
0.001
1.072
5.643
3.712
4.623
2.000
3.000
2.000
0
0
1.072
5.643
3.176
1.801
2.001
3.003
2.000
0
0
1.072
5.643
3.176
1.801
0
1.868
3.000
m21 0.5000
2.000
m31 0.0005
1.000
3.000
m 0.6300
0.5000 32
1.002
3.000
0.5000 ,
0.6870
50
得计算解为
x (0.4900, 0.05113, 0.3678)T x*.
这个例子告诉我们，在采用高斯消去法解方程组时，
(k )
小主元可能产生麻烦，故应避免采用绝对值小的主元素 akk
.
对一般矩阵来说，最好每一步选取系数矩阵(或消元后
的低阶矩阵)中绝对值最大的元素作为主元素，以使高斯消
去法具有较好的数值稳定性. 这就是全主元素消去法.
在选主元时要花费较多机器时间，目前主要使用的是
列主元消去法.
51
本节主要介绍列主元消去法，并假定(2.1)的 AR nn
为非奇异的.
设方程组(2.1)的增广矩阵为
a11
a21
B
a
n1
a12
a22
an 2
a1n
a1n
ann
b1
b2
.
bn
首先在 A的第一列中选取绝对值最大的元素作为主元素，
例如
ai1 ,1 max ai ,1 0,
1in
52
然后交换 B的第1行与第 i1 行，经第1次消元计算得
( A b) ( A( 2 ) b ( 2 ) ).
重复上述过程，设已完成第 k 1步的选主元素，交换两行
及消元计算， ( A b) 约化为
( A( k )
a11
b(k ) )
a12
a1k
a1n
a21
a2 k
akk
ank
a2 n
akn
ann
b1
b2
,
bk
bn
53
其中 A(k )的元素仍记为 aij , b ( k ) 的元素仍记为 bi .
第 k 步选主元素(在 A(k )右下角方阵的第1列内选)，
即确定 ik ，使
aik ,k max aik 0.
k in
交换 ( A( k ) b ( k ) ) 第 k 行与 ik 行的元素，再进行消元计算，
最后将原方程组化为
a11
a12
a22
a1n x1 b1
a2 n x2 b2
x b
ann
n n
.
54
回代求解
xn b / ann ,
n
xi (bi aij x j ) / aii
(i n 1, n 2,,1).
j i1
算法1(列主元素消去法)
设 Ax b . 本算法用具有行交换的列主元素消去法，
消元结果冲掉 A，乘数 mij冲掉 aij ，计算解 x冲掉常数项
b, 行列式存放在 det中.
1. det 1
55
2. 对于 k 1,2,,n 1
(1) 按列选主元
aik ,k max aik
k in
(2) 如果 ai
k
,k
0 ，则计算停止 (det A 0)
(3) 如果 ik k 则转(4)
交换行： akj ai , j ( j k , k 1,, n)
k
bk bik
det det
(4) 消元计算
, ,n
对于 i k 1
56
(a) aik mik aik / akk
(b) 对于 j k 1,,n
aij aij mik *akj .
(c)
(5)
bi bi mik *bk
det akk *det
3. 如果 ann 0，则计算停止 (det A 0)
4. 回代求解
(1) bn bn / ann
(2) 对于 i n 1,,2,1
bi (bi
n
a
j i 1
ij
* b j ) / aii
57
5. det ann *det
例3的解法2用的就是列主元素消去法.
列主元素消去法也可用矩阵运算描述：
L1 I1,i1 A(1) A( 2 ) ,
Lk I k ,ik A( k ) A( k 1) ,
L1 I1,i1 b (1) b ( 2 ) ,
Lk I k ,ik b( k ) b( k 1) .
（2.15）
其中 Lk 的元素满足 mik 1 (k 1,2,,n 1),
I k ,ik 是初等置换阵.
58
利用(2.15)得到
Ln1I n1,in1 L2 I 2,i2 L1I1,i1 A A( n ) U .
若记
~
P Ln1 I n1,in1 L2 I 2,i2 L1 I1,i1 .
(1)
( 2)
L
I
b
b
,
1 1,i1
1
~
~
(n )
P A U ,
Pb b ,
Lk I k ,ik A( k ) A( k ~1) , Lk I k ,ik b( k ) b( k 1) .
考虑 n 4 时的 P .
则有 L1 I1,i A(1) A( 2 ) ,
（2.15）
U A( 4) L3 I 3,i3 L2 I 2,i2 L1 I1,i1 A
L3 ( I 3,i3 L2 I 3,i3 )( I 3,i3 I 2,i2 L1 I 2,i2 I 3,i3 )
~ ~ ~
( I 3,i3 I 2,i2 I1,i1 ) A L3 L2 L1 PA,
（2.16）
59
其中
~
L1 I 3,i3 I 2,i2 L1 I 2,i2 I 3,i3 ,
~
L2 I 3,i3 L2 I 3,i3 ,
~
L3 L3 ,
P I 3,i3 I 2,i2 I1,i1 .
~
Lk (k 1,2,3) 为单位下三角阵，其元素的绝对值不超过1.
记
~ ~ ~
L L3 L2 L1 ,
1
由(2.16)得到
PA LU .
其中 P为排列矩阵， L 为单位下三角阵， U为上三角阵.
60
这说明对(2.1)应用列主元素消去法相当于对 ( A b)先
进行一系列行交换后对 PAx Pb再应用高斯消去法.
而在实际计算中只能在计算过程中做行的交换.
定理8 (列主元素的三角分解定理)如果 A为非奇异矩阵,
则存在排列矩阵 P使
PA LU ,
其中 L 为单位下三角阵，U 为上三角阵.
编程时， L 元素存放在数组 A的下三角部分， U元素
存放在 A上三角部分，由记录主行的整型数组 Ip(n) 可知
P 的情况.
61
5.3
矩阵三角分解法
5.3.1
直接三角分解法
将高斯消去法改写为紧凑形式，可以直接从矩阵 A的
元素得到计算 L, U元素的递推公式，而不需任何中间步骤，
这就是直接三角分解法.
一旦实现了矩阵 A的 LU分解，那么求解 Ax b 的问
①
Ly b,
求 y;
②
Ux y, 求 x.
62
1. 不选主元的三角分解法
设 A为非奇异矩阵，且有分解式
A LU ,
其中 L为单位下三角阵， U为上三角阵，即
1
l21
A
l
n1
1
ln 2
u11
1
L, U的元素可以由
u12
u 22
u1n
u2 n
.
u nn
（3.1）
n步直接计算定出，其中第 r 步
定出 U 的第 r 行和 L 的第 r 列元素.
63
由(3.1)有：
a1i u1i
得U的第1行元素；
(i 1,2,, n),
u1n
ai1 l1
li1 ai1 / u11 (iu112,u,12
n),
得
L的第1列元素
.
i1u11 ,
1
u 22 u 2 n
l21
A
.
（3.1）
设已经定出
的第1行到第
行元素与 L的第1列
r
1
U
到第 r 1ln列元素.
ln 2 1
u nn
1
由(3.1)，利用等式两边元素比较及当 r k 时，lrk 0,
有
n
故
r 1
a
lrk u ki ulrk u ki u u ri .
1 ri
k1 11
12
k 1
1
u 22
l21
A
l
n1 ln 2 1
u1n
u2 n
.
u nn
64
（3.1）
r 1
u ri ari lrk u ki
(i r , r 1, , n),
k 1
又由(3.1)有
n
r 1
k 1
k 1
air lik u kr lik u kr lir u rr .
用直接三角分解法解
u(要求
u1n
Ax
1
u
A 的所有顺序主子式
11 b
12
1
u 22 u 2 n
l21
都不为零)的计算公式.
A
.
),
① ln1u1i lna21i (i
1,2,
,
,
n
1 ,n), li1 ai1 / u11 (i 2u
nn
（3.1）
计算 U的第 r 行和 L 的第 r 列元素 (r 2,3, , n).
r 1
②
u ri ari lrk u ki
k 1
(i r , r 1, , n);
（3.2）
65
③
r 1
lir (air lik ukr ) / urr
(i r 1,, n;且r n);
k 1
（3.3）
求解 Ly b, Ux y 的计算公式；
④
y1 b1 ;
i1
yi bi lik yk
（3.4）
(i 2,3,, n);
k 1
⑤
xn y n / u nn ;
n
xi yi uik xk / uii
k i1
(i n 1, n 2,,1).
（3.5）
66
例5
用直接三角分解法解
1
2
3
解
2
5
1
3 x1 14
2 x2 18 .
5
x3 20
用分解公式(3.2)—(3.3)计算
u11 a11 1,
u12 a12 2,
l21 a21 / u11 2/11,
r 1
u13 a13 3,
l31 a31 / u11 3/1 3
u ri ari lrk u ki (i r , r 1, , n); （3.2）
u22
a22 l21uk12
5 22 1,
1
r 1
ukr)2
/ urr23
(i
r4
,1,, n;且r n);（3.3）
ulir23(a
air23
l21lu
ik13
k 1
l32 (a32 l31u12 ) / u22 (132) /1 5,
67
u33 a33 l31u13 l32u23 5 33 (5)(4) 24.
从而
1
A 2
3
0
1
5
0 1
0 0
0
1
2
1
0
3
4 LU .
24
按(3.4)求解 L y (14,18,20)T , 得
求解
y (14,10,72)T ,
x (1,2,3)T .
U x (14,10,72)T , 得
y1 b1 ;
（3.4）
i1
yi bi lik yk (i 2,3,, n);
k 1
68
再考虑存储问题
当 u ri计算好后 ari就不用了，故可将 u ri仍存放在 ari
的相应位置. 例如
a11
a21
A
a
31
a
41
a12
a22
a13
a23
a32
a33
a42
a43
a14
u11
a24
l21
a34
l
31
l
a44
41
u12
u 22
u13
u23
l32
u33
l42
l43
u14
u 24
u34
u 44
最后在存放 A的数组中得到 L, U 的元素.
直接分解法大约需要 n 3 / 3次乘除法，和高斯消去法
计算量基本相同.
69
如果已经有了 A的 LU分解，则求解具有相同系数的
方程组 Ax (b1 b2 bm ) 是相当方便的，每解一个方程
组 Ax b j 仅需要增加 n 2次乘除法运算.
矩阵 A的分解公式(3.2)，(3.3)又称为杜利特尔分解.
2. 选主元的三角分解法
r 1
从直接三角分解公式可看出当
urr 0时计算将中断，
u ri ari lrk u ki (i r , r 1, , n); （3.2）
k 1
或者当 urr绝对值很小时，按分解公式计算可能引起舍入误
r 1
lir (air lik ukr ) / urr (i r 1,, n;且r n); （3.3）
k 1
差的累积.
如果 A 非奇异，可通过交换 A的行实现矩阵 PA的 LU
分解.
70
采用与列主元消去法类似的方法，将直接三角分解法
修改为(部分)选主元的三角分解法.
设第 r 1步分解已完成，这时有
u11
l21
A lr 1,1
l
r1
l
n1
u12
u 22
lr 1, 2
lr 2
ln 2
u1, r 1
u 2 , r 1
u1r
u2 r
u r 1, r 1
lr , r 1
ln , r 1
u r 1, r
arr
anr
u1n
u2 n
u r 1, n .
arn
ann
71
第 r 步分解需用到(3.2)及(3.3)式，为了避免用小的数 urr作
除数，引进量
r 1
u
于是有ri
si rair
lik u kr (i r 1,, n).
1
ari lrkku1ki (i r , r 1, , n);
r 1
（3.2）
k 1
lir (air
r (1
urrlikuskrr ,) /lu
si (/is
i ,
r ,n1;,且
r,n).n); （3.3）
ir rr
r
k 1
取 max si si , 交换 A的 r 行与 ir行元素，将 si调到 (r , r )
r i n
r
r
位置(将 (i, j )位置的新元素仍记为 lij及 aij )， 于是有
lir 1
(i r 1,,n),
由此再进行第 r 步分解计算.
72
算法2
如果 A为非奇异矩阵，设 Ax b，其中 A为
非奇异矩阵.
1. 对于 r 1,2,, n
(1) 计算 si
r 1
air si air lik ukr
(i r , r 1, , n),
k 1
(2) 选主元 si max si ,
r
r i n
Ip (r ) ir
(3) 交换 A的 r 行与 ir 行元素
ari air ,i
(i 1,2, , n)
(4) 计算 U的第 r行元素， L 的第 r 列元素
73
arr urr sr
air lir si / urr air / arr
r 1
ari uri ari lrk uki
(i r 1,, n,
(i r 1,, n,
且r n)
且r n)
k 1
求解 Ly b, Ux y.
2. 对于 i 1,2,, n 1
(1) t Ip(i )
(2) 如果 i t 则转(3)
bi bt
74
(3) (继续循环)
i 1
3. bi bi lik bk
(i 2,3, , n)
k 1
4.
n
bn bn / u nn , bi bi uik bk / uii
k i 1
(i n 1, n 2, ,1)
利用算法2的结果 PA LU , 可以计算 A 的逆矩阵
A1 U 1 L1 P.
步骤：
(1) 计算上三角矩阵的逆阵 U 1 ;
(2) 计算 U 1 L1 ;
75
(3)
U 1 L1列(利用 Ip(n) 最后记录).
3
A1大约需要 n 次乘法运算.
76
5.3.2
平方根法
平方根法是利用对称正定矩阵的三角分解而得到的
求解对称正定方程组的一种有效方法.
设 A为对称矩阵，且 A的所有顺序主子式均不为零，
定理7知， A 可惟一分解为
1
l21
A
l
n1
1
ln 2
u11
1
u12
u 22
为了利用 A 的对称性，将 U再分解为
u1n
u2 n
.
u nn
77
u11
U
u 22
1
u nn
u12
u11
1
u 23
u 22
u1n
u11
u2 n
u 22
1
DU 0 ,
其中 D为对角阵，U 0为单位上三角阵. 于是
A LU LDU 0 .
（3.6）
又
T
T
A AT U 0 ( DL ),
由分解的惟一性, 得
U 0T L.
代入(4.6)得到对称矩阵 A 的分解式 A LDLT .
78
定理9 (对称阵的三角分解定理)
设 A 为 n阶对称
阵，且 A的所有顺序主子式均不为零，则 A可唯一分解
为
A LDLT ,
其中 D为对角阵， L 为单位下三角阵.
设 A为对称正定矩阵，则在分解式 A LDLT中，D
的对角元素 d i 均为正数.
事实上，由 A的对称正定性，有
d1 D1 0, di Di / Di1 0
(i 2,3,,n).
于是
79
d1
D
1
2
dn
dn
d1
d1
dn
1
2
D D ,
由定理6得到
1
2
1
2
1
2
1
2
A LDL LD D L ( LD )( LD )T L1 LT1 ,
T
T
1
2
L1 LD 为下三角矩阵.
80
定理10 (对称正定矩阵的三角分解或Cholesky分解)
如果 A为 n阶对称正定矩阵，则存在一个实的非奇异下三
角阵 L , 使 A LLT , 当限定 L 的对角元素为正时，这种
分解是唯一的.
可以用直接分解方法来确定计算 L元素的递推公式.
因为
l11
l21
A
l
n1
l22
ln 2
l11
lnn
l21
l22
ln1
ln 2
,
lnn
81
其中 lii 0(i 1,2,, n). 由矩阵乘法及 l jk 0(当j k时),
按等式两边对应元素相等，得
n
j 1
k 1
k 1
aij lik l jk lik l jk l jj lij ,
于是得到解对称正定方程组 Ax b的平方根法计算公式：
对于
l.
2.
j 1,2, , n
l jj
1
2
2
a jj l jk
,
k 1
j 1
j 1
lij aij lik l jk
/ l jj
k 1
（3.7）
(i j 1, , n);
82
求解 Ax b, 即求解两个三角形方程组：
(1)
(2)
Ly b,
求y;
求x.
LT x y,
3.
i 1
yi bi lik yk / lii
k 1
4.
n
xi bi lki xk / lii
k i 1
(i 1,2, , n).
（3.8）
(i n, n 1, ,1).
由计算公式(3.7)知
j
a jj l 2jk
( j 1,2, , n),
k 1
所以
83
l 2jk a jj max {a jj },
1 j n
于是
max {l 2jk } max {a jj }.
j ,k
1 j n
这个结果说明，分解过程中元素 l jk
l jj 恒为正数.
于是不选主元素的平方根法是一个数值稳定的方法.
当求出 L 的第 j 列元素时， LT
j 行元素亦算出.
所以平方根法约需 n 3 / 6 次乘除法，大约为一般直接分解
法计算量的一半.
84
由于 A为对称阵，因此在计算机实现时只需存储 A的
下三角部分.
下三角部分共需存储 n(n 1) / 2个元素，可按行主序
用一维数组存放，即
A(n (n 1) / 2) {a11, a21, a22 ,, an1 , an 2 ,, ann }.
矩阵元素 aij
A(i (i 1) / 2 j ), L 的
元素存放在 A的相应位置.
用平方根法解对称正定方程组时，需要用到开方运算.
为了避免开方，用定理9
85
A LDLT ,
即
1
l21
A
l
n1
1
ln 2
d1
1
d2
1
d n
l21
1
ln1
ln 2
.
1
由矩阵乘法，并注意 l jj 1, l jk 0 ( j k ), 得
n
n
aij ( LD) ik ( L ) kj lik d k l jk
T
k 1
k 1
j 1
lik d k l jk lij d j l jj .
k 1
86
于是得到计算 L 的元素及 D的对角元素公式：
对于 j 1,2,, n
l.
2.
j 1
lij
aij lik d k l jk
/dj
k 1
( j 1,2, , i 1);
i 1
d i aii lik2 d k .
（3.9）
k 1
为避免重复计算，引进
tij lij d j ,
由(3.9)得到按行计算 L, T 元素的公式：
d1 a11
87
对于
1.
i 2,3, , n.
j 1
tij aij tik l jk
( j 1,2, , i 1);
k 1
2.
3.
lij tij / d j
( j 1,2, , i 1);
（3.10）
i 1
d i aii tik lik
k 1
计算出 T LD的第 i 行元素 tij ( j 1,2, , i 1) 后，
存放在 A的第 i 行相应位置，然后再计算 L的第 i 行元素，
存放在 A的第 i 行.
D的对角元素存放在 A的相应位置.
88
例如
a11
a21
A
a
31
a
41
d1
l21
l
31
l
41
a22
对称
a32
a33
a42
a43
d2
l32
d3
l42
l43
d1
l21
l
31
t
a44
41
d2
l32
d3
t 42
t 43
a44
.
d 4
对称正定矩阵 A按 LDLT分解和按 LLT分解计算量差
LDLT分解不需要开方计算.
89
求解 Ly b, DLT x y 计算公式
4.
y1 b1 ;
i 1
y b
lik yk
i
i
k 1
5.
xn y n / d n ;
n
x y / d
lki xk
i
i
i
k i 1
(i 2, , n).
（3.11）
(i n 1, ,2,1).
计算公式(3.10)，(3.11)称为改进的平方根法.
j 1
tij aij tik l jk
( j 1,2, , i 1);
k 1
lij tij / d j
( j 1,2, , i 1);
90
i 1
d i aii tik lik
k 1
（3.10）
91
5.3.3
追赶法
实际问题中，通常要求解系数矩阵为对角占优的三对
b1
a2
c1
b2
c2
an 1
bn 1
an
x1
x2
cn 1 xn 1
bn xn
f1
f2
,
f n 1
f n
简记为 Ax f . 其中，当i j 1时，aij 0,
(a )
b1 c1 0;
（3.12）
且：
92
(b)
bi ai ci , ai ,ci 0
(c)
bn an 0.
(i 2,3,,n 1);
利用矩阵的直接三角分解法推导解三对角线方程组
（3.12)的计算公式.
由系数阵 A的特点，可以将 A分解为两个三角阵的乘积，
即
A LU ,
其中 L 为下三角矩阵， U为单位上三角矩阵.
93
设
b1
a2
A
c1
b2
1
2
c2
an 1
bn 1
an
2
n
其中 i , i , i为待定系数.
cn 1
bn
1
n
（3.13）
1
1
,
n 1
1
94
比较(3.13)两边得
b1 1 , c1 11 ,
ai i , bi i i1 i , (i 2,,n),
ci i i (i 2,, n 1)
（3.14）
由 1 b1 0， b1 c1 0，1 c1 / b1 , 得 0 1 1.
1
1 1
用归纳法可以证明
1
2 2
A
,
i ci 0 (i 1,2,, n 1), n1
n n
1
即 0 i 1, 从而由(4.14)可求出 i
（3.15）
95
i 1时由对角占优的条件，(3.15)是成立的.
现设(3.15)对 i 1成立，求证对 i 亦成立.
由归纳法假设 0 i 1 1, 又由(3.15)及 A 的对角占
优条件，有
（3.15）
i ci 0 (i 1,2,, n 1),
i bi ai i1 bi ai i1 bi ai ci 0,
也就是
0 i 1.
由(3.14)得到
i bi ai i1
(i 2,,n);
96
i ci /(bi ai i1 )
(i 2,3,,n 1).
确定了 { i },{ i },{ i } ，就实现了 A 的 LU分解.
求解 Ax f , 等价于求解两个三角形方程组：
(1) Ly f ,
求y;
(2) Ux y,
求x.
解三对角线方程组的追赶法公式：
1. 计算 { i } 的递推公式
1 c1 /b1
i ci /(bi ai i1 )
(i 2,3,,n 1);
97
2. 解 Ly f
y1 f1 / b1 ,
yi ( f i ai yi1 ) /(bi ai i1 )
(i 2,3,,n);
3. 解 Ux y
xn yn ,
xi yi i xi1
(i n 1,n 2,,2,1).
称计算系数 1 2 n 1及 y1 y2 yn
的过程为追的过程.
称计算方程组的解 xn xn1 x1 的过程为赶的
过程.
98
定理11 设有三对角线方程组 Ax f , 其中 A满足条件
(a)，(b)，(c)，则 A为非奇异矩阵且追赶法计算公式中的
{ i },{ i } 满足：
1° 0 i 1 （i 1,2,, n 1);
2° 0 ci bi ai i bi ai
(i 2,3,,n 1);
0 bn an n bn an .
追赶法公式实际上就是把高斯消去法用到求解三对角
线方程组上去的结果.
99
由于 A 特别简单，因此求解的计算公式也非常简单，
计算量也很小.
计算机实现时只需用三个一维数组分别存储 A的三条线
元素 {ai },{bi },{ci }，此外再用两组工作单元保存 { i }, { yi }
或 {xi }.
100
5.4
向量和矩阵的范数
101
5.4.1
向量范数
向量范数概念是三维欧氏空间中向量长度概念的推广，
在数值分析中起着重要作用.
定义2 设
n
x ( x1 , x2 , , xn )T , y ( y1 , y2 , , yn )T R （或
C n ）.
将实数
n
( x, y ) y x xi yi
T
i 1
n
（或复数 ( x, y ) y x xi yi ）
H
i 1
称为向量 x, y的数量积.
102
将非负实数
x
n 2
( x, x) xi
i 1
1
2
2
1
2
或
x
2
( x, x )
1
2
n
xi
i 1
2
1
2
称为向量 x 的欧氏范数 .
关于范数，成立如下定理.
定理12
n
n
x
,
y
R
（
或
C
), 则
设
1. ( x, x) 0,
当且仅且 x 0 时成立；
103
2. (x, y) ( x, y), 为实数,
(或 (x, y) ( x, y), 为复数）；
3. ( x, y ) ( y, x) (或 ( x, y ) ( y, x));
4. ( x1 x2 , y) ( x1 , y) ( x2 , y);
5. (Cauchy-Schwarz不等式)
( x, y) x 2 y 2 ,
等号当且仅当 x与 y 线性相关时成立；
6. 三角不等式
x y
2
x 2 y 2.
104
向量的欧式范数可以看成是对 R n中向量“大小”的一
种度量.
也可以用其他办法来度量向量的“大小”.
例如，对于 x ( x1 , x2 )T R 2 , 可以用一个 x 的函数
N ( x) max xi
i 1, 2
来度量 x 的“大小”，而且这种度量“大
小”的方法计算起来比欧氏范数方便.
一般要求度量向量“大小”的函数
N ( x)
满足正定性、
齐次性和三角不等式.
105
定义3 (向量的范数) 如果向量 x R n(或 C n )的某
个实值函数 N ( x) x ，满足条件：
1. x 0
( x 0 当且仅且 x 0 )
(正定条件)，
2. x x , R n (或 Cn ),
3.
x y x y
（4.1）
(三角不等式)，
则称 N (x)是 R n (或 C n )上的一个向量范数(或模).
由条件3
x x y y x y y .
106
y y x x y x x .
从而有
4.
x y x y .
（4.2）
几种常用的向量范数.
1. 向量的 -范数(最大范数)：
x
max xi .
1i n
2. 向量的1-范数：
n
x
1
xi .
i 1
107
3. 向量的2-范数：
x
1
2
2
n
1
2
( x, x) ( xi2 ) .
i 1
也称为向量 x的欧氏范数.
4. 向量的 p -范数：
n
x
p
( xi
p
)1/ p ,
i 1
其中 p [1, ) .
可以证明向量函数
N ( x) x
n 上向量的范数,
是
R
p
且容易说明上述三种范数是 p -范数的特殊情况.
108
例6 计算向量 x (1,2,3)T 的各种范数.
解
x 1 6,
定义4
x
3,
x 2 14.
设 {x (k ) } 为 R n中一向量序列，x*R n , 记
x ( k ) ( x1( k ) , x2( k ) ,, xn( k ) )T , x* ( x1* , x2* ,, xn* )T .
如果 lim xi( k ) xi* (i 1,2,, n), 则称 x ( k ) 收敛于向量 x *，
k
记为
lim x ( k ) x* .
k
109
定理13 （ N ( x)的连续性）设非负函数
N ( x) x
为 R n上任一向量范数，则 N (x) 是 x 的分量 x1 , x2 ,, xn
的连续函数.
证明
n
n
i 1
i 1
设 x xi ei , y yi ei ,
其中 ei (0,,1,0,,0),
只须证明当 x y时 N ( x) N ( y ). 事实上
N ( x) N ( y ) x y
n
( x
i 1
i
x y
yi )ei
110
n
xi yi ei
i 1
x y
n
e
i 1
i
,
即
N ( x) N ( y ) c x y
其中
c
0 (当x y时）
,
n
i 1
ei .
111
定理14 (向量范数的等价性) 设
x s, x
n上向量
为
R
t
的任意两种范数，则存在常数 c1 , c2 0,使得对一切 x R n
有
c1 x
s
x
t
c2 x s ,
证明 只要就
x
s
x
证明上式成立即可，即证明
存在常数 c1 , c2 0, 使
c1
x
x
t
c2 , 对一切xR n 且x 0.
考虑泛函
f ( x) x t 0,
xR n .
112
记 S {x x
1, x R n }, 则 S 是一个有界闭集.
由于 f (x) 为 S 上的连续函数，所以 f (x) 于 S 上达到
最大最小值，即存在 x, x S 使得
f ( x) min f ( x) c1 ,
xS
设 x R 且 x 0, 则
n
x
c1 f
x
f ( x) max f ( x) c2 .
xS
x
S , 从而有
x
c2 ,
（4.3）
显然 c1 , c2 0, 上式为
113
c1
x
x
c2 ,
t
即
c1 x
n
x
c
x
,
对一切
x
R
.
2
t
定理14不能推广到无穷维空间. 由定理14可得到结论：
如果在一种范数意义下向量序列收敛时，则在任何一种范
数意义下该向量序列均收敛.
114
定理15
lim x ( k ) x* lim x ( k ) x * 0,
k
k
其中‖·‖为向量的任一种范数.
证明
显然， lim x ( k ) x* lim x ( k ) x * 0,
k
k
而对于 R n 上任一种范数‖·‖，由定理14，存在常数 c1 , c2 0,
使
c1 x ( k ) x*
x ( k ) x* c2 x ( k ) x* ,
于是又有
lim x ( k ) x *
k
0 lim x ( k ) x * 0.
k
115
矩阵范数
5.4.2
向量范数概念可以推广到矩阵.
R
R 中的向量，
n2
nn
R 上的2范数
n2
R nn 中矩阵的一种范数
F ( A) A
2
a
i
,j
i
,
j
1
n
F
1
2
,
称为 A 的弗罗贝尼乌斯范数.
A
F
显然满足正定性、齐次性及三角不等式.
116
定义5 (矩阵的范数) 如果矩阵 A R nn
实值函数
1.
N ( A) A
,满足条件
A 0 ( A 0 Ax 0 ) (正定条件),
（4.4）
n
2. cA c A , c R (齐次条件);
3.
4.
A B A B
(三角不等式);
AB A B .
则称 N ( A)
R nn 上的一个矩阵范数(或模).
117
上面定义的 F ( A) A F
R nn上的一个矩阵范数.
由于在大多数与估计有关的问题中，矩阵和向量会同
时参与讨论，所以希望引进一种矩阵的范数，它和向量范
数相联系.
如要求对任何向量 x R n及 A R nn 都成立
Ax A x .
（4.5）
这时称矩阵范数和向量范数相容.
118
定义6 (矩阵的算子范数)
给出一种向量范数
x
v
(如
设 x R n, A R nn，
v 1,2 或∞)，相应地定义一
个矩阵的非负函数
A
可以验证
v
A
max
x0
Ax
x
v
（4.6）
.
v
满足定义4，所以
v
A
nn 上矩阵的一
是
R
v
个范数，称为 A 的算子范数，也称从属范数.
119
定理16
设
x
v
是 R n上一个向量范数，则
A
nn
是
R
v
上矩阵的范数，且满足相容条件
Ax
v
Av x
v
（4.7）
.
证明 由(4.6)相容性条件(4.7)是显然的.
现只验证定义5中条件4.
由(4.7)，有
Ax
v
ABx
A
.
A v vmax v Bx
.v A v B v x v （4.6）
x0
x v
当 x 0 时，有
ABx
x
v
v
A
v
B
v
.
120
故
AB
max
v
ABx
x0
显然这种矩阵范数
x
A
v
v
A v B
A
v
.
v
依赖于具体的向量范数
也就是说，给出一种具体的向量范数
到一种矩阵范数
v
x
v
x
v
.
，相应地就可得
.
121
定理17
1.
2.
3.
A
设 x R n , A R nn ，则
max
1i n
A 1 max
1 j n
A
2
n
a
ij
j 1
n
a
i 1
ij
max ( AT A)
(称为A的行范数 ),
(称为A的列范数),
(称为A的2 范数）.
其中 max ( AT A) 表示 AT A 的最大特征值.
122
证明 只就1，3给出证明，2同理.
1. 设 x ( x1 , , xn )T 0, 不妨设 A 0 . 记
t x
n
max aij ,
max xi ,
1in
1in
则
Ax
max
1in
n
a
j 1
ij
j 1
n
xj
max aij x j
i
j 1
n
t max aij .
i
j 1
这说明对任何非零 x R n ，有
Ax
x
.
（4.8）
123
接下来说明有一向量 x0 0 , 使
Ax0
x0
.
n
T
设 ai j ，取向量 x0 ( x1 , , xn ) , 其中
j 1
0
x j sgn( ai0 j ) ( j 1,2,, n).
显然
x0
1 ,且 Ax 的第 i 个分量为
0
0
n
n
a
j 1
i0 j
x j ai0 j
j 1
这说明
Ax0
max
1in
n
a
j 1
ij
n
x j ai0 j .
j 1
3. 由于对一切 xR n , Ax 2 ( Ax, Ax) ( AT Ax, x) 0,
2
从而 AT A 的特征值为非负实数，
1 2 n 0.
（4.9）
124
AT A为对称矩阵，设 u1 , u2 ,, un 为 A 的相应于(4.9)
的特征向量且 (ui , u j ) ij ，又设 x R n为任一非零向量，
于是有
n
x ci ui ,
i 1
其中 ci 为组合系数，则
n
Ax
x
2
2
2
2
T
( A Ax, x)
( x, x )
2
c
i 1
i 1
n
2
c
i
1.
i 1
另一方面，取 x u1 ，则上式等号成立，故
A
max
2
x0
Ax
x
2
2
1 max ( AT A) .
125
2
，计算 A的各种范数.
4
例7
1
设 A
3
解
A 1 max{ 1 3 , 2 4 } 6,
A
A
F
max{ 1 2 , 3 4 } 7,
12 (2) 2 (3) 2 42 5.477,
A 2 max ( AT A) 15 221 5.46.
126
对于复矩阵(即 A C nn )定理17中的第1,2项显然也成
立，3应改为
1/ 2
x H A H Ax
A 2 max
H
x 0
x x
max ( AH A) .
127
定理18
对任何 A R nn , 为任一种算子范数,则
( A) A
（对 A F也成立） （4.10）
反之，对任意实数 0,至少存在一种算子范数 ，使
A ( A) .
证明
设 是 A的任一特征值，
（4.11）
x为相应的特征向量，
则 Ax x ，由相容性条件 (4.7) 得
x x Ax A x ,
Ax v A v x
注意到 x 0 , 即得 A
v
.
（4.7）
128
( A).
定理19
如果 A R nn 为对称矩阵，则 A
定理20
如果 B 1，则 I B为非奇异矩阵，且
( I B)
1
1
,
1 B
2
（4.12）
其中‖·‖是指矩阵的算子范数.
证明
用反证法. 若 det( I B) 0，则 ( I B) x 0
有非零解，即存在 x0 0 ,使 Bx0 x0 , Bx 0 1 ,
x0
故 B 1，与假设矛盾.
又由 ( I B)( I B) 1 I，有
129
( I B) 1 I B( I B) 1 ,
从而
( I B ) 1 I B ( I B ) 1 ,
( I B)
1
1
.
1 B
130
5.5
误差分析
5.5.1
矩阵的条件数
考虑线性方程组
Ax b,
其中设 A为非奇异矩阵， x 为方程组的精确解.
由于 A(或 b )元素是测量得到的，或者是计算的结果，
在第一种情况 A (或 b )常带有某些观测误差.
在后一种情况 A (或 b )又包含有舍入误差.
因此我们处理的实际矩阵是 A A (或 b b ).
131
下面研究数据 A(或 b )的微小误差对解的影响.
即考虑估计 x y ，其中 y 是 ( A A) y b 的解.
例8
设有方程组
1
1
x1 2
1 1.0001
2
.
x
2
（5.1）
记为 Ax b, 它的精确解为 x (2,0)T .
现在考虑常数项的微小变化对方程组解的影响，即考
1 x1
2
1
1 1.0001 x 2.0001
,
2
（5.2）
132
也可表示为 A( x x) b b , 其中 b (0,0.0001)T ,
y x x, x 为(5.1)的解.
方程组(5.2)的解为 x x (1,1)T .
可以看到(5.1)的常数项 b 的第2个分量只有
1
的
10000
微小变化，方程组的解却变化很大.
1
1
x1 2 这样的方程组称为病
（5.1）
态方程组.
1 1.0001
2
.
x
1
1 x21 2
1 1.0001 x 2.0001
,
2
（5.2）
以下给出病态性的定义.
133
定义7
如果矩阵 A 或常数项 b 的微小变化，引起方
程组 Ax b解的巨大变化，则称此方程组为“病态”方
程组，矩阵 A称为“病态”矩阵(相对于方程组而言)，否
则称方程组为“良态”方程组， A 称为“良态”矩阵.
矩阵的“病态”性质是矩阵本身的特性，下面希望找出
刻画矩阵“病态”性质的量.
设有方程组
Ax b
（5.3）
其中 A为非奇异阵， x 为(5.3)的准确解.
134
设 A是精确的，b有误差 b，解为 x x , 则利用
A( x x) b b,
Ax b
有
从而
x A1b,
x A1 b .
（5.4）
由 Ax b
b A x,
1
x
A
b
于是由(5.4)及(5.5)，得
(设b 0).
（5.5）
135
定理21
设 A 是非奇异阵，Ax b 0，且
A( x x) b b,
则
x
x
A
1
A
b
b
.
上式给出了解的相对误差的上界，常数项 b 的相对误
差在解中可能放大 A1 A 倍.
现设 b是精确的，A 有微小误差(扰动) A , 解为 x x ,
则
( A A)( x x) b,
( A A)x (A) x. （5.6）
136
如果不对 A 加以限制， A A可能奇异，而
A A A( I A1A),
由定理20知，当 A1A 1时， ( I A1A) 1 存在.
由(5.6)式可推出
x ( I A1A) 1 A1 (A) x,
因此
x
A1 A x
1 A (A)
1
.
137
设 A1 A 1
x
x
A
1
A
A
1 A
1
A
A
A
.
（5.7）
A
138
定理22 设 A为非奇异矩阵，Ax b 0，且
( A A)( x x) b,
如果 A1 A 1 ，则(6.7)式成立.
如果 A充分小，且在条件 A1 A 1下，(5.7)式
说明矩阵 A 的相对误差 A / A 在解中可能放大 A1 A 倍.
量 A1 A 愈小，由 A (或 b )的相对误差引起的解的
相对误差就愈小；
量 A1 A 愈大，解的相对误差就可能愈大.
所以量 A1 A 实际上刻画了解对原始数据变化的
灵敏程度，即刻画了方程组的“病态”程度，
139
定义8
设 A为非奇异阵，称数 cond ( A)v A1
v
A
v
(v 1,2或)为矩阵 A的条件数.
矩阵的条件数与范数有关.
由上面讨论知，当 A 的条件数相对的大，即 cond ( A) 1
时， Ax b是“病态”的(即 A 是“病态”矩阵，或者说 A
是坏条件的，相对于解方程组).
当 A 的条件数相对的小，则 Ax b是“良态”的(或者说
A 是好条件的).
方程组病态性质是方程组本身的特性. A 的条件数愈大，
方程组的病态程度愈严重，也就愈难用一般的计算方法求得
140
比较准确的解.
常用的条件数有
(1) cond ( A) A1
A ;
(2) A 的谱条件数
cond ( A) 2 A
1
2
A
2
max ( AT A)
.
T
min ( A A)
当 A为对称矩阵时
1
cond ( A) 2
,
n
其中 1 , n为 A 的绝对值最大和绝对值最小的特征值.
141
条件数的性质：
1. 对任何非奇异矩阵 A，都有 cond ( A) v 1.
事实上，
cond ( A) v A1
v
A
v
A1 A 1;
v
2. 设 A为非奇异阵且 c 0(常数)，则
cond (cA) v cond ( A) v ;
3.
如果 A 为正交矩阵，则 cond ( A) 2 1 ;
如果 A 为非奇异矩阵， R 为正交矩阵，则
cond ( RA) 2 cond ( AR) 2 cond ( A) 2 .
142
例9
已知希尔伯特（Hilbert）矩阵
1
1
Hn 2
1
n
1
2
1
3
1
n 1
1
n
1
n 1 ,
1
2n 1
计算 H 3的条件数.
解
1
1
H3
2
1
3
1
2
1
3
1
4
1
3
1
,.
4
1
5
9
1
H 3 36
30
36
192
180
30
180 .
180
143
(1) 计算 H 3 条件数 cond ( H 3 )
H3
11 / 6,
H 31
408,
所以 cond ( H 3 ) 748.
同样可计算
cond ( H 6 ) 2.9107 , cond ( H 7 ) 9.85108.
当 n愈大时，H n 矩阵病态愈严重.
(2)
考虑方程组
H 3 x (11 / 6,
13 / 12,
47 / 60)T b,
设 H 3及 b 有微小误差(取3位有效数字)有
144
1.00
0.500
0.333
0.333 x1 x1 1.83
0.250 x2 x2 1.08 ,
x x 0.783
0.200
3
3
0.500
0.333
0.250
（5.8）
简记为 ( H 3 H 3 )( x x) b b.
方程组 H 3 x b 与(5.8)的精确解分别为：
x (1, 1, 1,)T ,
x x (1.089512538, 0.487967062, 1.491002798)T .
于是
x (0.0895, 0.5120, 0.4910)T ,
H 3
H3
0.1810 3 0.02%,
145
b
b
0.182%,
x
x
51.2%.
这就是说 H 3与 b 相对误差不超过 30%，而引起解的相
对误差超过 50%.
由上讨论，要判别一个矩阵是否病态需要计算条件数
cond ( A) A1 A .
而计算 A1 是比较费劲的，最好在计算中能发现病态
情况.
(1) 如果在 A 的三角约化时(尤其是用主元素消去法解
146
Ax b 时)出现小主元，对大多数矩阵来说，A是病态矩阵.
(2)
系数矩阵的行列式值相对说很小，或系数矩阵某
些行近似线性相关，这时 A可能病态.
(3)
系数矩阵 A元素间数量级相差很大，并且无一定
规则， A 可能病态.
用选主元素的消去法不能解决病态问题.
对于病态方程组可采用高精度的算术运算(采用双倍
字长进行运算),
或者采用预处理方法.
即将求解 Ax b 转化为一等价方程组
PAQy Pb;
y Q 1 x.
147
选择非奇异矩阵 P, Q使
cond ( PAQ) cond ( A).
一般选择 P, Q 为对角阵或者三角矩阵.
当矩阵 A 的元素大小不均时，对 A 的行(或列)引进适当
的比例因子(使矩阵 A 的所有行或列按∞-范数大体上有相同
的长度，使 A 的系数均衡)，对 A 的条件数是有影响的.
但这种方法不能保证 A 的条件数一定得到改善.
148
例10
设
1 10 4 x1 10 4
,
1
1 x2 2
（5.9）
计算 cond ( A) .
1 10 4
,
A
1
1
4
1
10
1
1
A 4
10 1 1
1
(110 4 ) 2
4
cond ( A)
10
.
4
10 1
这个矩阵的元素大小不均，在 A的第一行引进比例因子.
a1i 104 除第一个方程式，
如用 s1 max
1i 2
149
得 Ax b ，即
10 4
1
1 x1 1
,
1 x2 2
（5.10）
而
1
1
1
,
( A)
4
4
1 4 10
10
1 10 4 1
x10
1
,
1
1 x2 2
1
于是
cond ( A)
（5.9）
1
4.
4
1 10
当用列主元消去法解(5.9)时(计算到三位数字)，
150
1
( A b)
0
10 4
10 4
于是得到很坏的结果 :x2 1,
10 4
,
4
10
x1 0.
用列主元消去法解(5.10)，得到
1 2
1
1 1 2
( A b)
10 4 1 1
0 1 1
,
10 4 1 x1 1
,
（5.10）
1
x
2
1
1, x 1.
2 x
从而得到较好的计算解：
2
1
设 x 为方程组 Ax b的近似解，于是可计算 x 的剩余
向量 r b Ax ，当 r 很小时，x 是否为 Ax b 一个较好
的近似解呢?
151
定理23 (事后误差估计)
设 A 为非奇异矩阵， x 是
Ax b 0 的精确解.再设 x 是此方程组的近似解，r b Ax
则
xx
x
证明
cond ( A)
r
b
.
（5.11）
由 x x A1r ，得
（5.12）
x x A1 r ,
又有
b Ax A x ,
A
1
,
x
b
由(5.12)及(5.13)即得到(5.11).
（5.13）
152
(5.11)式说明，近似解 x 的精度(误差界)不仅依赖于剩
余 r 的“大小”，而且依赖于
A
的条件数.
当 A是病态时，即使有很小的剩余 r ，也不能保证 x
是高精度的近似解.
xx
x
cond ( A)
r
b
.
（5.11）
153
5.5.2
迭代改善法
设 Ax b，其中 A R nn 为非奇异矩阵，且为病态
方程组(但不过分病态).
本节研究的问题是，当求得方程组的近似解 x1后，
如何对其精度进行改善.
首先用选主元三角分解法实现分解计算
PA
LU ,
其中 P为置换阵， L为单位下三角阵， U为上三角阵，
且求得计算解 x1 .
154
现利用 x1 的剩余向量来提高 x1 的精度.
计算剩余向量
r1 b Ax1 ,
求解 Ad r1,
（5.14）
得到的解记为 d1 . 然后改善
x2 x1 d1.
（5.15）
如果(5.14)，(5.15)及解 Ad r1 的计算没有误差，则
Ax2 A( x1 d1 ) Ax1 Ad1
Ax1 r1 b,
说明 x2 就是 Ax b 的精确解.
155
但在实际计算中，由于有舍入误差，x2 只是方程组的
近似解，重复(5.14)，(5.15)的过程，就产生一近似解序列
{xk }，有时可能得到比较好的近似.
算法3 (迭代改善法)
设 Ax b，其中 A R nn 为非
奇异矩阵，且 Ax b为病态方程组(但不过分病态)，用选
主元分解法得到 PA
LU及计算解 x1 .
本算法用迭代改善法提高近似解 x1 的精度.
设计算机字长为 t ，用数组 A(n, n)保存 A元素，数组
C (n, n)保存三角矩阵 L 及 U，用 Ip(n) 记录行交换信息，
156
x(n) 存贮 x1 及 xk , r (n) 保存 rk 或 d k
.
1. 用选主元三角分解实行分解计算 PA LU 且求计
算解 x1 (用单精度)
2. 对于 k 1,2, , N 0
(1) 计算 rk b Axk (用原始 A 及双精度计算)
(2) 求解 LUdk P rk ,
(3) 如果 d k
(4)
/ xk
Ly P rk ,
(用单精度)
Udk y.
即
10t 则输出 k , xk , rk , 停机
改善 xk 1 xk d k (用单精度计算)
157
3.
输出迭代改善方法迭代 N 0次失败信息
当 Ax b不是过分病态时，迭代改善法是比较好的改进
近似解精度的一种方法，当 Ax b 非常病态时，{xk }可能
不收敛.
例11
用迭代改善法解
1.0303
0.99030
0.99030 x1 2.4944
0.95285 x2 2.3988
(记为Ax b) (这里 10, t 5 ，用5位浮点数运算).
158
解
精确解 x* (1.2240,1.2454)T (舍入到小数后
第4位).
cond ( A) A
A1
22000 4000.
首先实现分解计算 PA LU , 且求 x1
1
A
0.9118
0 1.0303
1
0
0.99030
LU ,
0.00099
计算解 x1 (1.2560,1.2121)T .
应用迭代改善法需要用原始矩阵 A且用双倍字长精度
计算剩余向量 r b Ax，其他计算用单精度.
159
计算结果如下表.
x1
r1
1.2560
5.7 10 7
1.2121 3.3715 10 5
d1
x2
r2
0.03220 1.2238 1.18 10 6
0.033502 1.2456 9 10 7
d2
2.285 10 4
2.365 10 4
x3 (1.2240,1.2454)T ,
r3 (0.68210 5 , 0.65910 5 )T ,
d 3 (0.271710 4 , 0.351510 4 )T .
如果 xk 需要更多的数位，迭代可以继续.
160
161