浙江大学研究生学位课程
《实用数值计算方法,1
第三章 线性代数方程组
3.1 问题概述
3.2 直接法
3.3 迭代法
3.4 稀疏矩阵
3.5 其他特殊形式的矩阵浙江大学研究生学位课程
《实用数值计算方法,2
3.1 问题概述
3.1.1 问题提出线性代数方程组
23
13
2211
22222121
11212111
bAx
bxaxaxa
bxaxaxa
bxaxaxa
MNMNMM
NN
NN
用矩阵形式表示：
33,,2
1
2
1
21
22221
11211
mnmnmm
n
n
b
b
b
b
x
x
x
x
aaa
aaa
aaa
A
系数矩阵 未知向量 右顶端浙江大学研究生学位课程
《实用数值计算方法,3
当 M=N时，如果 A非奇异，则方程组（ 3-1）
存在唯一解 。
3.1.2 矩阵的存储与结构
1,存储方式
a.满存方式 N2个实数
b.部分存储方式 非零元素个数稀疏矩阵、对称矩阵、块状矩阵
2,存储结构数组在计算机内存中总是一维存放的但是它的顺序在不同的高级语言中不一定相同。
浙江大学研究生学位课程
《实用数值计算方法,4
3.1.2
例如,FORTRAN语言中的矩阵是按列存放的：
,11a,21a?,1na,12a,22a?,2na,1na,2na?,nna
3,存储方式与存储结构是不同的两个概念
4,矩阵的逻辑维数与物理维数逻辑维数,实际参与计算的矩阵阶数物理维数,该矩阵可能出现的最大阶数例如,调用一个子程序，计算一个 4?4矩阵的转置，调用形式为：
CALL MATINV（ A,AI,NL,NP）
这里,NL是逻辑维数,=4。而 NP是物理维数，即数组 A的实际定义维数。
浙江大学研究生学位课程
《实用数值计算方法,5
3.1.2
假定 NP=6，则 4?4矩阵存放于 6?6矩阵中：
7
(1 )
5
(7 )
1
(1 3 )
2
(1 9 )
X
(2 5 )
X
(3 1 )
2
(2 )
3
(8 )
1
(1 4 )
1
(2 0 )
X
(2 6 )
X
(3 2 )
4
(3 )
9
(9 )
0
(1 5 )
5
(2 1 )
X
(2 7 )
X
(3 3 )
8
(4 )
1
(1 0 )
2
(1 6 )
4
(2 2 )
X
(2 8 )
X
(3 4 )
X
(5 )
X
(1 1 )
X
(1 7 )
X
(2 3 )
X
(2 9 )
X
(3 5 )
X
(6 )
X
(1 2 )
X
(1 8 )
X
(2 4 )
X
(3 0 )
X
(3 6 )
NL
NP
NL=4,NP=6 NP,NP
Dimension A(6,6)
A内存放 4?4矩阵
NL?NL求 A( i,j ) 1? i,j?4
但地址是：
（ j-1）?NP+i
浙江大学研究生学位课程
《实用数值计算方法,6
3.1.3 向量范数、矩阵范数定义 1,对于任一向量是,按照一定规则确定一个实数与它对应。 该实数记为,
若 满足下面三个性质：
那么实数 称为向量 x的范数。
设,则它常用的几种范数有：
nRx?
。对
，对任意实数
，
yxyxRy
xx
xxx
n
,3;2;0001
x
,,,,21 Tnxxxx
x
x
浙江大学研究生学位课程
《实用数值计算方法,7
3.1.3
63,,,m a x
53
43
21
1
211
2
1
1
222
2
2
12
n
n
i
in
n
i
in
xxxx
xxxxx
xxxxx
可以验证，以上定义的几种范数均满足三个范数的性质。它们的几何意义见图 3.1
从向量范数出发可以定义矩阵范数。
定义 2：设 A为 n?n阶矩阵。定义
Ax
x
Ax
A
nn Rx
x
Rx
x
10
m a xm a x
浙江大学研究生学位课程
《实用数值计算方法,8
3.1.3
5.5
12
4.7
1
2
x
x
x
321,,xxxx
2x
1x
3x
0.2
5.5
5.4
图 3.1 向量范数的几何意义浙江大学研究生学位课程
《实用数值计算方法,9
3.1.3
这样定义的矩阵范数具有性质：
。对有对
，对
BABARRB
xAAxRx
BABARRB
AAR
AAA
nn
n
nn
,5;,4;,3;2;00,01
显然，这样定义的矩阵范数与向量范数的定义方法有关。
前面三种常用的向量范数相应的矩阵范数是
的最大特征值是 AA
A
T
1
12 73
浙江大学研究生学位课程
《实用数值计算方法,10
3.1.3
93m a x
83m a x
1
1
1
1
1
n
j
ij
ni
n
i
ij
nj
aA
aA
另外，关于范数有一个很重要的等价定理：
定理：有穷线性空间上的一切范数都是等价的。即对任意两种范数有关系式：
,
有关的常数，是与其中
Mm
Mm
,
浙江大学研究生学位课程
《实用数值计算方法,11
3.1.4 线性代数方程组的性态线性方程组 (3-2) 的解完全由 A和 b确定。
在实际问题中,
由于各种原因,A和 b是有误差研究 A和 b的微小摄动对解 x的影响十分重要的这种影响的大小反映了问题的“性态”。
在（ 3-2）中，如果 A-1存在，则解可表示为
x=A-1b
又设 为 A,b的微小摄动，而 是由此而使 x产生的误差。即
bA,x?
bbAAAx 1
浙江大学研究生学位课程
《实用数值计算方法,12
3.1.3
示例 在 4位字长的计算机上解方程组
的。
）是“坏条件”象（微小摄动而得。可以想通过只是则方程组矛盾。事实上改成：如果把误差很大作为主元的消元法解得而取此方程组的真解是
103
103113
113
14 7.537 4.100 0.1
41.1512 2.400 0.3
103
50 0.3,95 1.9
:00 0.3
2.6,66 58.13
:
103
14 7.537 4.100 0.1
41.1512 7.400 0.3
21
21
21
11
21
21
21
xx
xx
xx
a
xx
xx
xx
浙江大学研究生学位课程
《实用数值计算方法,13
3.1.3
现给出一般情形的估计式设则可以证明以下估计式其中可见,k(A)近似表示了方程组求解的误差的相对放大率
11AA?
1232
1
b
b
A
AAk
x
x
11
1
AA
AAAk
浙江大学研究生学位课程
《实用数值计算方法,14
3.1.3
换言之,的大小表示问题的病态程度
k(A)称为计算问题（ 3-2）的条件数
Condition Number
现在再看前例：
2.7 8 42
6.1 4 25
501.5m a x
0.6000.200
6.8258.274
347.1000.1
127.4000.3
1
1
1
1
1
1
1
1
1
1
AAAkAC o n d
A
aA
A
A
n
i
ij
nj
1 AA
浙江大学研究生学位课程
《实用数值计算方法,15
3.1.3
可见计算对象（ 3-10）是病态的。
应当指出：
1,det(A)的值小未必会 引起 A病态例如：
det(A)=0.02，而 A是好条件的 。
2,严格来说，估计式只给出了好条件的充分条件，但不是必要条件。
3,由于数据误差在线性系统中引起的固有的不可靠的使得任何过分精度要求的企图都是徒劳的。
2.010.0
0.010.0A
浙江大学研究生学位课程
《实用数值计算方法,16
3.2 直接法
3.2.1 直接三角分解法 （ LU分解）
考虑
AX=b A非奇异在 Gauss 消去法中,每一步消元过程相当于对 A作一次初等变换 。即左乘一个初等变换矩阵 T。
第一步,
133
100
10
1
1
1
31
21
1
n
l
l
l
T 0
浙江大学研究生学位课程
《实用数值计算方法,17
3.2.1
第 k步
)143(
1
1
1
1
1
,1
nk
kk
k
l
l
T
0
0
0
到 k步以后 A化为当 k=n-1时,A化为上三角阵 U，即因此
)153(11 ATTT kk?
)163(121 UATTT nn?
)173(111211 LUUTTTA n?
浙江大学研究生学位课程
《实用数值计算方法,18
3.2.1
其中为下三角矩阵，称（ 3-17）为 A的三角分解。
由 Ti性质，可以知
1
1
1
1
,
,1
1
in
ii
i
l
l
T
)183(111211 nTTTL?
1
1
1
1
121
3231
21
1
1
1
2
1
1
nnnn
n
lll
ll
l
TTTL
0
0
浙江大学研究生学位课程
《实用数值计算方法,19
3.2.1
Gauss 消去法的基本步骤
4
3
5
311
612
294
3
2
1
x
x
x
75.2
5.0
5
5.225.10
55.00
294
:1
4
13,1
4
22
3
2
1
x
x
x
5.1
5.0
5
1000
55.00
294
:2
5.0
25.13
3
2
1
x
x
x
浙江大学研究生学位课程
《实用数值计算方法,20
3.2.1
所作的初等变换为
5.225.10
55.00
294
311
642
294
10
4
1
01
4
2
001
T1 * A = A1
1000
55.00
294
5.225.10
55.00
294
15.025.10
010
001
T2 * A1 = A2 u
例,设记
084
162
122
A
浙江大学研究生学位课程
《实用数值计算方法,21
3.2.1
则，它的 LU分解为,
112
011
001
110
010
001
102
011
001
200
040
122
240
040
122
110
10
1
240
040
122
084
162
122
102
11
1
1
2
1
1
122
11
21
21
TTL
ATA
ATA
uAAA
TT 左乘左乘浙江大学研究生学位课程
《实用数值计算方法,22
3.2.1
如果有了 LU分解。则解方程组就变得非常容易了。因为
bLY
UXY
bL U XAX
LUA
则令由于 L,U都是三角矩阵，则 Y,X 都可以很容易从上面两式中求得。
如果预先知道 A存在 LU分解，则我们可以直接把这种分解求出来。
浙江大学研究生学位课程
《实用数值计算方法,23
3.2.1
15.0105.1
5.05.055.0
05.04925
5.1
5.0
5
1000
55.00
294
5.15.05.25
4
1
4
5.05
4
2
3,5
4
3
5
15.2
4
1
01
4
2
001
1000
55.00
294
15.2
4
1
01
4
2
001
311
642
294
435
3
32
231
3
2
1
3
21
3
2
1
x
xx
xxx
x
x
x
y
yy
y
y
y
UL
A
b
T
解：
解：
前例：
浙江大学研究生学位课程
《实用数值计算方法,24
3.2.1
但是 Gauss 消去法并不是对任何非奇异矩阵都能顺利进行。 如我们有这样的结果：
任何非奇异的 n?n阶矩阵 A，总存在一个排列矩阵 P，使 PA能进行 LU分解。
上述结果 表明,非奇异矩阵总是能进行选主元的 Gauss消去法 。
01
10
A
浙江大学研究生学位课程
《实用数值计算方法,25
3.2.1
331
421
642
100
010
642
010
100
642
331
421
642
行变换选主元消元消去法：选主元的 G aus s
021 42?
浙江大学研究生学位课程
《实用数值计算方法,26
3.2.1
为使 LU 分解标准化。把 U 改换成 DU
是可行的。其中 D 是非奇异对角矩阵。而
U 则为单位上三角矩阵。如
100
010
2
111
200
040
002
200
040
122
称 A=LDU 为矩阵 A 的 LDU 分解。
定理,设
nk
aaa
aaa
A
nnaA
kkkk
k
k
ij
,,2,1
21
11211
阶矩阵。是一个浙江大学研究生学位课程
《实用数值计算方法,27
3.2.1
表示 A 的各阶主子矩阵。那么,A 存在唯一的 LDU 分解的 充分必要条件 是：
Ak （ k=1,2,….,n ） 都是非奇异。
如 A 为对称正定，或者对角占优，则
A 的各阶主子矩阵均非奇异。
现在我们直接从 A 求 LU 分解：
设 矩阵 A 有 LDU 分解。记 LD=L。 则
A=LU,L为下三角矩阵,U为单位上三角。
这种 LU分解称为 Crout 分解 。
浙江大学研究生学位课程
《实用数值计算方法,28
3.2.1
LUA
u
uu
uuu
llll
lll
ll
l
aaaa
aaaa
aaaa
aaaa
jn
nj
nj
nnnjnn
ijii
nnnjnn
inijii
nj
nj
1000
100
10
1
0
00
000
22
1112
21
21
2221
11
21
21
222221
111211
浙江大学研究生学位课程
《实用数值计算方法,29
3.2.1
行已得到。的前列和的前假定故时，
时，特别，
所以，
因为，
显然，
令：
11
)213(,,2
1
203,,2,1
1
193,,2,1,
,,2,1,1
,
11
1
1
1111
11111
,m i n
1
kukL
nj
l
a
u
ula
i
nilula
j
njiula
LUA
njiu
uUlL
j
j
jj
iii
ji
r
rjirij
ii
ijij
浙江大学研究生学位课程
《实用数值计算方法,30
3.2.1
233,,1
223,,
,,
193,1
1
1
1
1
1
1
nkj
lulau
nkiulal
nkiulla
u
kk
k
r
rjkrkjkj
k
r
rkirikik
k
r
rkirikik
kk
类似地，有，
故有所以根据由于关于 L,U元素的存放从（ 3-22），（ 3-23）可以看出。 A的元素 aij在计算出 或 以后就不再有用了 。
ijl iju
浙江大学研究生学位课程
《实用数值计算方法,31
3.2.1
故 L,U的非零元素便可以存放在矩阵 A中的相应位置上了。
最后,A所存放的元素是：
nnnn
n
n
lll
ull
uul
21
22221
11211
容易证明：
1,即是 Gauss消去法 中各次的主元素。 nn
ll,,11?
浙江大学研究生学位课程
《实用数值计算方法,32
3.2.1
2,如果 A 的各阶主子矩阵均非奇异，
上述过程可以一直进行到底。
矩阵 A 除了以上介绍的 LU 分解以外，
还有一种重要的分解方法可以用于求解线性方程组。即 QR 分解。
定理,任意矩阵 A，总是以分解成正交矩阵 Q 与上三角矩阵 R 的乘积：
如果 A 是非奇异的，则在规定 R的对角元素的符号下，分解式是唯一的。
如果,则 Q 是正交矩阵
RQA
TQQ1
浙江大学研究生学位课程
《实用数值计算方法,33
3.2.1
常见的有 Householder 正交三角分解有了正交三角分解以后，解方程组就变得非常简单：
bQRX
bQ RX
bAX
T?
用 Householder 变换进行分解。从而求解线性方程组的过程是非常稳定的。也不必选主元。特别在方程性质不太好时，
更显其优越性，但计算量大。
Householder 正交三角分解还有其他用途。
浙江大学研究生学位课程
《实用数值计算方法,34
3.2.1
PROGRAM EXAMPLE
DIMENSION A(5,5),B(3),IND(3),C(3)
DATA B /5.0,3.0,4.0/
A(1,1)=4.0
A(3,3)=3.0
N=3
NP=5
CALL LUDCMP(A,N,NP,IND,D)
CALL LUBKSB(A,N,NP,IND,B)
WRITE(*,*) B
C(1)=1.
C(2)=2.
C(3)=3.
CALL LUBKSB(A,N,NP,IND,C)
WRITE (*,*) C
STOP
END
浙江大学研究生学位课程
《实用数值计算方法,35
3.2.1
15.0
5.0
05.0
4
3
5
00000
00000
00105.225.0
0055.05.0
00294
00000
00000
00311
00642
00294
L U B K S B
L U B K S B
L U D C M P
B
A
的变化：
前同的变化：
3
2
1
0.1DIN D 的值：
浙江大学研究生学位课程
《实用数值计算方法,36
3.2.1
程序手稿原程序文件目标程序执行程序输出结果编辑（ WS）
编译（ F77）
连结（ LINK）
执行（程序名）
库修改图 3.2 程序调试流程图浙江大学研究生学位课程
《实用数值计算方法,37
3.2.2 迭代改进由前面讨论，我们知道，由于各种原因，
特别是方程组本身的病态，将会导致解与真解之间的误差不能达到令人满意的程度。有各种处理病态方程组的方法。这里介绍一种简单实用的方法。它还可以用于一般方程组的高精度求解。 该方法的基本思想是利用迭代逐步改善解的精度 （关键在于在迭代过程中有些运算需用双精度进行）
迭代改善过程的框图如下：
浙江大学研究生学位课程
《实用数值计算方法,38
3.2.2
LUAA?的三角分解作
bxL U xx 011,求三角方程组的解
AXb用双精度计算残向量：
L UZZ,求出解的改善量
Zxx 1求改善解：
Z
1
1 xx
xx1?
结束是否图 3.3 迭代改善过程浙江大学研究生学位课程
《实用数值计算方法,39
3.2.2
由上图可知，迭代改进的计算，只要在开始时作一次三角分解。以后只是反复求解三角方程组就可以获得解的改善。
关于解的改善的程度有以下结论：
假设在浮点数的尾数是 t 位（二进制）的计算机上用消去法计算。 A 的条件数为一般可设则：经过一次迭代，解的精度近似地改善了 t-p 位。因此迭代次数 k 满足 k(t-p)?t
就够了。 或近似地用估算式：
pACo n d 2? tp0
1
1
2lo g Z
xtk
浙江大学研究生学位课程
《实用数值计算方法,40
3.2.3 奇异值分解 （ SVD）
解方程组时，常常会出现奇异或十分接近于奇异的情况。用通常的 Gauss消去法或者 LDU分解法难以得到满意的解。对这类问题，我们有一个有效的方法，它叫 奇异值分解法 Singular Value Decomposition
它不但可以诊断出问题所在，而且有时能给出一个有用的数值结果。
SVD方法基于下面这个线性代数定理，
由于定理本身不是我们研究的重点，这里不加证明。
浙江大学研究生学位课程
《实用数值计算方法,41
3.2.3
定理,对任意 M?N矩阵 A,这里 M?N，
它都可以分解成下列形式：
A=UWVT
其中,U是 M?N列正交矩阵；
V是 N?N 正交矩阵；
W是 N?N对角矩阵。且对角元素非负。
即：
Nlk
uu
V
w
w
w
UA
kl
M
i
ilik
T
n
,1
1
2
1
且
浙江大学研究生学位课程
《实用数值计算方法,42
3.2.3
IVVUU
Nlkvv
TT
N
j
kljljk
用矩阵表示：
,1,
1
由于 V是方阵，故它还是行正交即不管矩阵 A是否奇异。 SVD分解都能做到。而且它“几乎”是唯一的。即允许
U,V的列进行变换并进行 W相应交换的前提下是唯一的。
如果 A是 N?N方阵，则 U,V,W都是方阵。且分解式可写为：
.IVV T?
浙江大学研究生学位课程
《实用数值计算方法,43
3.2.3
T
j
T
j
U
w
di agVA
Vwdi agUA
11
由此但是，如果其中有一些 或十分接近于 0，即 A奇异或接近奇异，那么 A-1
就难以求得。因此,SVD分解可以通过判断 的大小来分析 A的性态。
当 A是奇异时,我们就考虑
0?jw
jw
jj wwAC o n d m inm a x?
243 bAx nn
在某种意义下的解。
浙江大学研究生学位课程
《实用数值计算方法,44
3.2.3
由代数知识，有当 rank(A)<rank(A,b)时,
（ 3-24）为矛盾方程组。无解。
当 rank(A)=rank(A,b)时,
（ 3-24）有无穷多解。
在这些解中。有一个 最小长度解在某些场合是有意义的。
在求这个解以前，先介绍几个重要的概念。
1.（ 3-24）定义了一个从向量空间 x到向量空间 b的线性映射 A。
2,如果Ａ奇异，则有,使得
m in )( 2?x
0?x
浙江大学研究生学位课程
《实用数值计算方法,45
3.2.3
0?AX253?
所有满足（ 3-25）的 x全体称 Ａ的零空间 。
零空间的维数称为零度 。
３,所有存在原象 x,使得（ 3-24）成立的向量
b的全体称为 Ａ的域 。 域的维数称为Ａ的秩 。
４,如果Ａ非奇异，则零空间只有零向量；
零度＝０；Ａ的秩＝Ｎ；
如果Ａ奇异，则零度 >0; 秩 <N;
但可以证明,零度＋秩＝Ｎ 。
现在再看看ＳＶＤ分解的意义。
浙江大学研究生学位课程
《实用数值计算方法,46
3.2.3
ＳＶＤ分解事实上分别构造了Ａ的零空间和域的各一组正交基。
由所有 对应的Ｕ的第 j列向量是张成域的一组正交基 。
由所有 对应的Ｖ的第 j列向量是张成零空间的一组正交组 。
现在考虑（ 3-24）的解。如果 的域，
则方程有无穷解。因为零空间的任何向量的
0?jw
0?jw
Ab?
bxV
w
w
w
uuuu
T
n
jnj
1
21
,,,,,
浙江大学研究生学位课程
《实用数值计算方法,47
3.2.3
buxvw
uxvw
xvw
xvw
xvw
xvw
uuu
xV
w
w
w
w
uuuxA
k
w
nk
k
T
kk
n
i
i
T
ii
T
nn
T
jj
T
T
nj
T
n
j
n
0
1
1
22
11
1
2
1
21
,,,,
,,,
故任意的域中向量 b可以由 线性表示，
即这样的 张成域的一组基。
0?kk wu
0?kk wu
（ 1）
浙江大学研究生学位课程
《实用数值计算方法,48
0
1
2211
2
1
2
1
,0,0
),,2,1(
0
0
,,2,10
0
0
0
0
0
k
w
nk
kk
jjj
nn
k
j
T
j
j
T
kk
T
n
T
T
n
T
T
vx
vjw
vvvX
nkvx
vxxv
w
nkxvw
xv
xv
xv
w
w
w
XWV
UU W VA
Ax
Ax
于是内积即得上式两边与的下标对应于的线性组合时表示成由此若令向量正交与即时当故即是非奇异阵，得及由的零空间，即设
3.2.3
（ 2）
浙江大学研究生学位课程
《实用数值计算方法,49
3.2.3
xxxxxx
bbUU
bU
w
di agUW
bU
w
V di agU W V
xU W VAx
w
wbU
w
di agVx
bAxxx
bbAxxAxxA
bxAx
xAx
T
T
j
T
j
T
T
j
j
T
j
现在要证明：。之和：的向量及零空间中表示为特解因为任何解向量都可以的所有解向量中长度最小式还是事实上，
的解。式是问题故则时，取当令的解仍然是方程所以的解时，是当是零空间中的向量。即设解上，仍是方程的解：解的线性组合加在某个
''
'
''
''
)243(
243
1
1
0
1
0,
1
0
,0
0,1 0,0{ jjjj ww令对角阵
0jw ijubAb?的域，是属于因为浙江大学研究生学位课程
《实用数值计算方法,50
3.2.3
的解。是所有解中，长度最小中的特解的长度最小；或者说，，换言之，当
，故综合以上两点：
时，自然当上式长度更小，又才能使分量时，最好取这些下标所以取分量时，由于当因为
x
xx
xxV
xvw
xvj
bu
w
w
xv
xv
xv
bu
w
bu
w
bu
w
xVbUW
xVxV
xxVxxx
T
T
jj
T
j
T
j
j
j
n
T
n
T
T
T
n
n
T
T
TT
TT
T
0
00
00
,0
,0
1
0
1
1
1
'
''
'
'
'
'
22
'
11
2
2
1
1
'1
'
''
浙江大学研究生学位课程
《实用数值计算方法,51
3.2.3
的大小控制。误差值被“零化”；多少大小的奇异值可以较小。但应注意掌握倒可能会使问题求最小长度解此时，把方程作为奇异可能很大差分解直接求解将会使误分解或用时，
件数对病态方程组，即使条最小的解。
即使向量长度最小二乘意义下的解，
式给出了无解，但的域时，当说明：
r
r
AXbr
SV DLU
w
w
Axbr
Ab
i
i
2
1
:
1
m i n
m a x
.2
243.1
浙江大学研究生学位课程
《实用数值计算方法,52
3.2.3
以上讨论用图表示如下：
x
b
A
bAx?
零空间的解dAx?
解的 SVDdAx?
的域A
d
'c
的解'cAx?
最小二乘解的 SV DcAx?零空间?c
图 3.4 解空间示意图浙江大学研究生学位课程
《实用数值计算方法,53
3.2.3
T
v
T
T
VWU
U W VASV DA
bAX
bAr ankAr ank
bA
bAX
3
3
1
6
2
0
3
1
6
1
2
1
3
1
6
1
2
1
000
030
002
3
2
23
1
2
1
3
1
23
4
0
3
2
23
1
2
1
,
3,,2
0
0
1
,
23
1
3
1
6
1
6
1
32
1
2
3
22
3
2
3
2
23
1
3
1
6
1
6
1
32
1
2
，，
其中分解：进行把最小长度解。现求最小二乘意义下的为矛盾方程组，无解所以因为，
是矛盾方程组。例：下列线性代数方程
1u 2u
浙江大学研究生学位课程
《实用数值计算方法,54
3.2.3
32
1
9
1
9
1
36
1
22
1
0
3
1
2
1
2
01
23
1
23
4
23
1
2
1
0
2
1
2
3
1
6
2
0
1
2
2211
3
3
21
3
bUVx
ukukAx
Rx
kvA
uu
v
Ar ankA
T
：原方程组的最小长度解即有
，，且域的一组基为：域的维数为
，且零空间的一组基为零度为奇异，且由此，
浙江大学研究生学位课程
《实用数值计算方法,55
3.3 迭代法迭代法适用于解高阶稀疏非病态方程组。它只需要存储非零元素。
但对有些问题，迭代可能发散或收敛很慢。
3.3.1 Jacobi 迭代法设有方程组其中 A是 N?N非奇异矩阵。故有唯一解。
把（ 3-26）改写成迭代形式
bAX263?
浙江大学研究生学位课程
《实用数值计算方法,56
273 gBXX
于是设 X(0)是一个任意的初始迭代向量。
代入（ 3-27）的右边，得到新的向量记为 X(1)：
一般地有：
我们得到向量序列如果 收敛于 X,即：
则向量 X 是（ 3-27）的解。
现讨论解的求法及收敛性
gBXX )0()1(
)283(,2,11)(kgBXX kk
,1,0,?kX k
kX
XX k
klim
3.3.1
浙江大学研究生学位课程
《实用数值计算方法,57
3.3.1
设矩阵 A的对角元素 。记则 D非奇异。将 A表示成和式其中：
0?iia
nnaaad ia gD,,,2211
293 uL CCDA
方法的定义是：唯一确定。述矩阵迭代法都可以用上我们会看到，
J ac obiCCD
SO RSGJ ac obi
a
aa
aaa
C
aa
aa
a
C
uL
n
n
n
u
nn
L
,,
,,
0000
000
00
0
,
0
00
000
0000
3
223
11312
21
3231
21
浙江大学研究生学位课程
《实用数值计算方法,58
303,,2,1,
1
1
nibxaxa
n
ij
j
i
k
jij
k
iii?
3.3.1
利用（ 3-29）式的矩阵符号。 Jacobi 法可以表示成：
这里 B 是 Jacobi 迭代矩阵。其定义为
gBxx kk 1
333,,,
323,,,
313
21
1
21
11
T
n
Tk
n
kkk
uL
bbbDg
xxxx
ADICCDB
Jacobi 方法收敛的充要条件是：
迭代矩阵 B的谱半径 1)(?B?
浙江大学研究生学位课程
《实用数值计算方法,59
3.3.1
利用这个判断标准。
，故较多的计算可以较方便检验由这是由于法收敛方法的一个充分条件：给出为此，的条件较难检验由于的特征值。是这里模最大者：谱半径定义为特征值的
B
BB
J ac obiB
J ac obi
B
Bni
B
i
i
ni
1
,1
,,2,1
m a x
1
1
,,m a x
i
iiiB
即为设单位长度的特征向量
B
BBB
B
i
iiiii
iii
于是则浙江大学研究生学位课程
《实用数值计算方法,60
3.3.1
现在讨论一下 Jacobi方法的收敛速度的快慢及影响速度的因素。
的精确解。是方程组其中，
方法成立，则定理：若
273
353
1
343
1
1
1
01
x
xx
B
B
xx
xx
B
B
xx
J ac obiB
kkk
k
k
B
BI
1
11
示：其证明作为思考题，提浙江大学研究生学位课程
《实用数值计算方法,61
3.3.1
（ 3-34）式可作为误差估计式，并且说明 ||B||越小,x(k) 收敛越快。（ 3-35）式说明了只要 ||B||不是很小，当相邻两次迭代向量 x(k) 和 x(k-1) 很接近时，解 x* 和 x(k) 也是很接近的。因此，可以用 || x(k)- x(k-1)|| 作为迭代终止的判别式。这个结论对逐次逼近法也成立。
特别，应该指出，以上讨论的收敛性及收敛速度均是对一般迭代式（ 3-27）而言。因此，其结论对以后要讨论的 G--S迭代法和 SOR法也成立。
浙江大学研究生学位课程
《实用数值计算方法,62
3.3.1
Jacobi 迭代法的计算步骤：
)1
,2,13;2;1
1
1
0
11
akk
xxb
gBxxa
k
x
bDgADIB
kk
kk
，转则结束，否则若计算对和初始迭代向量给出精度小量
，计算
浙江大学研究生学位课程
《实用数值计算方法,63
n
k
n
k
n
k
n
k
nn
kk
k
nn
kk
k
i
k
j
k
j
k
i
k
i
n
k
n
k
n
k
n
k
nn
kk
k
nn
kk
gxbxbx
gxbxbx
gxbxbx
xxijx
xx
gxbxbx
gxbxbx
gxbxbx
J ac obi
0
0
0
:1,,1
0
0
0
2211
2
1
21212
1
1
1
1
2121
1
1
1
22
1
11
2
1
2
1
1212
1
1
1
1
2121
来求代替故可考虑用已求出。时，由上可知当计算迭代格式：仔细写出
3.3.2 Gauss--Seidel 迭代法 。
G--S迭代法的基本思想来自 Jacobi
迭代法。只是迭代矩阵 B不同。
浙江大学研究生学位课程
《实用数值计算方法,64
3.3.2
这就是 Gauss--Seidel 迭代格式，它也可写为：
若用（ 3-29）式记号，则 G--S迭代式可表示为：
或其中
),,1(
1
1
1
1
ni
gxbxbx i
n
ij
k
jij
i
j
k
jij
k
i
bxCxCD kukL 1
3631' hxLx kk
UL CDUCDL
bDLIhULIL
11
111 373'
，
，
浙江大学研究生学位课程
《实用数值计算方法,65
有关收敛性和收敛速度的讨论完全与 Jacobi 方法相同，只是迭代阵为 L'。
3.3.3 逐次超松弛（ SOR） 法
SOR法定义为,
法可写成：的矩阵记号利用法。法退化为时，当
。是实数，称为松弛因子其中
S O R
S e id e lG a us sS O R
nix
gxbxbx
k
i
i
n
ij
k
jij
i
j
k
jij
k
i
,293
1
,,2,11
1
1
1
1
1
浙江大学研究生学位课程
《实用数值计算方法,66
3.3.3
1)(
,
393
1
383
1
11
11
1
1
1
1
L
SO RJ ac obi
CDUCDL
bDLIh
IULIL
hxLx
Dx
bxCxCDx
UL
kk
k
k
L
k
L
k
法收敛的充要条件是法类似，与其中：
或浙江大学研究生学位课程
《实用数值计算方法,67
2.3.13
1)(1
2
:
403
40320
1)(
L
SO R
SO R
A
SO R
L
opt
opt
，可以证明对某些特殊类型的矩阵
。为称为最佳松弛因子，记因子法收敛速度最快的松弛使收敛的充分条件。也是是对称正定矩阵，则另外，若收敛的必要条件是
，因此可以证明曹志浩等求根》参考《矩阵计算与方程
浙江大学研究生学位课程
《实用数值计算方法,68
3.3.3
三种迭代法的比较：
1,一般情况下,J方法与 G--S方法比较并无优劣。收敛情况与速度均不一定。
例：
011
202
110
022
101
220
19.09.0
9.019.0
9.09.01
B
B
A
法收敛方法不收敛
SG
J
法不收敛方法收敛
SGL
JB
,1'
,1
1'
,1
L
B
浙江大学研究生学位课程
《实用数值计算方法,69
3.3.3
2,但是具有相容次序的矩阵，在相同精度要求下，迭代次数分别为：
Jacobi 方法,1154
G--S 方法,578
SOR 方法,59
可见 对这类矩阵,G--S 法比 J方法快一倍，
而 SOR 法的收敛速度可提高一个数量级。
最后介绍一类对于三种方法都收敛的矩阵。 先定义可约矩阵和对角占优矩阵：
opt
浙江大学研究生学位课程
《实用数值计算方法,70
3.3.3
(1),称 n阶方阵为 可约矩阵,如果存在 n阶排列阵 P，使得其中 A11为 r阶方阵,A22为（ n-r）阶方阵不是可约矩阵就是 不可约矩阵
(2),称 n 阶方阵 A是 对角占优矩阵,如果若上式严格不等号成立，则称 A为 严格对角占优矩阵 。
4130
22
1211
A
AAPAP T
4231
1
niaa
n
ij
j
ijii,
浙江大学研究生学位课程
《实用数值计算方法,71
3.3.3
(3),称 A是 不可约对角占优矩阵 。如果 A是不可约的且对角占优，而（ 3-42）式中至少对一个 i 有严格不等号成立。
我们有如下结论,
1,若 A是严格对角占优或不可约对角占优矩阵，则 A非奇异。
2,若线性方程组系数 A是上述矩阵，则：
1） J方程和 G--S方法都收敛。
2）当 SOR方法收敛。10
浙江大学研究生学位课程
《实用数值计算方法,72
3.4 稀疏矩阵前面介绍了 直接法 和 迭代法 。但是实践中如何选择计算方法是一个很复杂的问题。需要考虑的因素很多，主要有：
1.方程组的性质； （病态、对角占优、正定等）
2.相类似问题出现的次数；
3.方程组的阶数；
4.计算机的容量、字长、运算速度等；
5.系数矩阵的结构； （稀疏、三对角等）
6.对结果的精度要求。
浙江大学研究生学位课程
《实用数值计算方法,73
3.4
迭代法 是解超高阶线性代数方程组的重要方法之一，但是它也有缺点：
不适合解病态问题；
精度不高；
收敛速度依赖于加速因子的选取；
收敛性不能保证；
所需乘除运算量大。
在这些方面，相对来说，直接法运算步骤一定，运算次数有限，精度高。对解大型稀疏矩阵问题可以利用和保持系数矩阵的稀疏性来降低存储容量和缩短计算时间。
浙江大学研究生学位课程
《实用数值计算方法,74
3.4.3
但是，为了保持系数矩阵的稀疏性。必须使用特殊的计算方法和高级程序设计技术。
近年来在这方面有较大的发展。其优越性大大超过了迭代法。
在这一节里，不详细介绍计算方法，
只给出一些常用的标准稀疏矩阵形式。同时给出一种一般稀疏矩阵的存储结构。这一节里还给出一个有效的计算公式。
浙江大学研究生学位课程
《实用数值计算方法,75
3.4.1 稀疏矩阵的结构和存储从本质上讲，解稀疏问题的方法与一般 Gausss消去法 LU分解法并没有什么不同，只不过前者更注重于运算过程中对零元素填入量的控制。
稀疏系数矩阵问题的直接法必定依赖于矩阵的稀疏形状，下面列举的是常见的标准的稀疏矩阵。
1.稀疏矩阵的结构浙江大学研究生学位课程
《实用数值计算方法,76
3.4.1
带对角块三角阵块三角阵单边块对角双边块对角
e
a
c d
b
单边块三角
f
图 3.5
浙江大学研究生学位课程
《实用数值计算方法,77
3.4.1
双边带三角单边带对角双边带对角其他其他
k
g
i j
h
图 3.6 稀疏矩阵的形式浙江大学研究生学位课程
《实用数值计算方法,78
3.4.1
2,稀疏矩阵的存储一般矩阵的存储方式有两种：满存和部分存储。稀疏矩阵则采用后一种方式。
稀疏矩阵的种类不同，所采用的存储方法也不一样，但是不论采用什么方法，它的标准有两条：
1）节省存储空间 ；
2）存取方便，即节省存取时间。
存储方法的选择一般要根据矩阵的种类和对程序的要求在上述两个标准之间权衡。
浙江大学研究生学位课程
《实用数值计算方法,79
3.4.1
(1) 带对角矩阵的存储
imjkmjib
mji
a
mmjia
ki
ij
ij
1
0
,0
,
，这里当当为半带宽。称当图 3.7 带对角矩阵的存储方式
nnijaA )( )12()( mnijbB
浙江大学研究生学位课程
《实用数值计算方法,80
3.4.1
(2)块三角型矩阵的存储
H1
i1
iL
il+1
N1
H2 N2
Hl+1
aij
Hl Nl
i
1
i
2
i
l - 1
i
l
N
1
N
2
N
L
bb
ID2?L
b
1
b
2
b
M
1
2
1
11
ll
l
lll
HHN
iiH
11 lM NMB
1
1
0
,
0,
llk
l
ij
ll
ijnnijnn
ijib
ij
a
iii
ijaaA
则设时当
1
2
1
1
ll
l
h
ll
h
ijii
iiii
Nk
这里：
iL
图 3.8 块三角型矩阵的存储形式浙江大学研究生学位课程
《实用数值计算方法,81
3.4.1
(3) 一般稀疏矩阵的存储
B,一维存储非零元素（按行）
（ m）
（ m）
JA,B中相应元素所在的列号 （ A中）
IA,每行第一个非零元素在 B中的下标。
k
ij b
elsea 0
jkJA
iIAkiIA
k
）
）
使得如果
2
11
:,
（ n）
浙江大学研究生学位课程
《实用数值计算方法,82
3.4.1
例如：
8706
0050
4003
0210
A
1 2 3 4 5 6 7 8
2 3 1 4 2 1 3 4
1 3 5 6
IA
JA
B
图 3.9 矩阵 A非零元素存放 的一般形式浙江大学研究生学位课程
《实用数值计算方法,83
3.4.1
(4) 对称矩阵的存储因为故只要存储即可。
任何形状的矩阵，包括稀疏的和非稀疏的，只要是对称的。
均可节省将近一半的存储量。
njiaa jiij,,2,1,
njiija ij,,2,1,,
浙江大学研究生学位课程
《实用数值计算方法,84
3.4.2 Sherman-Morrison & Woodburg 公式我们知道，矩阵求逆问题等价于解一个线性方程组。如果我们已经得到矩阵 A 的逆矩阵 A-1。 现在想在 A上作个小小变动。如：一个元素、一行、一列，
或者部分元素，那么是否可以不再做复杂的求逆运算，而只是在原来的 A-1的基础上作些适当的运算就可得到变动后的矩阵的逆呢？回答是肯定的。
假定矩阵变化为：
433 VUAA
浙江大学研究生学位课程
《实用数值计算方法,85
3.4.2
附带复习一下内积，外积的定义：
T
nnnn
n
n
jiij
nn
T
n
i
ii
T
nn
uv
vuvuvu
vuvuvu
vuvuvu
vu
njivuA
Auvvu
uvvuvuvu
v
v
v
v
u
u
u
u
21
22212
12111
1
2
1
2
1
,,2,1,
,,
,
即其中：
外积：
内积：
内积、外积均有结合律
（叉积）
（点积） (或,u?v)
浙江大学研究生学位课程
《实用数值计算方法,86
3.4.2
这里 u?v 是外积，它表示 A的变动。
Sherman--Morrison公式给出了一个计算的方法：
1 vuA
443
1
1
11
1
2111
1111
1
1
11
AvuA
A
AvuAA
AvuAvuAvuAI
AvuAIvuA
这里从（ 3-44）（ 3-45）式可以看出，给了
A-1和 u,v只需作矩阵的乘积运算和一个向量的加法：
4531 uAv?
浙江大学研究生学位课程
《实用数值计算方法,87
3.4.2
473
1
463,,
11
11
Wz
AA
zvvAWuAz
T
则：
令：
Sherman--Morrison公式可以直接应用于一类稀疏线性方程组上。如果我们已经有一个有效的算法得到矩阵 A的逆（如 A
是三对角阵。或者其他标准稀疏形状的矩阵），那么要解 A的稍作变动的矩阵决定的方程组，只要运用一次或多次
Sherman--Morrison 公式即可。
但对有些稀疏问题，这个公式并不能直接运用。理由很简单，因为 A-1不可能在机器内部全部都储存着。在实际运用中这个公式还有变化。
A
浙江大学研究生学位课程
《实用数值计算方法,88
3.4.2
变化一,求解求解线性系统问题：
首先可以对 A用有效的方法解两个辅助问题：
得到了向量 y和 z后，（ 3-48）的解既是：
483 bXvuA
493, uAzbAy
5031
z
zv
yvyx
浙江大学研究生学位课程
《实用数值计算方法,89
3.4.2
z
zv
yv
yx
uAzbAy
uAzbAy
zbAv
bzWbzWbWz
zv
bWz
bA
b
zv
Wz
A
bvuAX
TT
1
,
,
1
1
11
1
1
1
1
有时故当由于证明浙江大学研究生学位课程
《实用数值计算方法,90
3.4.18
变化二,即 Woodbury 公式当在标准稀疏矩阵上添加的元素不只是一行或一列时，我们就不能简单地重复上述过程。因为新的 A-1 并没有存储下来。
这时要用 Woodbury 公式：
.,
513
111111
NPPNVUNNA
AVUAVIUAAVUA TTT
矩阵，是阵，是这里，
（ 3-51）式主要是计算因此可以看到，（ 3-51）式把一个 N?N阵求逆问题转化为一个 P? P阵的有求逆问题。由于 P<<N,故问题得以大大简化。
11 )( PPT UAVI
浙江大学研究生学位课程
《实用数值计算方法,91
3.4.2
。的右边式子也同样计算左乘用均为单位阵即可：
右乘和左乘的右边用只要证公式：证明
513
513
11
1
1
111
1
1
111
1
1
111
T
TT
TTTT
TTTT
TTT
T
UVA
IUVAUVAI
VUAVIUAVIUAUVAI
UVAIVUAVIUAUVAI
UVAAVUAVIUAA
b
VUA
W ood bur y
Woodbury公式和连续地运用 Sherman--
Morrison公式的联系在于：
如果记
PP vvVuuuU,,,,,,121
523
1
P
k
kk
T VUAVUA
则有浙江大学研究生学位课程
《实用数值计算方法,92
3.4.2
上式表明，如果已经求得 A-1，计算的方法有两种：
（ 1）利用（ 3-51），计算一个 P?P 的逆矩阵。
（ 2）利用（ 3-47），执行 P次。
如果 A-1没有存储，则解下列问题的步骤是
(1) 解 P个辅助方程组：
得矩阵
1 TUVA
533
1
bxvuA
P
k
kk
543,,2,1 PiuAz ii?
PzzzZ,,,21
浙江大学研究生学位课程
《实用数值计算方法,93
3.4.21
(2) 计算 P?P矩阵的逆
(3) 再解辅助方程组
(4) 得方程组（ 3-53）的解因为：
5531ZVIH T
563 bAy
573 yVHZyx T
bAVHUAbA
bAVUAVIUAA
bVUAx
T
TT
T
111
11111
1
yVHZy T
浙江大学研究生学位课程
《实用数值计算方法,94
2.5 特殊矩阵在很多问题中，我们会遇到特殊类型的矩阵。对于这些具有特殊性质的矩阵，
不仅存储方式上应予研究并很好利用其性质以节省存储，在计算方法上也要利用它们的好的性质以减少运算量和提高算法稳定性。
在这一节里，我们研究几种特殊矩阵，
对称、正定、带状，和三对角矩阵。
由第二节的定理知，若的各阶主子式
nnijaA
。则 L D UA,0
浙江大学研究生学位课程
《实用数值计算方法,95
3.5
同样我们可以证明，若对称矩阵Ａ的各阶主子式 其中
Ｌ是单位下三角阵，Ｄ为对角矩阵。
这个分解叫做对称矩阵的 分解。
对一般非奇异对称矩阵，则存在排列矩阵Ｐ，使故可通过选主元方式达到 分解。
解线性问题 可转化为解三个简单的线性问题：
。则 TL D LA,0
TLDL
TLDL
TLDLPA?
bAX?
zxL
yDz
PbLy
T?
浙江大学研究生学位课程
《实用数值计算方法,96
3.5
如果Ａ是正定矩阵，则Ａ是对称的并且Ａ的各阶主子式不等于０。因此，
存在且故可令：
从而这里Ｇ是下三角矩阵。这种分解叫做对称正定矩阵的平方根分解。
583 TL D LA
niddddd ia gD iinn,10,,,,2211
nndddd ia gD,,,~ 2211
593
~~~~
T
TTT
GG
DLDLLDDLLD LA
浙江大学研究生学位课程
《实用数值计算方法,97
3.5
ijgga
g
g
gggg
ijgga
gag
ggga
G
ggg
gg
g
aaa
aaa
aaa
GGA
j
k
jkikij
jj
ij
jjij
j
k
jkik
j
k
jkikij
i
k
ikiiii
ii
i
k
ik
i
k
ikii
T
nnnnnnnn
n
n
T
1
1
1
1
1
2
1
1
1
2
2
1
1
2
1
2
21
2221
11
21
22221
11211
1
而，
故，
即：
因为
浙江大学研究生学位课程
《实用数值计算方法,98
3.5
平方根分解的算法如下：
个算法是稳定的。能够得到控制，因此这量的平方根分解法的中间知，由所以的元素有如下关系的对角元素和因为执行对
ij
iiij
i
j
ijii
j
k
jkikijjjijij
i
k
ikiiiiii
gA
ijag
ga
GA
ijggaggaii
gagai
ni
613
613
1,,1)
603)
,,2,1
1
2
1
1
1
1
1
2
浙江大学研究生学位课程
《实用数值计算方法,99
3.5
这个算法的另一个优越性是 不用选主元 。
但是从（ 3-60）知道，这个算法在求要进行平方根运算。
为了避免开方运算，有时我们仍要用
ＬＤＬ Ｔ 分解。
在 3.5 中，引入了带宽为 2m+1的带状矩阵。他们满足：
当带状矩阵当 m<<n 时才有意义 。
带的宽度依赖于方程组未知量的次序和 方程的顺序。
nig ii,1?
0 ijamji 时，
浙江大学研究生学位课程
《实用数值计算方法,100
3.5
例如,m=1 的三对角方程的矩阵
623
1
11121
343332
232221
1211
nnnn
nnnnnn
aa
aaa
aaa
aaa
aa
A
若Ａ的相邻任意两行交换一下，则带宽从
３增加到５。
对带状矩阵，有以下定理：
633,,1,0
,,1,0
,,
12
mjjjiu
jmjmjil
uUlL
m
ij
ij
ijij
仅对仅对三角带状矩阵，即：
是则有ＬＵ分解的带状矩阵定理：如果带宽为浙江大学研究生学位课程
《实用数值计算方法,101
3.5
我们在上述定理中假定ＬＵ的分解存在。事实上未必总是如此，可能需要交换方程。这就可能要增加带的宽度。但不可能增大２倍。
但有些实际问题中遇到的带型矩阵性质是比较好的，并不需要选主元。
三对角矩阵（ m=1） 在带状矩阵中占有重要的位置。若（ 3-62）中的各阶主子式
０，则可用 追赶法 很方便求解 Ax＝ b。
所幸的是实际应用中，有很大一类矩阵是对角占优矩阵,因而各主子式?０。
浙江大学研究生学位课程
《实用数值计算方法,102
第三章 习题
分解法求解。分解及用态；试分析以上方程组的性设有线性代数方程组２
。有极小值，并求极小值取何值时，，求：当
１１
２
设Ａ＝１
SV DLU
xxxx
xxxx
xxxx
xxxx
AC ond
2
1
9
1
7
1
6
1
5
1
4
1
8
1
6
1
5
1
4
1
3
1
7
1
5
1
4
1
3
1
2
1
6
1
4
1
3
1
2
1
.
.
4321
4321
4321
4321
浙江大学研究生学位课程
《实用数值计算方法,103
次数，并讨论之。
。试比较迭代，精度控制在解迭代法求解。取初始的因子为以及取松弛请分别用设有线性代数方程组所需的乘除法次数。
求解线性代数方程组阶下三角矩阵，试计算是设矩阵
4
54
543
432
321
21
100
46.1
,,
12
02
12
02
12
.4
.3
x
SO R
Se i de lG aus sJ ac obi
xx
xxx
xxx
xxx
xx
bLx
nL
浙江大学研究生学位课程
《实用数值计算方法,104
解，并解释解的含义。时分别计算当取的维数及基。线性变换的零空间及域定义的的性态及由试分析分解式为的其中系数矩阵阶代数方程组：设有
SV D
b
b
AbAx
V
di agW
U
U W VA
SV DA
bAx
T
T
T
T
3334.14714.00
,1112
1
5774.05774.05774.0
8165.04082.04082.0
07071.07071.0
1001
6667.03333.06667.0
2357.09428.02357.0
7071.007071.0
33.5
2
1
浙江大学研究生学位课程
《实用数值计算方法,105
。的逆矩阵非奇异矩阵程组的程序来计算一个性代数方试说明如何利用求解线亦是对称正定。是对称正定，则若亦是对称矩阵试证明：
约化为消去法一步，
，经为对称矩阵，且设
1
2
2
2
12111
11
.7
2
1
0
0.6
AA
AA
A
A
Aa
AA
AG aus s
aA
《实用数值计算方法,1
第三章 线性代数方程组
3.1 问题概述
3.2 直接法
3.3 迭代法
3.4 稀疏矩阵
3.5 其他特殊形式的矩阵浙江大学研究生学位课程
《实用数值计算方法,2
3.1 问题概述
3.1.1 问题提出线性代数方程组
23
13
2211
22222121
11212111
bAx
bxaxaxa
bxaxaxa
bxaxaxa
MNMNMM
NN
NN
用矩阵形式表示：
33,,2
1
2
1
21
22221
11211
mnmnmm
n
n
b
b
b
b
x
x
x
x
aaa
aaa
aaa
A
系数矩阵 未知向量 右顶端浙江大学研究生学位课程
《实用数值计算方法,3
当 M=N时，如果 A非奇异，则方程组（ 3-1）
存在唯一解 。
3.1.2 矩阵的存储与结构
1,存储方式
a.满存方式 N2个实数
b.部分存储方式 非零元素个数稀疏矩阵、对称矩阵、块状矩阵
2,存储结构数组在计算机内存中总是一维存放的但是它的顺序在不同的高级语言中不一定相同。
浙江大学研究生学位课程
《实用数值计算方法,4
3.1.2
例如,FORTRAN语言中的矩阵是按列存放的：
,11a,21a?,1na,12a,22a?,2na,1na,2na?,nna
3,存储方式与存储结构是不同的两个概念
4,矩阵的逻辑维数与物理维数逻辑维数,实际参与计算的矩阵阶数物理维数,该矩阵可能出现的最大阶数例如,调用一个子程序，计算一个 4?4矩阵的转置，调用形式为：
CALL MATINV（ A,AI,NL,NP）
这里,NL是逻辑维数,=4。而 NP是物理维数，即数组 A的实际定义维数。
浙江大学研究生学位课程
《实用数值计算方法,5
3.1.2
假定 NP=6，则 4?4矩阵存放于 6?6矩阵中：
7
(1 )
5
(7 )
1
(1 3 )
2
(1 9 )
X
(2 5 )
X
(3 1 )
2
(2 )
3
(8 )
1
(1 4 )
1
(2 0 )
X
(2 6 )
X
(3 2 )
4
(3 )
9
(9 )
0
(1 5 )
5
(2 1 )
X
(2 7 )
X
(3 3 )
8
(4 )
1
(1 0 )
2
(1 6 )
4
(2 2 )
X
(2 8 )
X
(3 4 )
X
(5 )
X
(1 1 )
X
(1 7 )
X
(2 3 )
X
(2 9 )
X
(3 5 )
X
(6 )
X
(1 2 )
X
(1 8 )
X
(2 4 )
X
(3 0 )
X
(3 6 )
NL
NP
NL=4,NP=6 NP,NP
Dimension A(6,6)
A内存放 4?4矩阵
NL?NL求 A( i,j ) 1? i,j?4
但地址是：
（ j-1）?NP+i
浙江大学研究生学位课程
《实用数值计算方法,6
3.1.3 向量范数、矩阵范数定义 1,对于任一向量是,按照一定规则确定一个实数与它对应。 该实数记为,
若 满足下面三个性质：
那么实数 称为向量 x的范数。
设,则它常用的几种范数有：
nRx?
。对
，对任意实数
，
yxyxRy
xx
xxx
n
,3;2;0001
x
,,,,21 Tnxxxx
x
x
浙江大学研究生学位课程
《实用数值计算方法,7
3.1.3
63,,,m a x
53
43
21
1
211
2
1
1
222
2
2
12
n
n
i
in
n
i
in
xxxx
xxxxx
xxxxx
可以验证，以上定义的几种范数均满足三个范数的性质。它们的几何意义见图 3.1
从向量范数出发可以定义矩阵范数。
定义 2：设 A为 n?n阶矩阵。定义
Ax
x
Ax
A
nn Rx
x
Rx
x
10
m a xm a x
浙江大学研究生学位课程
《实用数值计算方法,8
3.1.3
5.5
12
4.7
1
2
x
x
x
321,,xxxx
2x
1x
3x
0.2
5.5
5.4
图 3.1 向量范数的几何意义浙江大学研究生学位课程
《实用数值计算方法,9
3.1.3
这样定义的矩阵范数具有性质：
。对有对
，对
BABARRB
xAAxRx
BABARRB
AAR
AAA
nn
n
nn
,5;,4;,3;2;00,01
显然，这样定义的矩阵范数与向量范数的定义方法有关。
前面三种常用的向量范数相应的矩阵范数是
的最大特征值是 AA
A
T
1
12 73
浙江大学研究生学位课程
《实用数值计算方法,10
3.1.3
93m a x
83m a x
1
1
1
1
1
n
j
ij
ni
n
i
ij
nj
aA
aA
另外，关于范数有一个很重要的等价定理：
定理：有穷线性空间上的一切范数都是等价的。即对任意两种范数有关系式：
,
有关的常数，是与其中
Mm
Mm
,
浙江大学研究生学位课程
《实用数值计算方法,11
3.1.4 线性代数方程组的性态线性方程组 (3-2) 的解完全由 A和 b确定。
在实际问题中,
由于各种原因,A和 b是有误差研究 A和 b的微小摄动对解 x的影响十分重要的这种影响的大小反映了问题的“性态”。
在（ 3-2）中，如果 A-1存在，则解可表示为
x=A-1b
又设 为 A,b的微小摄动，而 是由此而使 x产生的误差。即
bA,x?
bbAAAx 1
浙江大学研究生学位课程
《实用数值计算方法,12
3.1.3
示例 在 4位字长的计算机上解方程组
的。
）是“坏条件”象（微小摄动而得。可以想通过只是则方程组矛盾。事实上改成：如果把误差很大作为主元的消元法解得而取此方程组的真解是
103
103113
113
14 7.537 4.100 0.1
41.1512 2.400 0.3
103
50 0.3,95 1.9
:00 0.3
2.6,66 58.13
:
103
14 7.537 4.100 0.1
41.1512 7.400 0.3
21
21
21
11
21
21
21
xx
xx
xx
a
xx
xx
xx
浙江大学研究生学位课程
《实用数值计算方法,13
3.1.3
现给出一般情形的估计式设则可以证明以下估计式其中可见,k(A)近似表示了方程组求解的误差的相对放大率
11AA?
1232
1
b
b
A
AAk
x
x
11
1
AA
AAAk
浙江大学研究生学位课程
《实用数值计算方法,14
3.1.3
换言之,的大小表示问题的病态程度
k(A)称为计算问题（ 3-2）的条件数
Condition Number
现在再看前例：
2.7 8 42
6.1 4 25
501.5m a x
0.6000.200
6.8258.274
347.1000.1
127.4000.3
1
1
1
1
1
1
1
1
1
1
AAAkAC o n d
A
aA
A
A
n
i
ij
nj
1 AA
浙江大学研究生学位课程
《实用数值计算方法,15
3.1.3
可见计算对象（ 3-10）是病态的。
应当指出：
1,det(A)的值小未必会 引起 A病态例如：
det(A)=0.02，而 A是好条件的 。
2,严格来说，估计式只给出了好条件的充分条件，但不是必要条件。
3,由于数据误差在线性系统中引起的固有的不可靠的使得任何过分精度要求的企图都是徒劳的。
2.010.0
0.010.0A
浙江大学研究生学位课程
《实用数值计算方法,16
3.2 直接法
3.2.1 直接三角分解法 （ LU分解）
考虑
AX=b A非奇异在 Gauss 消去法中,每一步消元过程相当于对 A作一次初等变换 。即左乘一个初等变换矩阵 T。
第一步,
133
100
10
1
1
1
31
21
1
n
l
l
l
T 0
浙江大学研究生学位课程
《实用数值计算方法,17
3.2.1
第 k步
)143(
1
1
1
1
1
,1
nk
kk
k
l
l
T
0
0
0
到 k步以后 A化为当 k=n-1时,A化为上三角阵 U，即因此
)153(11 ATTT kk?
)163(121 UATTT nn?
)173(111211 LUUTTTA n?
浙江大学研究生学位课程
《实用数值计算方法,18
3.2.1
其中为下三角矩阵，称（ 3-17）为 A的三角分解。
由 Ti性质，可以知
1
1
1
1
,
,1
1
in
ii
i
l
l
T
)183(111211 nTTTL?
1
1
1
1
121
3231
21
1
1
1
2
1
1
nnnn
n
lll
ll
l
TTTL
0
0
浙江大学研究生学位课程
《实用数值计算方法,19
3.2.1
Gauss 消去法的基本步骤
4
3
5
311
612
294
3
2
1
x
x
x
75.2
5.0
5
5.225.10
55.00
294
:1
4
13,1
4
22
3
2
1
x
x
x
5.1
5.0
5
1000
55.00
294
:2
5.0
25.13
3
2
1
x
x
x
浙江大学研究生学位课程
《实用数值计算方法,20
3.2.1
所作的初等变换为
5.225.10
55.00
294
311
642
294
10
4
1
01
4
2
001
T1 * A = A1
1000
55.00
294
5.225.10
55.00
294
15.025.10
010
001
T2 * A1 = A2 u
例,设记
084
162
122
A
浙江大学研究生学位课程
《实用数值计算方法,21
3.2.1
则，它的 LU分解为,
112
011
001
110
010
001
102
011
001
200
040
122
240
040
122
110
10
1
240
040
122
084
162
122
102
11
1
1
2
1
1
122
11
21
21
TTL
ATA
ATA
uAAA
TT 左乘左乘浙江大学研究生学位课程
《实用数值计算方法,22
3.2.1
如果有了 LU分解。则解方程组就变得非常容易了。因为
bLY
UXY
bL U XAX
LUA
则令由于 L,U都是三角矩阵，则 Y,X 都可以很容易从上面两式中求得。
如果预先知道 A存在 LU分解，则我们可以直接把这种分解求出来。
浙江大学研究生学位课程
《实用数值计算方法,23
3.2.1
15.0105.1
5.05.055.0
05.04925
5.1
5.0
5
1000
55.00
294
5.15.05.25
4
1
4
5.05
4
2
3,5
4
3
5
15.2
4
1
01
4
2
001
1000
55.00
294
15.2
4
1
01
4
2
001
311
642
294
435
3
32
231
3
2
1
3
21
3
2
1
x
xx
xxx
x
x
x
y
yy
y
y
y
UL
A
b
T
解：
解：
前例：
浙江大学研究生学位课程
《实用数值计算方法,24
3.2.1
但是 Gauss 消去法并不是对任何非奇异矩阵都能顺利进行。 如我们有这样的结果：
任何非奇异的 n?n阶矩阵 A，总存在一个排列矩阵 P，使 PA能进行 LU分解。
上述结果 表明,非奇异矩阵总是能进行选主元的 Gauss消去法 。
01
10
A
浙江大学研究生学位课程
《实用数值计算方法,25
3.2.1
331
421
642
100
010
642
010
100
642
331
421
642
行变换选主元消元消去法：选主元的 G aus s
021 42?
浙江大学研究生学位课程
《实用数值计算方法,26
3.2.1
为使 LU 分解标准化。把 U 改换成 DU
是可行的。其中 D 是非奇异对角矩阵。而
U 则为单位上三角矩阵。如
100
010
2
111
200
040
002
200
040
122
称 A=LDU 为矩阵 A 的 LDU 分解。
定理,设
nk
aaa
aaa
A
nnaA
kkkk
k
k
ij
,,2,1
21
11211
阶矩阵。是一个浙江大学研究生学位课程
《实用数值计算方法,27
3.2.1
表示 A 的各阶主子矩阵。那么,A 存在唯一的 LDU 分解的 充分必要条件 是：
Ak （ k=1,2,….,n ） 都是非奇异。
如 A 为对称正定，或者对角占优，则
A 的各阶主子矩阵均非奇异。
现在我们直接从 A 求 LU 分解：
设 矩阵 A 有 LDU 分解。记 LD=L。 则
A=LU,L为下三角矩阵,U为单位上三角。
这种 LU分解称为 Crout 分解 。
浙江大学研究生学位课程
《实用数值计算方法,28
3.2.1
LUA
u
uu
uuu
llll
lll
ll
l
aaaa
aaaa
aaaa
aaaa
jn
nj
nj
nnnjnn
ijii
nnnjnn
inijii
nj
nj
1000
100
10
1
0
00
000
22
1112
21
21
2221
11
21
21
222221
111211
浙江大学研究生学位课程
《实用数值计算方法,29
3.2.1
行已得到。的前列和的前假定故时，
时，特别，
所以，
因为，
显然，
令：
11
)213(,,2
1
203,,2,1
1
193,,2,1,
,,2,1,1
,
11
1
1
1111
11111
,m i n
1
kukL
nj
l
a
u
ula
i
nilula
j
njiula
LUA
njiu
uUlL
j
j
jj
iii
ji
r
rjirij
ii
ijij
浙江大学研究生学位课程
《实用数值计算方法,30
3.2.1
233,,1
223,,
,,
193,1
1
1
1
1
1
1
nkj
lulau
nkiulal
nkiulla
u
kk
k
r
rjkrkjkj
k
r
rkirikik
k
r
rkirikik
kk
类似地，有，
故有所以根据由于关于 L,U元素的存放从（ 3-22），（ 3-23）可以看出。 A的元素 aij在计算出 或 以后就不再有用了 。
ijl iju
浙江大学研究生学位课程
《实用数值计算方法,31
3.2.1
故 L,U的非零元素便可以存放在矩阵 A中的相应位置上了。
最后,A所存放的元素是：
nnnn
n
n
lll
ull
uul
21
22221
11211
容易证明：
1,即是 Gauss消去法 中各次的主元素。 nn
ll,,11?
浙江大学研究生学位课程
《实用数值计算方法,32
3.2.1
2,如果 A 的各阶主子矩阵均非奇异，
上述过程可以一直进行到底。
矩阵 A 除了以上介绍的 LU 分解以外，
还有一种重要的分解方法可以用于求解线性方程组。即 QR 分解。
定理,任意矩阵 A，总是以分解成正交矩阵 Q 与上三角矩阵 R 的乘积：
如果 A 是非奇异的，则在规定 R的对角元素的符号下，分解式是唯一的。
如果,则 Q 是正交矩阵
RQA
TQQ1
浙江大学研究生学位课程
《实用数值计算方法,33
3.2.1
常见的有 Householder 正交三角分解有了正交三角分解以后，解方程组就变得非常简单：
bQRX
bQ RX
bAX
T?
用 Householder 变换进行分解。从而求解线性方程组的过程是非常稳定的。也不必选主元。特别在方程性质不太好时，
更显其优越性，但计算量大。
Householder 正交三角分解还有其他用途。
浙江大学研究生学位课程
《实用数值计算方法,34
3.2.1
PROGRAM EXAMPLE
DIMENSION A(5,5),B(3),IND(3),C(3)
DATA B /5.0,3.0,4.0/
A(1,1)=4.0
A(3,3)=3.0
N=3
NP=5
CALL LUDCMP(A,N,NP,IND,D)
CALL LUBKSB(A,N,NP,IND,B)
WRITE(*,*) B
C(1)=1.
C(2)=2.
C(3)=3.
CALL LUBKSB(A,N,NP,IND,C)
WRITE (*,*) C
STOP
END
浙江大学研究生学位课程
《实用数值计算方法,35
3.2.1
15.0
5.0
05.0
4
3
5
00000
00000
00105.225.0
0055.05.0
00294
00000
00000
00311
00642
00294
L U B K S B
L U B K S B
L U D C M P
B
A
的变化：
前同的变化：
3
2
1
0.1DIN D 的值：
浙江大学研究生学位课程
《实用数值计算方法,36
3.2.1
程序手稿原程序文件目标程序执行程序输出结果编辑（ WS）
编译（ F77）
连结（ LINK）
执行（程序名）
库修改图 3.2 程序调试流程图浙江大学研究生学位课程
《实用数值计算方法,37
3.2.2 迭代改进由前面讨论，我们知道，由于各种原因，
特别是方程组本身的病态，将会导致解与真解之间的误差不能达到令人满意的程度。有各种处理病态方程组的方法。这里介绍一种简单实用的方法。它还可以用于一般方程组的高精度求解。 该方法的基本思想是利用迭代逐步改善解的精度 （关键在于在迭代过程中有些运算需用双精度进行）
迭代改善过程的框图如下：
浙江大学研究生学位课程
《实用数值计算方法,38
3.2.2
LUAA?的三角分解作
bxL U xx 011,求三角方程组的解
AXb用双精度计算残向量：
L UZZ,求出解的改善量
Zxx 1求改善解：
Z
1
1 xx
xx1?
结束是否图 3.3 迭代改善过程浙江大学研究生学位课程
《实用数值计算方法,39
3.2.2
由上图可知，迭代改进的计算，只要在开始时作一次三角分解。以后只是反复求解三角方程组就可以获得解的改善。
关于解的改善的程度有以下结论：
假设在浮点数的尾数是 t 位（二进制）的计算机上用消去法计算。 A 的条件数为一般可设则：经过一次迭代，解的精度近似地改善了 t-p 位。因此迭代次数 k 满足 k(t-p)?t
就够了。 或近似地用估算式：
pACo n d 2? tp0
1
1
2lo g Z
xtk
浙江大学研究生学位课程
《实用数值计算方法,40
3.2.3 奇异值分解 （ SVD）
解方程组时，常常会出现奇异或十分接近于奇异的情况。用通常的 Gauss消去法或者 LDU分解法难以得到满意的解。对这类问题，我们有一个有效的方法，它叫 奇异值分解法 Singular Value Decomposition
它不但可以诊断出问题所在，而且有时能给出一个有用的数值结果。
SVD方法基于下面这个线性代数定理，
由于定理本身不是我们研究的重点，这里不加证明。
浙江大学研究生学位课程
《实用数值计算方法,41
3.2.3
定理,对任意 M?N矩阵 A,这里 M?N，
它都可以分解成下列形式：
A=UWVT
其中,U是 M?N列正交矩阵；
V是 N?N 正交矩阵；
W是 N?N对角矩阵。且对角元素非负。
即：
Nlk
uu
V
w
w
w
UA
kl
M
i
ilik
T
n
,1
1
2
1
且
浙江大学研究生学位课程
《实用数值计算方法,42
3.2.3
IVVUU
Nlkvv
TT
N
j
kljljk
用矩阵表示：
,1,
1
由于 V是方阵，故它还是行正交即不管矩阵 A是否奇异。 SVD分解都能做到。而且它“几乎”是唯一的。即允许
U,V的列进行变换并进行 W相应交换的前提下是唯一的。
如果 A是 N?N方阵，则 U,V,W都是方阵。且分解式可写为：
.IVV T?
浙江大学研究生学位课程
《实用数值计算方法,43
3.2.3
T
j
T
j
U
w
di agVA
Vwdi agUA
11
由此但是，如果其中有一些 或十分接近于 0，即 A奇异或接近奇异，那么 A-1
就难以求得。因此,SVD分解可以通过判断 的大小来分析 A的性态。
当 A是奇异时,我们就考虑
0?jw
jw
jj wwAC o n d m inm a x?
243 bAx nn
在某种意义下的解。
浙江大学研究生学位课程
《实用数值计算方法,44
3.2.3
由代数知识，有当 rank(A)<rank(A,b)时,
（ 3-24）为矛盾方程组。无解。
当 rank(A)=rank(A,b)时,
（ 3-24）有无穷多解。
在这些解中。有一个 最小长度解在某些场合是有意义的。
在求这个解以前，先介绍几个重要的概念。
1.（ 3-24）定义了一个从向量空间 x到向量空间 b的线性映射 A。
2,如果Ａ奇异，则有,使得
m in )( 2?x
0?x
浙江大学研究生学位课程
《实用数值计算方法,45
3.2.3
0?AX253?
所有满足（ 3-25）的 x全体称 Ａ的零空间 。
零空间的维数称为零度 。
３,所有存在原象 x,使得（ 3-24）成立的向量
b的全体称为 Ａ的域 。 域的维数称为Ａ的秩 。
４,如果Ａ非奇异，则零空间只有零向量；
零度＝０；Ａ的秩＝Ｎ；
如果Ａ奇异，则零度 >0; 秩 <N;
但可以证明,零度＋秩＝Ｎ 。
现在再看看ＳＶＤ分解的意义。
浙江大学研究生学位课程
《实用数值计算方法,46
3.2.3
ＳＶＤ分解事实上分别构造了Ａ的零空间和域的各一组正交基。
由所有 对应的Ｕ的第 j列向量是张成域的一组正交基 。
由所有 对应的Ｖ的第 j列向量是张成零空间的一组正交组 。
现在考虑（ 3-24）的解。如果 的域，
则方程有无穷解。因为零空间的任何向量的
0?jw
0?jw
Ab?
bxV
w
w
w
uuuu
T
n
jnj
1
21
,,,,,
浙江大学研究生学位课程
《实用数值计算方法,47
3.2.3
buxvw
uxvw
xvw
xvw
xvw
xvw
uuu
xV
w
w
w
w
uuuxA
k
w
nk
k
T
kk
n
i
i
T
ii
T
nn
T
jj
T
T
nj
T
n
j
n
0
1
1
22
11
1
2
1
21
,,,,
,,,
故任意的域中向量 b可以由 线性表示，
即这样的 张成域的一组基。
0?kk wu
0?kk wu
（ 1）
浙江大学研究生学位课程
《实用数值计算方法,48
0
1
2211
2
1
2
1
,0,0
),,2,1(
0
0
,,2,10
0
0
0
0
0
k
w
nk
kk
jjj
nn
k
j
T
j
j
T
kk
T
n
T
T
n
T
T
vx
vjw
vvvX
nkvx
vxxv
w
nkxvw
xv
xv
xv
w
w
w
XWV
UU W VA
Ax
Ax
于是内积即得上式两边与的下标对应于的线性组合时表示成由此若令向量正交与即时当故即是非奇异阵，得及由的零空间，即设
3.2.3
（ 2）
浙江大学研究生学位课程
《实用数值计算方法,49
3.2.3
xxxxxx
bbUU
bU
w
di agUW
bU
w
V di agU W V
xU W VAx
w
wbU
w
di agVx
bAxxx
bbAxxAxxA
bxAx
xAx
T
T
j
T
j
T
T
j
j
T
j
现在要证明：。之和：的向量及零空间中表示为特解因为任何解向量都可以的所有解向量中长度最小式还是事实上，
的解。式是问题故则时，取当令的解仍然是方程所以的解时，是当是零空间中的向量。即设解上，仍是方程的解：解的线性组合加在某个
''
'
''
''
)243(
243
1
1
0
1
0,
1
0
,0
0,1 0,0{ jjjj ww令对角阵
0jw ijubAb?的域，是属于因为浙江大学研究生学位课程
《实用数值计算方法,50
3.2.3
的解。是所有解中，长度最小中的特解的长度最小；或者说，，换言之，当
，故综合以上两点：
时，自然当上式长度更小，又才能使分量时，最好取这些下标所以取分量时，由于当因为
x
xx
xxV
xvw
xvj
bu
w
w
xv
xv
xv
bu
w
bu
w
bu
w
xVbUW
xVxV
xxVxxx
T
T
jj
T
j
T
j
j
j
n
T
n
T
T
T
n
n
T
T
TT
TT
T
0
00
00
,0
,0
1
0
1
1
1
'
''
'
'
'
'
22
'
11
2
2
1
1
'1
'
''
浙江大学研究生学位课程
《实用数值计算方法,51
3.2.3
的大小控制。误差值被“零化”；多少大小的奇异值可以较小。但应注意掌握倒可能会使问题求最小长度解此时，把方程作为奇异可能很大差分解直接求解将会使误分解或用时，
件数对病态方程组，即使条最小的解。
即使向量长度最小二乘意义下的解，
式给出了无解，但的域时，当说明：
r
r
AXbr
SV DLU
w
w
Axbr
Ab
i
i
2
1
:
1
m i n
m a x
.2
243.1
浙江大学研究生学位课程
《实用数值计算方法,52
3.2.3
以上讨论用图表示如下：
x
b
A
bAx?
零空间的解dAx?
解的 SVDdAx?
的域A
d
'c
的解'cAx?
最小二乘解的 SV DcAx?零空间?c
图 3.4 解空间示意图浙江大学研究生学位课程
《实用数值计算方法,53
3.2.3
T
v
T
T
VWU
U W VASV DA
bAX
bAr ankAr ank
bA
bAX
3
3
1
6
2
0
3
1
6
1
2
1
3
1
6
1
2
1
000
030
002
3
2
23
1
2
1
3
1
23
4
0
3
2
23
1
2
1
,
3,,2
0
0
1
,
23
1
3
1
6
1
6
1
32
1
2
3
22
3
2
3
2
23
1
3
1
6
1
6
1
32
1
2
，，
其中分解：进行把最小长度解。现求最小二乘意义下的为矛盾方程组，无解所以因为，
是矛盾方程组。例：下列线性代数方程
1u 2u
浙江大学研究生学位课程
《实用数值计算方法,54
3.2.3
32
1
9
1
9
1
36
1
22
1
0
3
1
2
1
2
01
23
1
23
4
23
1
2
1
0
2
1
2
3
1
6
2
0
1
2
2211
3
3
21
3
bUVx
ukukAx
Rx
kvA
uu
v
Ar ankA
T
：原方程组的最小长度解即有
，，且域的一组基为：域的维数为
，且零空间的一组基为零度为奇异，且由此，
浙江大学研究生学位课程
《实用数值计算方法,55
3.3 迭代法迭代法适用于解高阶稀疏非病态方程组。它只需要存储非零元素。
但对有些问题，迭代可能发散或收敛很慢。
3.3.1 Jacobi 迭代法设有方程组其中 A是 N?N非奇异矩阵。故有唯一解。
把（ 3-26）改写成迭代形式
bAX263?
浙江大学研究生学位课程
《实用数值计算方法,56
273 gBXX
于是设 X(0)是一个任意的初始迭代向量。
代入（ 3-27）的右边，得到新的向量记为 X(1)：
一般地有：
我们得到向量序列如果 收敛于 X,即：
则向量 X 是（ 3-27）的解。
现讨论解的求法及收敛性
gBXX )0()1(
)283(,2,11)(kgBXX kk
,1,0,?kX k
kX
XX k
klim
3.3.1
浙江大学研究生学位课程
《实用数值计算方法,57
3.3.1
设矩阵 A的对角元素 。记则 D非奇异。将 A表示成和式其中：
0?iia
nnaaad ia gD,,,2211
293 uL CCDA
方法的定义是：唯一确定。述矩阵迭代法都可以用上我们会看到，
J ac obiCCD
SO RSGJ ac obi
a
aa
aaa
C
aa
aa
a
C
uL
n
n
n
u
nn
L
,,
,,
0000
000
00
0
,
0
00
000
0000
3
223
11312
21
3231
21
浙江大学研究生学位课程
《实用数值计算方法,58
303,,2,1,
1
1
nibxaxa
n
ij
j
i
k
jij
k
iii?
3.3.1
利用（ 3-29）式的矩阵符号。 Jacobi 法可以表示成：
这里 B 是 Jacobi 迭代矩阵。其定义为
gBxx kk 1
333,,,
323,,,
313
21
1
21
11
T
n
Tk
n
kkk
uL
bbbDg
xxxx
ADICCDB
Jacobi 方法收敛的充要条件是：
迭代矩阵 B的谱半径 1)(?B?
浙江大学研究生学位课程
《实用数值计算方法,59
3.3.1
利用这个判断标准。
，故较多的计算可以较方便检验由这是由于法收敛方法的一个充分条件：给出为此，的条件较难检验由于的特征值。是这里模最大者：谱半径定义为特征值的
B
BB
J ac obiB
J ac obi
B
Bni
B
i
i
ni
1
,1
,,2,1
m a x
1
1
,,m a x
i
iiiB
即为设单位长度的特征向量
B
BBB
B
i
iiiii
iii
于是则浙江大学研究生学位课程
《实用数值计算方法,60
3.3.1
现在讨论一下 Jacobi方法的收敛速度的快慢及影响速度的因素。
的精确解。是方程组其中，
方法成立，则定理：若
273
353
1
343
1
1
1
01
x
xx
B
B
xx
xx
B
B
xx
J ac obiB
kkk
k
k
B
BI
1
11
示：其证明作为思考题，提浙江大学研究生学位课程
《实用数值计算方法,61
3.3.1
（ 3-34）式可作为误差估计式，并且说明 ||B||越小,x(k) 收敛越快。（ 3-35）式说明了只要 ||B||不是很小，当相邻两次迭代向量 x(k) 和 x(k-1) 很接近时，解 x* 和 x(k) 也是很接近的。因此，可以用 || x(k)- x(k-1)|| 作为迭代终止的判别式。这个结论对逐次逼近法也成立。
特别，应该指出，以上讨论的收敛性及收敛速度均是对一般迭代式（ 3-27）而言。因此，其结论对以后要讨论的 G--S迭代法和 SOR法也成立。
浙江大学研究生学位课程
《实用数值计算方法,62
3.3.1
Jacobi 迭代法的计算步骤：
)1
,2,13;2;1
1
1
0
11
akk
xxb
gBxxa
k
x
bDgADIB
kk
kk
，转则结束，否则若计算对和初始迭代向量给出精度小量
，计算
浙江大学研究生学位课程
《实用数值计算方法,63
n
k
n
k
n
k
n
k
nn
kk
k
nn
kk
k
i
k
j
k
j
k
i
k
i
n
k
n
k
n
k
n
k
nn
kk
k
nn
kk
gxbxbx
gxbxbx
gxbxbx
xxijx
xx
gxbxbx
gxbxbx
gxbxbx
J ac obi
0
0
0
:1,,1
0
0
0
2211
2
1
21212
1
1
1
1
2121
1
1
1
22
1
11
2
1
2
1
1212
1
1
1
1
2121
来求代替故可考虑用已求出。时，由上可知当计算迭代格式：仔细写出
3.3.2 Gauss--Seidel 迭代法 。
G--S迭代法的基本思想来自 Jacobi
迭代法。只是迭代矩阵 B不同。
浙江大学研究生学位课程
《实用数值计算方法,64
3.3.2
这就是 Gauss--Seidel 迭代格式，它也可写为：
若用（ 3-29）式记号，则 G--S迭代式可表示为：
或其中
),,1(
1
1
1
1
ni
gxbxbx i
n
ij
k
jij
i
j
k
jij
k
i
bxCxCD kukL 1
3631' hxLx kk
UL CDUCDL
bDLIhULIL
11
111 373'
，
，
浙江大学研究生学位课程
《实用数值计算方法,65
有关收敛性和收敛速度的讨论完全与 Jacobi 方法相同，只是迭代阵为 L'。
3.3.3 逐次超松弛（ SOR） 法
SOR法定义为,
法可写成：的矩阵记号利用法。法退化为时，当
。是实数，称为松弛因子其中
S O R
S e id e lG a us sS O R
nix
gxbxbx
k
i
i
n
ij
k
jij
i
j
k
jij
k
i
,293
1
,,2,11
1
1
1
1
1
浙江大学研究生学位课程
《实用数值计算方法,66
3.3.3
1)(
,
393
1
383
1
11
11
1
1
1
1
L
SO RJ ac obi
CDUCDL
bDLIh
IULIL
hxLx
Dx
bxCxCDx
UL
kk
k
k
L
k
L
k
法收敛的充要条件是法类似，与其中：
或浙江大学研究生学位课程
《实用数值计算方法,67
2.3.13
1)(1
2
:
403
40320
1)(
L
SO R
SO R
A
SO R
L
opt
opt
，可以证明对某些特殊类型的矩阵
。为称为最佳松弛因子，记因子法收敛速度最快的松弛使收敛的充分条件。也是是对称正定矩阵，则另外，若收敛的必要条件是
，因此可以证明曹志浩等求根》参考《矩阵计算与方程
浙江大学研究生学位课程
《实用数值计算方法,68
3.3.3
三种迭代法的比较：
1,一般情况下,J方法与 G--S方法比较并无优劣。收敛情况与速度均不一定。
例：
011
202
110
022
101
220
19.09.0
9.019.0
9.09.01
B
B
A
法收敛方法不收敛
SG
J
法不收敛方法收敛
SGL
JB
,1'
,1
1'
,1
L
B
浙江大学研究生学位课程
《实用数值计算方法,69
3.3.3
2,但是具有相容次序的矩阵，在相同精度要求下，迭代次数分别为：
Jacobi 方法,1154
G--S 方法,578
SOR 方法,59
可见 对这类矩阵,G--S 法比 J方法快一倍，
而 SOR 法的收敛速度可提高一个数量级。
最后介绍一类对于三种方法都收敛的矩阵。 先定义可约矩阵和对角占优矩阵：
opt
浙江大学研究生学位课程
《实用数值计算方法,70
3.3.3
(1),称 n阶方阵为 可约矩阵,如果存在 n阶排列阵 P，使得其中 A11为 r阶方阵,A22为（ n-r）阶方阵不是可约矩阵就是 不可约矩阵
(2),称 n 阶方阵 A是 对角占优矩阵,如果若上式严格不等号成立，则称 A为 严格对角占优矩阵 。
4130
22
1211
A
AAPAP T
4231
1
niaa
n
ij
j
ijii,
浙江大学研究生学位课程
《实用数值计算方法,71
3.3.3
(3),称 A是 不可约对角占优矩阵 。如果 A是不可约的且对角占优，而（ 3-42）式中至少对一个 i 有严格不等号成立。
我们有如下结论,
1,若 A是严格对角占优或不可约对角占优矩阵，则 A非奇异。
2,若线性方程组系数 A是上述矩阵，则：
1） J方程和 G--S方法都收敛。
2）当 SOR方法收敛。10
浙江大学研究生学位课程
《实用数值计算方法,72
3.4 稀疏矩阵前面介绍了 直接法 和 迭代法 。但是实践中如何选择计算方法是一个很复杂的问题。需要考虑的因素很多，主要有：
1.方程组的性质； （病态、对角占优、正定等）
2.相类似问题出现的次数；
3.方程组的阶数；
4.计算机的容量、字长、运算速度等；
5.系数矩阵的结构； （稀疏、三对角等）
6.对结果的精度要求。
浙江大学研究生学位课程
《实用数值计算方法,73
3.4
迭代法 是解超高阶线性代数方程组的重要方法之一，但是它也有缺点：
不适合解病态问题；
精度不高；
收敛速度依赖于加速因子的选取；
收敛性不能保证；
所需乘除运算量大。
在这些方面，相对来说，直接法运算步骤一定，运算次数有限，精度高。对解大型稀疏矩阵问题可以利用和保持系数矩阵的稀疏性来降低存储容量和缩短计算时间。
浙江大学研究生学位课程
《实用数值计算方法,74
3.4.3
但是，为了保持系数矩阵的稀疏性。必须使用特殊的计算方法和高级程序设计技术。
近年来在这方面有较大的发展。其优越性大大超过了迭代法。
在这一节里，不详细介绍计算方法，
只给出一些常用的标准稀疏矩阵形式。同时给出一种一般稀疏矩阵的存储结构。这一节里还给出一个有效的计算公式。
浙江大学研究生学位课程
《实用数值计算方法,75
3.4.1 稀疏矩阵的结构和存储从本质上讲，解稀疏问题的方法与一般 Gausss消去法 LU分解法并没有什么不同，只不过前者更注重于运算过程中对零元素填入量的控制。
稀疏系数矩阵问题的直接法必定依赖于矩阵的稀疏形状，下面列举的是常见的标准的稀疏矩阵。
1.稀疏矩阵的结构浙江大学研究生学位课程
《实用数值计算方法,76
3.4.1
带对角块三角阵块三角阵单边块对角双边块对角
e
a
c d
b
单边块三角
f
图 3.5
浙江大学研究生学位课程
《实用数值计算方法,77
3.4.1
双边带三角单边带对角双边带对角其他其他
k
g
i j
h
图 3.6 稀疏矩阵的形式浙江大学研究生学位课程
《实用数值计算方法,78
3.4.1
2,稀疏矩阵的存储一般矩阵的存储方式有两种：满存和部分存储。稀疏矩阵则采用后一种方式。
稀疏矩阵的种类不同，所采用的存储方法也不一样，但是不论采用什么方法，它的标准有两条：
1）节省存储空间 ；
2）存取方便，即节省存取时间。
存储方法的选择一般要根据矩阵的种类和对程序的要求在上述两个标准之间权衡。
浙江大学研究生学位课程
《实用数值计算方法,79
3.4.1
(1) 带对角矩阵的存储
imjkmjib
mji
a
mmjia
ki
ij
ij
1
0
,0
,
，这里当当为半带宽。称当图 3.7 带对角矩阵的存储方式
nnijaA )( )12()( mnijbB
浙江大学研究生学位课程
《实用数值计算方法,80
3.4.1
(2)块三角型矩阵的存储
H1
i1
iL
il+1
N1
H2 N2
Hl+1
aij
Hl Nl
i
1
i
2
i
l - 1
i
l
N
1
N
2
N
L
bb
ID2?L
b
1
b
2
b
M
1
2
1
11
ll
l
lll
HHN
iiH
11 lM NMB
1
1
0
,
0,
llk
l
ij
ll
ijnnijnn
ijib
ij
a
iii
ijaaA
则设时当
1
2
1
1
ll
l
h
ll
h
ijii
iiii
Nk
这里：
iL
图 3.8 块三角型矩阵的存储形式浙江大学研究生学位课程
《实用数值计算方法,81
3.4.1
(3) 一般稀疏矩阵的存储
B,一维存储非零元素（按行）
（ m）
（ m）
JA,B中相应元素所在的列号 （ A中）
IA,每行第一个非零元素在 B中的下标。
k
ij b
elsea 0
jkJA
iIAkiIA
k
）
）
使得如果
2
11
:,
（ n）
浙江大学研究生学位课程
《实用数值计算方法,82
3.4.1
例如：
8706
0050
4003
0210
A
1 2 3 4 5 6 7 8
2 3 1 4 2 1 3 4
1 3 5 6
IA
JA
B
图 3.9 矩阵 A非零元素存放 的一般形式浙江大学研究生学位课程
《实用数值计算方法,83
3.4.1
(4) 对称矩阵的存储因为故只要存储即可。
任何形状的矩阵，包括稀疏的和非稀疏的，只要是对称的。
均可节省将近一半的存储量。
njiaa jiij,,2,1,
njiija ij,,2,1,,
浙江大学研究生学位课程
《实用数值计算方法,84
3.4.2 Sherman-Morrison & Woodburg 公式我们知道，矩阵求逆问题等价于解一个线性方程组。如果我们已经得到矩阵 A 的逆矩阵 A-1。 现在想在 A上作个小小变动。如：一个元素、一行、一列，
或者部分元素，那么是否可以不再做复杂的求逆运算，而只是在原来的 A-1的基础上作些适当的运算就可得到变动后的矩阵的逆呢？回答是肯定的。
假定矩阵变化为：
433 VUAA
浙江大学研究生学位课程
《实用数值计算方法,85
3.4.2
附带复习一下内积，外积的定义：
T
nnnn
n
n
jiij
nn
T
n
i
ii
T
nn
uv
vuvuvu
vuvuvu
vuvuvu
vu
njivuA
Auvvu
uvvuvuvu
v
v
v
v
u
u
u
u
21
22212
12111
1
2
1
2
1
,,2,1,
,,
,
即其中：
外积：
内积：
内积、外积均有结合律
（叉积）
（点积） (或,u?v)
浙江大学研究生学位课程
《实用数值计算方法,86
3.4.2
这里 u?v 是外积，它表示 A的变动。
Sherman--Morrison公式给出了一个计算的方法：
1 vuA
443
1
1
11
1
2111
1111
1
1
11
AvuA
A
AvuAA
AvuAvuAvuAI
AvuAIvuA
这里从（ 3-44）（ 3-45）式可以看出，给了
A-1和 u,v只需作矩阵的乘积运算和一个向量的加法：
4531 uAv?
浙江大学研究生学位课程
《实用数值计算方法,87
3.4.2
473
1
463,,
11
11
Wz
AA
zvvAWuAz
T
则：
令：
Sherman--Morrison公式可以直接应用于一类稀疏线性方程组上。如果我们已经有一个有效的算法得到矩阵 A的逆（如 A
是三对角阵。或者其他标准稀疏形状的矩阵），那么要解 A的稍作变动的矩阵决定的方程组，只要运用一次或多次
Sherman--Morrison 公式即可。
但对有些稀疏问题，这个公式并不能直接运用。理由很简单，因为 A-1不可能在机器内部全部都储存着。在实际运用中这个公式还有变化。
A
浙江大学研究生学位课程
《实用数值计算方法,88
3.4.2
变化一,求解求解线性系统问题：
首先可以对 A用有效的方法解两个辅助问题：
得到了向量 y和 z后，（ 3-48）的解既是：
483 bXvuA
493, uAzbAy
5031
z
zv
yvyx
浙江大学研究生学位课程
《实用数值计算方法,89
3.4.2
z
zv
yv
yx
uAzbAy
uAzbAy
zbAv
bzWbzWbWz
zv
bWz
bA
b
zv
Wz
A
bvuAX
TT
1
,
,
1
1
11
1
1
1
1
有时故当由于证明浙江大学研究生学位课程
《实用数值计算方法,90
3.4.18
变化二,即 Woodbury 公式当在标准稀疏矩阵上添加的元素不只是一行或一列时，我们就不能简单地重复上述过程。因为新的 A-1 并没有存储下来。
这时要用 Woodbury 公式：
.,
513
111111
NPPNVUNNA
AVUAVIUAAVUA TTT
矩阵，是阵，是这里，
（ 3-51）式主要是计算因此可以看到，（ 3-51）式把一个 N?N阵求逆问题转化为一个 P? P阵的有求逆问题。由于 P<<N,故问题得以大大简化。
11 )( PPT UAVI
浙江大学研究生学位课程
《实用数值计算方法,91
3.4.2
。的右边式子也同样计算左乘用均为单位阵即可：
右乘和左乘的右边用只要证公式：证明
513
513
11
1
1
111
1
1
111
1
1
111
T
TT
TTTT
TTTT
TTT
T
UVA
IUVAUVAI
VUAVIUAVIUAUVAI
UVAIVUAVIUAUVAI
UVAAVUAVIUAA
b
VUA
W ood bur y
Woodbury公式和连续地运用 Sherman--
Morrison公式的联系在于：
如果记
PP vvVuuuU,,,,,,121
523
1
P
k
kk
T VUAVUA
则有浙江大学研究生学位课程
《实用数值计算方法,92
3.4.2
上式表明，如果已经求得 A-1，计算的方法有两种：
（ 1）利用（ 3-51），计算一个 P?P 的逆矩阵。
（ 2）利用（ 3-47），执行 P次。
如果 A-1没有存储，则解下列问题的步骤是
(1) 解 P个辅助方程组：
得矩阵
1 TUVA
533
1
bxvuA
P
k
kk
543,,2,1 PiuAz ii?
PzzzZ,,,21
浙江大学研究生学位课程
《实用数值计算方法,93
3.4.21
(2) 计算 P?P矩阵的逆
(3) 再解辅助方程组
(4) 得方程组（ 3-53）的解因为：
5531ZVIH T
563 bAy
573 yVHZyx T
bAVHUAbA
bAVUAVIUAA
bVUAx
T
TT
T
111
11111
1
yVHZy T
浙江大学研究生学位课程
《实用数值计算方法,94
2.5 特殊矩阵在很多问题中，我们会遇到特殊类型的矩阵。对于这些具有特殊性质的矩阵，
不仅存储方式上应予研究并很好利用其性质以节省存储，在计算方法上也要利用它们的好的性质以减少运算量和提高算法稳定性。
在这一节里，我们研究几种特殊矩阵，
对称、正定、带状，和三对角矩阵。
由第二节的定理知，若的各阶主子式
nnijaA
。则 L D UA,0
浙江大学研究生学位课程
《实用数值计算方法,95
3.5
同样我们可以证明，若对称矩阵Ａ的各阶主子式 其中
Ｌ是单位下三角阵，Ｄ为对角矩阵。
这个分解叫做对称矩阵的 分解。
对一般非奇异对称矩阵，则存在排列矩阵Ｐ，使故可通过选主元方式达到 分解。
解线性问题 可转化为解三个简单的线性问题：
。则 TL D LA,0
TLDL
TLDL
TLDLPA?
bAX?
zxL
yDz
PbLy
T?
浙江大学研究生学位课程
《实用数值计算方法,96
3.5
如果Ａ是正定矩阵，则Ａ是对称的并且Ａ的各阶主子式不等于０。因此，
存在且故可令：
从而这里Ｇ是下三角矩阵。这种分解叫做对称正定矩阵的平方根分解。
583 TL D LA
niddddd ia gD iinn,10,,,,2211
nndddd ia gD,,,~ 2211
593
~~~~
T
TTT
GG
DLDLLDDLLD LA
浙江大学研究生学位课程
《实用数值计算方法,97
3.5
ijgga
g
g
gggg
ijgga
gag
ggga
G
ggg
gg
g
aaa
aaa
aaa
GGA
j
k
jkikij
jj
ij
jjij
j
k
jkik
j
k
jkikij
i
k
ikiiii
ii
i
k
ik
i
k
ikii
T
nnnnnnnn
n
n
T
1
1
1
1
1
2
1
1
1
2
2
1
1
2
1
2
21
2221
11
21
22221
11211
1
而，
故，
即：
因为
浙江大学研究生学位课程
《实用数值计算方法,98
3.5
平方根分解的算法如下：
个算法是稳定的。能够得到控制，因此这量的平方根分解法的中间知，由所以的元素有如下关系的对角元素和因为执行对
ij
iiij
i
j
ijii
j
k
jkikijjjijij
i
k
ikiiiiii
gA
ijag
ga
GA
ijggaggaii
gagai
ni
613
613
1,,1)
603)
,,2,1
1
2
1
1
1
1
1
2
浙江大学研究生学位课程
《实用数值计算方法,99
3.5
这个算法的另一个优越性是 不用选主元 。
但是从（ 3-60）知道，这个算法在求要进行平方根运算。
为了避免开方运算，有时我们仍要用
ＬＤＬ Ｔ 分解。
在 3.5 中，引入了带宽为 2m+1的带状矩阵。他们满足：
当带状矩阵当 m<<n 时才有意义 。
带的宽度依赖于方程组未知量的次序和 方程的顺序。
nig ii,1?
0 ijamji 时，
浙江大学研究生学位课程
《实用数值计算方法,100
3.5
例如,m=1 的三对角方程的矩阵
623
1
11121
343332
232221
1211
nnnn
nnnnnn
aa
aaa
aaa
aaa
aa
A
若Ａ的相邻任意两行交换一下，则带宽从
３增加到５。
对带状矩阵，有以下定理：
633,,1,0
,,1,0
,,
12
mjjjiu
jmjmjil
uUlL
m
ij
ij
ijij
仅对仅对三角带状矩阵，即：
是则有ＬＵ分解的带状矩阵定理：如果带宽为浙江大学研究生学位课程
《实用数值计算方法,101
3.5
我们在上述定理中假定ＬＵ的分解存在。事实上未必总是如此，可能需要交换方程。这就可能要增加带的宽度。但不可能增大２倍。
但有些实际问题中遇到的带型矩阵性质是比较好的，并不需要选主元。
三对角矩阵（ m=1） 在带状矩阵中占有重要的位置。若（ 3-62）中的各阶主子式
０，则可用 追赶法 很方便求解 Ax＝ b。
所幸的是实际应用中，有很大一类矩阵是对角占优矩阵,因而各主子式?０。
浙江大学研究生学位课程
《实用数值计算方法,102
第三章 习题
分解法求解。分解及用态；试分析以上方程组的性设有线性代数方程组２
。有极小值，并求极小值取何值时，，求：当
１１
２
设Ａ＝１
SV DLU
xxxx
xxxx
xxxx
xxxx
AC ond
2
1
9
1
7
1
6
1
5
1
4
1
8
1
6
1
5
1
4
1
3
1
7
1
5
1
4
1
3
1
2
1
6
1
4
1
3
1
2
1
.
.
4321
4321
4321
4321
浙江大学研究生学位课程
《实用数值计算方法,103
次数，并讨论之。
。试比较迭代，精度控制在解迭代法求解。取初始的因子为以及取松弛请分别用设有线性代数方程组所需的乘除法次数。
求解线性代数方程组阶下三角矩阵，试计算是设矩阵
4
54
543
432
321
21
100
46.1
,,
12
02
12
02
12
.4
.3
x
SO R
Se i de lG aus sJ ac obi
xx
xxx
xxx
xxx
xx
bLx
nL
浙江大学研究生学位课程
《实用数值计算方法,104
解，并解释解的含义。时分别计算当取的维数及基。线性变换的零空间及域定义的的性态及由试分析分解式为的其中系数矩阵阶代数方程组：设有
SV D
b
b
AbAx
V
di agW
U
U W VA
SV DA
bAx
T
T
T
T
3334.14714.00
,1112
1
5774.05774.05774.0
8165.04082.04082.0
07071.07071.0
1001
6667.03333.06667.0
2357.09428.02357.0
7071.007071.0
33.5
2
1
浙江大学研究生学位课程
《实用数值计算方法,105
。的逆矩阵非奇异矩阵程组的程序来计算一个性代数方试说明如何利用求解线亦是对称正定。是对称正定，则若亦是对称矩阵试证明：
约化为消去法一步，
，经为对称矩阵，且设
1
2
2
2
12111
11
.7
2
1
0
0.6
AA
AA
A
A
Aa
AA
AG aus s
aA
