高等工程数学第二章 线性代数方程组数值解法

第2章 线性代数方程组数值解法 研究
n 阶线性方程组
Ax b =
的数值解法.
()
ij A a =是
n n
⨯矩阵且非奇异，
12(,,
,)T
n x x x x =,
12(,,
,)T n b b b b =
两类数值方法：
(1) 直接法：通过有限次的算术运算，若计算过程中没有舍入误差，可以求出精确解的方法.
Ax b
Gx d ==等价变换
G 通常是对角矩阵、三角矩阵或者是一些结构简单的矩阵的乘积.
(2) 迭代法：用某种极限过程去逐次逼近方程组的解的方法.
(1)()i i Ax b
x Bx k x Bx k +==+−−−−−→=+等价变换
建立迭代格式
，0,1,
i =
一、向量范数与矩阵范数 1. 向量范数
【定义】 若对n
K 上任一向量x ，对应一个非负实数
x ，对任意,n
x y R ∈及K α∈，满足如下条件
(向量范数三公理) (1) 非负性：0x ≥，且0x =的充要条件是0x =；
(2)
齐次性：
x x
αα=；
(3)
三角不等式：
x y x y
+≤+.
则称x
为向量x
的范数.
常用的向量范数： (1) 1—范数
11
n
i
i x x ==∑
(2) 2—范数
12
2
21
()
n
i i x x ==∑
(3) ∞—范数
1max i
i n
x
x ∞
≤≤=
(4) 一般的p —范数
11
()
p
n
p
i p
i x
x ==∑
2. 矩阵范数
【定义】 若n n
K ⨯上任一矩阵
()ij n n A a ⨯=
，对应一个非负实数A ，对任意的,n n
A B K ⨯∈和
K α∈，满足如下条件(矩阵范数公理)：
(1) 非负性：0A ≥，且0A =的充要条件是0A =；
(2)
齐次性：
A A
αα=；
(3)
三角不等式：A B A B +≤+；
(4)
乘法不等式：
AB A B
≤.
则称A
为矩阵A
的范数.
矩阵范数与向量范数是相容的：
Ax A x
≤
向量范数产生的从属范数或算子范数：
10max max
x x Ax
A Ax x
=≠==
常见从属范数：
(1) 1—范数
111
max ||
n
ij j n
i A a ≤≤==∑
(2) ∞—范数
11
max ||
n
ij i n
j A a ∞≤≤==∑
(3) 2—范数
2A =谱半径1()max ||H i i n A A ρλ≤≤=，i
λ
为
H A A 的特征值.H A 为A 的共轭转置. 注：矩阵
A 的谱半径不超过A 的任一范数，即
()A A ρ≤
范数等价性定理：
,s t x x
为
n R 上向量的任意两种范数，则存在常数12,0c c >，使得
12,n
s t s c x x c x x R ≤≤ ∀∈.
注：矩阵范数有同样的结论. 【定理2.1】
是任一向量范数，向量序列()k x 收敛于向量*
x 的充要条件是
()*
0,k x x k -→ →∞
二、 Gauss 消去法 1．顺序Gauss 消去法 将方程
Ax b =写成如下形式
11112211,121122222,11122,1
n n n n n n n n nn n n n a x a x a x a a x a x a x a a x a x a x a ++++++=⎧⎪
+++=⎪⎨⎪
⎪+++=
⎩
其中记
,1,1,2,
,.
i n i a b i n +==
消元过程：
第一次消元：设110a ≠，由第2,3,
,n 个方程减去第一个方程乘以1
111/(2,3,
,)i i m a a i n ==，则将方程
组中第一个未知数1x
消去，得到同解方程
11112211,1
(1)
(1)(1)
22222,1
(1)(1)(1)22,1
n n n n n n n nn n n n a x a x a x a a x a x a a x a x a ++++++=⎧⎪ ++=⎪⎨⎪
⎪ +
+=
⎩
其中，
(1)
11,2,3,,;2,3,,,1ij ij i j a a m a i n j n n =-==+. 1111/i i m a a =，2,3,,i n =.
第二次消元：设(1)
220a ≠，.由第2,3,,n 个方程减去方程组中的第2个方程乘以
(1)(1)
2222/(3,4,
,)i i m a a i n ==，则将方程组第2个未知数2x 消去，得到同解方程
11112213311,1
(1)(1)(1)(1)2222322,1(2)(2)(2)33333,1(2)(2)(2)33,1n n n n n n n n n n nn n n n a x a x a x a x a a x a a x a a x a x a a x a x a ++++++++=⎧⎪ +++=⎪⎪ ++=⎨⎪⎪
⎪ ++=
⎩
其中
(2)(1)(1)22, 3,4,
,; 3,4,
,,1ij ij i j a a m a i n j n n =-==+. (1)
(1)2222/i i m a a =，
3,4,,i n =.
经过
1n -次消元后，原方程组变成等价方程组
11112213311,1(1)(1)(1)(1)2222322,1(2)(2)(2)33333,1(1)(1),1
n n n n n n n n n n n nn n n n a x a x a x a x a a x a a x a a x a x a a x a +++--+++++=⎧⎪ +++=⎪⎪ ++=⎨⎪⎪
⎪ =
⎩
其中
()(1)(1), 1,2,,k k k ij ij ik ij a a m a i k k n --=-=++， 1,2,
,,1j k k n n =+++.
(1)(1)/k k ik ik kk
m a a --=，1,2,,i k k n =++；1,2,,1k n =-.
回代过程：
(1)(1)
,1(1)(1)(1),1,,1/[]/,1,2,,2,1.
n n n n n m n i i i i
i n i j j i j j i x a a x a a x a i n n --+---+=+⎧=⎪⎨=-=--⎪⎩
∑
计算量：按常规把乘除法的计算次数合在一起作为Gauss 消去法总的计算量，而略去加减法的计算次数. 在消去过程中，对固定的消去次数(1,2,,1)k k n =-，有：
除法
(1)(1)
,,/,1,1,,k k ik i k k k m a a i k k n --= =++ 共计
n k -次；
乘法
(1),,1
,2,,;1,2,
,,1k ik k j m a i k k n j k k n n - =++ =+++
共计()(1)n k n k --+次.
因此，消去过程总的计算量为
1
3
11[()(1)]3n k M n k n k n k n
-==--++-≈∑ 回代过程的乘除法计算次数为2
1()2n n +.与消去法计算量相比可以略去不计.所以， Gauss 消去法总的计算量大
约为3
13n .
2. Gauss-Jordan 消去法
Gauss-Jordan 消去法是Gauss 消去法的一种变形.此方法的第一次消元过程同Gauss 消去法一样，得到
(1)(1)
(1)(1)
11112213311,1(1)(1)
(1)(1)22223322,1(1)(1)
(1)(1)32233333,1(1)(1)(1)(1)2233,1,,,
,
n n n n n n n n n n
n nn n n n a x a x a x a x a a x a x a x a a x a x a x a a x a x a x a ++++⎧++++=⎪ +++=⎪ +++=⎨ +++
=⎪⎪⎪
⎪⎩
其中，(1)
11,2,,,1j
j a a j n n ==+. 第二次消元：设(1)
220a ≠，由第1,3,4,
,n 个方程减去第2个方程乘以(1)(1)
2222/(1,3,4,
,)i i m a a i n ==,则得
到同解方程组
(2)(2)(2)
11113311,1(1)(2)
(2)(2)22223322,1(2)
(2)(2)33333,1(2)(2)33,1,,,
n n n n n n n n n n
nn n n n a x a x a x a a x a x a x a a x a x a a x a x a +++++ +++= +++= ++= ++
=(2),
⎧⎪⎪⎪⎨⎪⎪
⎪⎩
继续类似的过程，在第k 次消元时，设
(1)k kk a -，将第i
个方程减去第
k 个方程乘以(1)(1)
/k k ik ik kk m a a --=，
这里
1,3,4,
1,1,,i k k n =-+.经过1n -次消元，得到
(2)
1111,
1(1)(2)2222,1(2)(2)
33,1,
,,n n n n n a x a a x a a x a +++⎧ =⎪ =⎪⎪ ⎨⎪⎪
⎪ =
⎩
其中
()(1)(1),1,2,,1,1,,k k k ij ij ik kj a a m a i k k n --=-=-+；
1,2,,,1; 1,2,,1j n n k n =+=-.
此时，求解回代过程为
(1)(1),1/,1,2,,n i i i n ii
x a a i n --+= =
经统计，总的计算量约为
3
1
2M n ≈次乘除法. 从表面上看Gauss-Jordan 消去法似乎比Gauss 消去法好，但从计算量上看Gauss －Jordan 消去法明显比Gauss
消去法的计算量要大，这说明用Gauss-Jordan 消去法解线性方程组并不可取.但用此方法求矩阵的逆却很方便. 3．列选主元Gauss 消去法
在介绍Gauss 消去法时，始终假设(1)0k kk a -≠，称(1)k kk a -为主元.若(1)
0k kk
a -=，显然消去过程无法进行.
实际上，既使
(1)0k kk
a -≠，但
(1)k kk
a -很小时，用它作除数对实际计算结果也是很不利的.称这样的(1)k kk
a -为小主元.
【例2.2】设计算机可保证10位有效数字，用消元法解方程
11
12120.3100.7,0.9,x x x x -⎧⨯+=⎪⎨ +=⎪⎩
【解】经过第一次消元：第2个方程减去第1个方程乘以212111/m a a =得
1112(1)(1)222230.3100.7x x a x a -⎧⨯+=⎪⎨ =⎪⎩
其中
(1)1222222111/0.333333333310a a a a =-=-⨯，
(1)123323211113(/)0.233333333310a a a a a =-⋅=-⨯
于是解得
(1)(1)223221/0.7000000000,0.0000000000,x a a x ⎧==⎪⎨
=⎪⎩
而真解为
120.2,0.7x x = =
注：造成结果失真的主要因素是主元素11a
太小，而且在消元过程中作了分母，为避免这个情况发生，应在消元之前，作行交换.
【定义】 若 (1)(1)
||max ||
k k k r k ik k i n
a a --≤≤=，则称
(1)||
k k r k a - 为列主元素. k r 行为主元素行，这时可将第 k r
行
与第
k 行进行交换，使
(1)||
k k r k a - 位于交换后的等价方程组的 (1)
k kk a - 位置，然后再施实消去法，这种方
法称为列选主元Gauss 消去法或部分主元Gauss 消去法.
【例2.3】 应用列选主元Gauss 消去法解上述方程. 【解】 因为
2111a a >，所以先交换第1行与第2行，得
1211120.9,0.3100.7,x x x x -⎧+=⎪⎨⨯+=⎪⎩ 然后再应用Gauss 消去法，得到消元后的方程组为
1220.9,0.7.x x x ⎧+=⎨
=⎩
回代求解，可以得到正确的结果.即1
20.2,0.7x x = =.
三、三角分解法 设方程组
Ax b =的系数矩阵A 的顺序主子式不为零.即
1112121222110,1,2,
,.
k
k k k k kk
a a a a a a k n a a a ∆=
≠=
在Gauss 消去法中，第一次消元时，相当于用单位下三角阵
21113111101
0010n m L m m -⎡⎤⎢⎥- ⎢⎥⎢⎥=- ⎢⎥ ⎢⎥⎢⎥- ⎢⎥⎣⎦，
左乘方程组
Ax b =，得
11A x b =，
其中
11121(1)(1)1
22211(1)200n n n nn a a a a a A L a a -(1)
⎡⎤
⎢⎥ ⎢⎥==
⎢⎥ ⎢⎥
⎢⎥ ⎣⎦，
1(1)(1)111,11,1,1(,,,)T
n n n n b L b a a a -+++==.
第二次消元时，相当于用单位下三角阵
1232210101001n L m m - ⎡⎤
⎢⎥ ⎢⎥⎢⎥= - ⎢⎥ ⎢⎥⎢⎥ - ⎢⎥⎣⎦0，
左乘方程组
11A x b =，得
22A x b =
其中
11121(1)(1)22211(2)(2)221333(2)(2)300000n n n n nn a a a a a A L L A a a a a --⎡⎤ ⎢⎥ ⎢⎥⎢⎥== ⎢⎥⎢⎥ ⎢⎥ ⎢⎥⎣⎦，
11(1)(2)
(2)2211,12,13,1,1(,,,
,).T
n n n n n b L L b a a a a --++++==
经过
1n -次消元，最后得到等价方程组
11n n A x b --=
其中
11121(1)22211
11
112
21(1)
n n n n n n nn a a a a a A L L L L A a (1)--------⎡⎤
⎢⎥ ⎢⎥==⎢⎥
⎢⎥
⎢⎥ ⎣⎦
11
11(1)
(1)112
221,12,1,1(,,
,)n T
n n n n n n n b L L L L b a a a --------+++==
注意到
1n A -是一个上三角阵，记
11
11
112
21n n n U A L L L L A -------==
则
12
1()n A L L L U LU -==
其中，121n L L L L -=.
不难验证
21313212_1111
n n nn m L m m m m m ⎡⎤
⎢⎥ ⎢⎥
⎢⎥= ⎢⎥ ⎢⎥⎢⎥ 1 ⎢⎥⎣⎦
是单位下三角阵.
于是解线性方程组Ax b =，就转化为解方程 LUx b =，若令Ux y =
就得到一个与 Ax b =等价的方程组
Ly b Ux y =⎧⎨
=⎩
【定理2.2】 若 A 为 n 阶方阵，且 A 的所有顺序主子式0k ∆≠，1,2,
,k n =.则存在唯一的
一个单位下三角矩阵 L 和一个上三角矩阵 U ，使A LU =.
在上述过程中，若不假设A 的顺序主子式都不为零，只假设A 非奇异，那么Gauss 消去法将不可避免要应用两行对换的初等变换.
第一次消元，将第1行与第1r
行交换，相当于将方程组
Ax b =左乘矩阵11r P
： 1111r r P Ax P b
=
经第一次消元得
1111
1111r r L P Ax L P b
--=
即系数矩阵为
11111r A L P A
-=，
其中
1
10111r P
⎡⎢ ⎢
1= 1 0 1 ⎣0 0
⎤
⎥⎥⎢⎥⎢
⎥⎢⎥⎢⎥⎢
⎥⎢⎥⎢⎥⎢⎥⎢⎥⎦
1 列 1r
列 类似地，经
1n -次消元，有
12
111
1
111,22,11n n n n n r n n r r A L P L P L P A
----------=.
如果预先知道每一个(1,2,,1)
i
ir P i n =-，则在消元之前就全部作交换，得 12
11,2,1,n n n r n r r A P P P A PA
----==，
其中，
12
1
1,2,1,n n n r n r r P P P P ----=.即原方程变为
PAx Pb =
然后再消元，相当于对PA 做三角分解
PA LU =
由以上讨论，可得结论 【定理2.3】 若
A 非奇异，则一定存在排列矩阵 P ，使得 PA 被分解为一个单位下三角阵和一个上三角
1 行
1行
r
阵的乘积，即
PA LU =成立.
这时，原方程组
Ax b = 等价于 PAx Pb =，即等价于求解
LUx Pb =
令
Ux y =
则
Ly Pb =
实际求解时，先解方程组
Ly Pb =,再根据 y 求解 Ux y =,即得原方程组Ax b =的解. 这种求解
方法称为三角分解法.
常用三角分解方法有以下几种. 1．Doolittle 分解方法 假设系数矩阵
A 不需要进行行交换，且三角分解是唯一的. 记
21121110n n l L l l ⎡⎤⎢⎥ ⎢⎥=⎢⎥ ⎢⎥ ⎢⎥⎣⎦,
11121222n n nn u u u u u U u ⎡⎤⎢⎥ ⎢⎥=⎢⎥ ⎢⎥ 0 ⎣⎦ 于是有
11121111212222122221
12111110n n n n n n n n nn a a a u u u u u a a a l l l a a a ⎡⎤ ⎡⎤⎢⎥⎢⎥ ⎢⎥⎢⎥=⎢⎥⎢⎥ ⎢⎥⎢⎥ ⎢
⎥⎢⎥ ⎣⎦⎣⎦nn u ⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥
0 ⎣⎦
从前面讨论A 的LU 分解过程可看出,L 、U 的元素都是用有关的(1)
k ij a -来表示的，而它们的计算较麻烦.
现在给出直接从系数矩阵A ，通过比较等式的两边逐步把L 和U 构造出来的方法，而不必利用Gauss 消去法
的中间结果(1)
k ij a -.
计算步骤： (1) 由
L 阵的第1行分别乘U 阵的各列，先算出U 阵的第1行元素 11,1,2,,j j u a j n = =.
然后，由
L 阵的各行分别去乘U 阵的第1列，算出L 阵的第1列元素
1111/,2,3,
,i i l a a i n = =.
(2)现假设已经算出
U 阵的前1r -行元素，L 阵的前1r -列元素，下面来算U 阵的第r 行元素，L 阵的
第r 列元素.
由
L 阵的第r 行分别乘U 阵的第j 列(,1,,)j r r n =+，得
1
1
r ij rk kj rj
k a l u u -==+∑
所以，得
U 阵的第r 行元素
1
1
,,1,
,r rj rj rk kj k u a l u j r r n
-==- =+∑.
再由
L 阵的第i 行(1,2,,)i r r n =++分别去乘U 阵的第r 列,得
11
r ir ik kr ir rr
k a l u l u -==+∑,
所以，得
L 阵的第r 列元素
1
1
[]/,1,2,
,.
r ir ir ik kr rr k l a l u u i r r n -==- =++∑
取
1,2,,r n =逐步计算，就可完成三角分解A LU =；
(3)解与
Ax b = 等价的方程组
Ly b Ux y =⎧⎨
=⎩
逐次用向前代入过程先解
Ly b = 得
1111
,2,3,
,.
i i i ij j j y b y b l y i n -==⎧⎪⎨
=- =⎪⎩
∑
然后再用逐次向后回代过程解
Ux y =得
1
/,()/,1,2,
,2,1.
n n nn n i i ij j ii j i x y u x y u x u i n n =+=⎧⎪⎨
=- =--⎪⎩
∑
2．Crout 分解方法 仍假设系数矩阵
A 不需要进行行交换，且三角分解是唯一的.即ˆA L
=ˆU .与Doolittle 分解方法的区别在
111212122211n n n n nn a a a a a a a a a ⎡⎤ ⎢⎥ ⎢⎥=⎢⎥ ⎢⎥⎢⎥ ⎣⎦1122ˆˆl l ⎡⎤ 0⎢⎥ ⎢⎥⎢⎥ ⎢⎥⎢
⎥⎣⎦122ˆ1ˆ10n u u ⎡⎤⎢⎥ ⎢⎥⎢⎥ ⎢⎥ 1 ⎣⎦ 比较两边，则可推导出与Doolittle 分解方法类似的公式，不过Crout 分解方法是先算ˆL 的第r 列，然后再算ˆU
的第
r 行.
3．Cholesky 分解方法
若 A 为对称正定矩阵，则有 ˆT U L =，即
1
122
()()T
T T A LDL LD LD LL ===
其中
L 为下三角阵. 进一步展开为
1121111211112122221222221212n n n n n n nn n n nn a a a l l l l a a a l l l l l l l a a a ⎡⎤⎡⎤ ⎢⎥⎢⎥ 0 ⎢⎥⎢⎥=⎢⎥⎢⎥ ⎢⎥⎢⎥ ⎢⎥ ⎢
⎥⎣⎦⎣⎦0nn l ⎡⎤⎢
⎥⎢⎥⎢⎥⎢⎥⎢⎥
⎣⎦ 比较两边对应元素，容易得到
12
1
2
1
()
r rr rr rk k l a l -==-∑ ，
1
1
()/r ir ir ik rk rr
k l a l l l -==-∑ 1,2,,;1,2,,.r n i r r n ==++
Cholesky 分解的优点：不用选主元. 由
2
1
r
rr rk k a l ==∑ 可以看出
||1,2,,.rk l k r ≤=
这表明中间量rk l
得以控制，因此不会产生由中间量放大使计算不稳定的现象. Cholesky 分解的缺点：需要作开方运算. 改进的Cholesky 分解： 改为使用分解
T A LDL =
即
111211211212212221211
11n n n n n n n n nn a a a d l l l d a a a l l d a a a ⎡⎤ 1 ⎡⎤⎡⎤⎢⎥⎢
⎥⎢⎥ 1 1 ⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥ ⎢⎥⎢
⎥⎢⎥ ⎢⎥ ⎣⎦⎣⎦⎣⎦2n l ⎡⎤⎢⎥ ⎢
⎥⎢⎥ ⎢⎥ 1⎣⎦
其中
21ˆl 1ˆn l 2ˆn l ˆnn l 1ˆn u
1
211
1()/r r rr rk k k r ir ir ik k rk r
k d a l d l a l d l d
-=-=⎧=-⎪⎪⎨⎪=-⎪⎩∑∑，1,2,
,;1,2,,.r n i r r n ==++ Cholesky 分解方法或平方根法：应用Cholesky 分解可将
Ax b =分解为两个三角形方程组
T Ly b L x y ⎧= ⎪⎨= ⎪⎩
分别可解得
111111/,()/.i i i ik k ii k y b l y b l y l i n -=⎧=⎪⎨
=-, =2,3,
,
⎪⎩
∑
和
1/,()/1,.n n nn n i i ki k ii k i x y l x y l x l i n n =+⎧=⎪⎨
=-, =--2,,2,1
⎪⎩∑
改进的Cholesky 分解方法或改进的平方根法：应用改进的Cholesky 分解，将方程组Ax b =分解为下面两个
方程组
1,,T Ly b L x D y -= ⎧⎨= ⎩
同理可解得
1111,,2,3,
,.
i i i ik k k y b y b l y i n ==⎧=⎪⎨
=- =⎪⎩
∑
和
1/,
/,1,2,
,2,1.
n n n n i i i ki k k i x y d x y d l x i n n =+⎧=⎪⎨
=- =--⎪⎩
∑
4．解三对角方程组的追赶法
若
()ij n n A a ⨯=满足
1||||,1,2,
,.
n
ii ij j j i
a a i n =≠> =∑
则称
A 为严格对角占优矩阵.若A 满足
1||||,1,2,
,.
n
ii ij j j i
a a i n =≠≥ =∑
且其中至少有一个严格不等式成立，则称A 为弱对角占优矩阵.
现在考虑
Ax d = 的求解，即
11112222211111n n n n n n n n n b c x d a b c x d a b c x d d a b x -----⎡⎤⎡⎤⎡⎤ ⎢⎥⎢⎥⎢⎥ ⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥ = ⎢⎥⎢⎥⎢⎥ ⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥
⎢⎥⎢⎥⎢⎥⎣
⎦⎣⎦⎣⎦ 系数矩阵
A 满足条件
11||||0,
||||||,,0,2,3,, 1.
||||0,
i i i i i n n b c b a c a c i n b a ⎧>>⎪
≥+ ≠=-⎨⎪
>>⎩
采用Crout 分解方法
11112222221111n n n n n n n b c a b c a b c a b βαβγαγα---⎡⎤ ⎡⎤⎢⎥ 1 ⎢⎥⎢⎥ ⎢⎥⎢⎥ = ⎢⎥⎢⎥ ⎢⎥ ⎢⎥ ⎢⎥⎢⎥⎣⎦ ⎣⎦1n β-⎡⎤
⎢⎥⎢⎥⎢⎥ ⎢⎥
1 ⎢⎥
⎢⎥ 1 ⎣⎦
其中，
,,i i i αβγ为待定系数.比较上式两边可得到
111111,;
,,2,3,,;
,2,3,
, 1.
i i i i i i i i i b c a b i n c i n ααβγγβααβ-= == =+ == =-
进而可导出
1111111,2,3,,.,/,,2,3,,./(),2,3,
, 1.
i i i i i i i
i i i i a i n b c b b i n c b i n γαβααββαβ--⎧= =⎪
= =⎪⎨
=- =⎪⎪=- =-⎩
由此可看出，真正需要计算的是(1,2,,1)
i n β=-，而
i α可由,i i b a 和1i β-产生.
因此，实现了
A 的Crout 分解后，求解Ax d =就等价于解方程组
Ly d
Ux y =⎧⎨
=⎩
从而得到解三对角方程组的追赶法公式： (1) 计算
i β的递推公式：
1111/,/(),2,3,
, 1.
i i i i i c b c b i n ββαβ-⎧=⎪⎨=- =-⎪⎩
(2) 解方程组
Ly d =：
11111/()/(),2,3,
,.i i i i i i i y d b y d a y b a i n β--⎧=⎪⎨=-- =⎪⎩
(3) 解方程组
Ux y =：
1,1,2,
,2,1.n n i i i i x y x y x i n n β+⎧=⎪⎨=- =--⎪⎩
追赶法的乘除法次数是
66
n -次.将计算
121
n βββ-→→→及
12n y y y →→
→的过程称之为“追”的过程，将计算方程组Ax d =的解
121n n x x x x -→→
→→的过程称之为“赶”的过程.
四、迭代法 将
Ax b =
改写为一个等价的方程组 x Bx k =+
建立迭代公式 (1)
(),0,1,2,.
i i x Bx k i +=+ =
称矩阵
B 为迭代矩阵.
【定义】 如果对固定的矩阵
B
及向量
k
，对任意初始猜值向量
(0)
x ，迭代公式
(1)()i i +()i
()*
lim i i x x →+∞
=
成立，其中*x 是一确定的向量，它不依赖于(0)x 的选取.则称
此迭代公式是收敛的，否则称为发散的.
如果迭代收敛，则应有
**,x Bx k =+
1. 收敛性
()()*,0,1,2,
i i x x i ε=- =
为第i
步迭代的误差向量.则有
(1)(1)*()*()(),0,1,2,.
x x B x x B i εε++=-=-==
所以，容易推出
()(0),0,1,2,,i i B i εε= =
其中，
(0)
(0)
*
x
x
ε
=-为初始猜值的误差向量.
设n n
B K ⨯∈，lim 0
i i B →+∞
=⇔ ()1B ρ<.
迭代法收敛基本定理： 下面三个命题是等价的 (1) 迭代法
(1)()
i i x Bx k +=+收敛；
(2)
()1B ρ<；
(3) 至少存在一种矩阵的从属范数
⋅
，使
1B <
注：当条件()1B ρ<难以检验时，用1B 或
B ∞等容易求出的范数，检验
11B <或
1
B
∞
<来
作为收敛的充分条件较为方便.
常用迭代法如下. 2．Jacob 迭代 考察线性方程组Ax b =，设A 为非奇异的n 阶方阵，且对角线元素0ii a ≠(1,2,,)i n =.此时，
可将矩阵
A 写成如下形式
A D L U =++, 1122(,,
,)nn D diag a a a =，
21313212000
n n a L a a a a ⎡⎤⎢⎥ ⎢⎥
⎢⎥= ⎢⎥ ⎢⎥⎢⎥ 0 ⎢⎥⎣⎦,12131232000n n a a a a a U ⎡⎤ ⎢⎥ ⎢⎥⎢⎥= 0 ⎢⎥
⎢⎥⎢⎥ ⎢⎥⎣⎦,
建立Jacobi 迭代公式
(1)1()1(),i i x D L U x D b +--=-++
迭代矩阵
11()J B D L U I D A --=-+=-
J B 的具体元素为
112
111122122221200n n J n n nn nn a a a a a a B a a a a a a ⎡⎤ - -⎢⎥⎢⎥⎢⎥- - ⎢⎥=⎢⎥⎢⎥ ⎢⎥⎢⎥- - 0 ⎢⎥⎣⎦ Jacobi 迭代法的分量形式如下
1(1)()()11
1(),j n i i i j
j jm m jm m m m j jj x
b a x a x a -+==+=--∑∑
1,2,,;0,1,2,.j n i = =
3．Gauss-Seidel 迭代
容易看出，在Jacobi 迭代法中，每次迭代用的是前一次迭代的全部分量()(1,2,,)i j
x j n =.实际上，在
计算(1)i j x +时，最新的分量(1)(1)(1)121,,,i i i j x x x +++-已经算出，但没有被利用.事实上，如果Jacobi 迭代收敛，最新算出的分量一般都比前一次旧的分量更加逼近精确解，因此，若在求
(1)
i j x
+时，利用刚刚计算出的
新分量(1)(1)(1)
1
21,,,i i i j x x x
+++-，对Jacobi 迭代加以修改，可得迭代公式
1(1)
(1)()11
1(),j n
i i i j
j jm m jm m m m j jj x
b a x a x a -++==+=--∑∑
1,2,,;0,1,2,.j n i = =
矩阵形式
(1)1()1()(),0,1,2,.i i x D L Ux D L b i +--=-++-+=
1()G B D L U -=--+
注：
（1）两种迭代法均收敛时，Gauss-Seidt 迭代收敛速度更快一些.
（2）但也有这样的方程组，对Jacobi 迭代法收敛，而对Gauss-Seidel 迭代法却是发散的. 【例2.4】 分别用Jacobi 迭代法和Gauss-Seidel 迭代法求解下面的方程组
121232342,46,4 2.
x x x x x x x ⎧- =⎪
-+-=⎨⎪
-+=⎩
初始猜值取0(0,0,0)x =. 【解】 Jacobi 迭代公式为
(1)()1
2(1)()()
213(1)()321(2),4
1(6),0,1,2,
41(2),4i i i i i i i x x x x x i x x +++⎧=+⎪⎪
⎪=++=⎨⎪
⎪=+⎪⎩
迭代计算4次的结果如下 (1)(2)(3)(4)(0.5,1.5,0.5),(0.875,1.75,0.875),(0.938,1.938,0.938),(0.984,1.969,0.984).T T T T x x x x ====
Gauss-Seidel 迭代公式为
(1)()1
2(1)(1)()
213(1)(1)321(2),4
1(6),0,1,2,
41(2),4i i i i i i i x x x x x i x x +++++⎧=+⎪⎪
⎪=++=⎨⎪
⎪=+⎪⎩
迭代计算4次的结果如下
(1)(2)(3)(4)(0.5,1.625,0.9063),(0.9063,1.9532,0.9883),(0.9883,2.0,0.9985),(0.9985,1.999,0.9998).
T T T T x x x x ====
从这个例子可以看到，两种迭代法作出的向量序列
()
{}i x 逐步逼近方程组的精确解
*(1,2,1)T x =,
而且Gauss-Seidel 迭代法收敛速度较快.一般情况下，当这两种迭代法均收敛时，Gauss-Seidt 迭代收敛速度更
3．超松弛迭代法
为了加快迭代的收敛速度，可将Gauss-Seidel 迭代公式改写成
1(1)()(1)()
11(),
j n
i i i i j
j
j jm m jm m m m j
jj x
x b a x a x a -++===+--∑∑ 1,2,,;0,1,2,.j n i = =
并记
1(1)(1)()11(),j n
i i i j
j jm m jm m m m j
jj r
b a x a x a -++===--∑∑
称 (1)i j r + 为 1i + 步迭代的第 j 个分量的误差向量.当迭代收敛时，显然有所有的误差向量
(1)0(),1,2,,.i j r i j n +→→∞=
为了获得更快的迭代公式，引入因子R ω∈，对误差向量 (1)
i j r + 加以修正，得超松弛迭代法（简称SOR 方法）
(1)()(1),0,1,2,.i i i j j j x x r i ω++=+ =
即
1
(1)
()(1)()
1
(),j n
i i i i j
j
j jm m
jm m m m j
jj
x
x b a x
a x a ω
-++===+
--∑∑
1,2,,;0,1,2,.j n i = =
适当选取因子
ω，可望比Gauss-Seidel 迭代法收敛得更快.称ω为松弛因子.特别当1ω=时，SOR 方法就是
Gauss-Seidel 迭代法.
写成矩阵向量形式
(1)1()1()[(1)](),j i x D L D U x D L b ωωωωω+--=+--++
0,1,2,.i =
迭代矩阵为
1()[(1)].B D L D U ωωωω-=+--
实际计算时，大部分是由计算经验或通过试算法来确定opt ω的近似值.所谓试算法就是从同一初始向量出发，
取不同的松驰因子
ω
迭代相同次数(注意：迭代次数不应太少)，然后比较其相应的误差向量
()()i i r b Ax =-(或()(1)i i x x --)，
并取使其范数最小的松弛因子ω作为最佳松弛因子opt ω的近似值.
实践证明，此方法虽然简单，但往往是行之有效的. 4．迭代收敛其它判别方法：
用迭代法收敛基本定理来判断收敛性时，当
n 较大时，迭代矩阵的谱半径计算比较困难，因此，人们试图建立
直接利用矩阵元素的条件来判别迭代法的收敛定理. （1） 若方程组
Ax b =中的系数矩阵A 是对称正定阵，则
Gauss-Seidel 迭代法收敛. 对于SOR 方法，当
02ω<< 时迭代收敛
（2）若
A 为严格对角占优阵，则解方程组 Ax b = 的Jacobi 迭代法，Gauss －Seidel 迭代法均收敛. 对于
SOR 方法，当
01ω<< 时迭代收敛.
【例2.5】 设线性方程组为
121221,32,x x x x ⎧+=-⎪⎨+=⎪⎩
建立收敛的Jacobi 迭代公式和Gauss －Seidel 迭代公式. 【解】 对方程组直接建立迭代公式，其Jacobi 迭代矩阵为
0230J B -⎡⎤=⎢⎥
- ⎣⎦，
显见谱半径()1J B ρ=
>，故Jacobi 迭代公式发散.
同理Gauss －Seidel 迭代矩阵为
0206G B -⎡⎤=⎢⎥ ⎣⎦，
谱半径()61
G B ρ=>，故Gauss －Seidel 选代公式也发散. 若交换原方程组两个方程的次序，得一等价方程组
121232,21,x x x x ⎧+=⎪⎨+=-⎪⎩
其系数矩阵显然对角占优，故对这一等价方程组建立的Jacobi 迭代公式，Gauss －Seidel 迭代公式皆收敛. （3）SOR 方法收敛的必要条件是 02ω<<
【定理2.5】 如果A 是对称正定阵，且02ω<<，则解Ax b =的SOR 方法收敛.
注：当
(0,2)ω∈ 时，并不是对任意类型的矩阵A ，解线性方程组Ax b =的SOR 方法都是收敛的.
当SOR 方法收敛时，通常希望选择一个最佳的值opt ω使SOR 方法的收敛速度最快.然而遗憾的是，目前尚无确
定最佳超松弛因子
opt ω的一般理论结果.实际计算时，大部分是由计算经验或通过试算法来确定opt ω的近似值.
所谓试算法就是从同一初始向量出发，取不同的松驰因子
ω迭代相同次数(注意：迭代次数不应太少)，然后比
较其相应的误差向量()()i i r b Ax =-(或()(1)i i x x --)，并取使其范数最小的松弛因子ω作为最佳松
弛因子opt ω的近似值.实践证明，此方法虽然简单，但往往是行之有效的.
【例2.6】 求解线性方程组
Ax b =，其中
10.3000900.308980.30009100.4669110.274710.30898A - -- -0.46691 0= - -- 00.274711(5.32088,6.07624,8.80455,2.67600).T b ⎡⎤⎢⎥
⎢
⎥⎢⎥⎢⎥ - ⎣⎦ =-
分别利用Jacobi 迭代法，Gauss －Seidel 迭代法，SOR 迭代法求解. 【解】其结果列入下表中，方程组精确解(五位有效数字)为
*(8.4877,6.4275, 4.7028,4.0066).T x =-
Jacobi 迭代法计算结果
i
()1i x
()2i x
()3i x
()4i x
()
2||||i r
0 0
12.3095 1 5.3209 6.0762 －8.8046 2.6760 5.3609 2 7.9711 3.5621 －5.2324 1.9014
3.6318
20 8.4872 6.4263 －4.7035 4.0041 0.0041 21
8.4860
6.4271 －4.7050 4.0063 0.0028
Gauss-Seidel 迭代法计算结果
i
()1i x
()2
i x
()3
i x
()4
i x
()
2||
||i
r
0 0
12.3095 1 5.3209 7.6730 －5.2220 2.8855 3.6202 2 8.5150 6.1933 －5.1201 3.9004 0.4909 8 8.4832 6.4228 －4.7064 4.0043 0.0078 9
8.4855
6.4252
－4.7055
4.0055
0.0038
SOR 迭代法计算结果（
1.16ω=）
i
()1i x
()2
i x
()3
i x
()4
i x
()
2||
||i r
0 0
12.3095 1 6.1722 9.1970 －5.2320 3.6492 3.6659 2 9.6941 6.1177 －4.8999 4.4335 1.3313 6 8.4842 6.4253 －4.7005 4.4047 0.0051 7
8.4868
6.4288
－4.7031
4.0065
0.0016
计算结果表明，若求出精确到小数点后两位的近似解，Jacobi 迭代法需要21次，Gauss －Seidel 迭代法需要9次，而SOR 迭代法(选松弛因子 1.16ω=)仅需要7次，起到加速作用.
5．误差分析 【定理
2.6】设 *x 是方程 Ax b = 的惟一解，
v ⋅ 是某一种向量范数，若对应的迭代矩阵其范数
1
v B <，则迭代法
(1)
()
,0,1,2,.i i x
Bx k i +=+ = 收敛，且产生向量序列
(){}i x 满足
()*()(1)||||||||||||1||||i i i v
v v
v
B x x x x B --≤
--
()
*
(1)(0)||||||||||||1||||i i v
v v
v
B x x x x B -≤--
【证明】 由迭代收敛基本定理的(3)知，迭代法
(1)(),0,1,2,.i i x Bx k i +=+ =
收敛到方程的解*
x .于是，由迭代公式立即得到
(1)*()*(1)()()(1)(),().i i i i i i x x B x x x x B x x ++--=--=-
为书写方便把
v 范数中v 略去，有估计式
(1)*()*||||||||||||,i i x x B x x +-≤⋅-
(1)()()(1)||||||||||||.i i i i x x B x x +--≤⋅-
再利用向量范数不等式
||||||||||||x y x y -≥-
于是得第一个不等式
()(1)(1)()()*(1)*()*||||||||||||||||||||(1||||)||||,
i i i i i i i B x x x x x x x x B x x -++ -≥-≥--- ≥--
再反复递推即第二个不等式.
注：
（1）若事先给出误差精度
ε，利用第二个不等式可得到迭代次数的估计
(1)(0)
(1||||)ln ln ||||||||v v v B i B x x ε⎡⎤->⎢⎥-⎣⎦ （2）在
||||v B 不太接近1的情况下，由第一个不等式，可用
()(1)||||i i v x x ε--<
作为控制迭代终止的条件，并取 ()
i x 作为方程组 Ax b = 的近似解.但是在
||||v B 很接近1时，此方法并不可靠.一般可取
1,2,v =∞或F .
【例2.7】 用Jacobi 迭代法解方程组
123123123202324,812,231530.
x x x x x x x x x ⎧++=⎪
++=⎨⎪
-+=⎩
问Jacobi 迭代是否收敛?若收敛，取(0)(0,0,0)T x =，需要迭代多少次，才能保证各分量的误差绝对值小于
610-?
【解】
Jacobi 迭代的分量公式为
(1)()()
1
23(1)()()
213(1)()()3
121(2423)20
1(12),0,1,2,
81(3022),15i i i i i i i i i x x x x x x i x x x +++⎧=--⎪⎪
⎪=-- =⎨⎪
⎪=-+⎪⎩
Jacobi 迭代矩阵
J B 为
130102011
08821
0155J B ⎡⎤ - -⎢⎥⎢⎥⎢⎥
=- -⎢⎥⎢⎥⎢⎥- ⎢⎥⎣⎦，
由
5251||||max ,,1
208153J B ∞⎧⎫==<⎨⎬⎩⎭知，Jacobi 迭代收敛. 因设(0)(0,0,0)T
x =，用迭代公式计算一次得
(1)
(1)(1)1
2363,, 2.
52x x x = = =
而
(1)(0)|||| 2.x x ∞-=
于是有
6
110(1)13ln ln 13.23i -⎡⎤⋅-⎢⎥>=⎢⎥
⎢⎥
⎣⎦
所以，要保证各分量误差绝对值小于6
10-，需要迭代14次.
【例2.8】 用Gauss －Seidel 迭代法解例2.11中的方程组，问迭代是否收敛?若收敛，取(0)(0,0,0)T
x =，需要迭代多少次，才能保证各分量误差的绝对值小于6
10-?
【解】 Gauss －Seidel 迭代矩阵
G B 为
102403601()03025524000G B D L U - - ⎡⎤⎢
⎥
=-+= -⎢
⎥⎢⎥ 38 -3⎣⎦
显然
1
||||1
4G B =<，所以迭代收敛. Gauss －Seidel 迭代分量公式为
(1)()()
123(1)(1)()
213(1)(1)(1)3121(2423),201(12),0,1,2,
81(3022),15i i i i i i i i i x x x x x x i x x x ++++++⎧=--⎪⎪
⎪=-- =⎨⎪
⎪=-+⎪⎩
因取(0)(0,0,0)T x =，故迭代一次得
(1)(1)(1)1231.2, 1.35, 2.11x x x = = =
于是有
(1)(0)|||| 2.11x x ∞-=，
计算得
6
110(1)14ln
ln 10.2.114i -⎡⎤⋅-⎢⎥>=⎢⎥⎢⎥⎣⎦
所在，要保证各分量误差绝对值小于6
10-，需要迭代11次.。