第 5章 解线性方程组的直接法
实际中,存在大量的解线性方程组的问题。很多数值方
法到最后也会涉及到线性方程组的求解问题:如样条插值的
M和 m关系式,曲线拟合的法方程,方程组的 Newton迭代等
问题。
?
?
?
?
?
???
???
nnnnn
nn
bxaxa
bxaxa
?
?
?
11
11111
0)det( ?A
对线性方程组:
或者,bAx?
我们有 Gram法则:当且仅当 时,有唯一的解,而且解为:
?
?
?
?
?
?
?
?
?
?
???
??
??
nnninnin
nii
i
i
i
aabaa
aabaa
DAD
D
D
x
??
?????
??
111
11111111
d et),d et (,
但 Gram法则不能用于计算方程组的解,如 n= 100,1033次 /秒的计算机要算 10120年
解线性方程组的方法可以分为 2类:
① 直接法,准确,可靠,理论上得到的解是精确的
② 迭代法,速度快,但有误差
本章讲解直接法
5.1 消元法
我们知道,下面有 3种方程的解我们可以直接求出:
niabxaaad i a gA
ii
i
inn,,1,),,,( 2211 ?? ????
① n次运算
ni
l
xlb
x
lll
ll
l
A
ii
i
j
jiji
i
nnnn
,,1,
1
1
21
2221
11
?
?
???
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
② ( n+ 1) n/2次运算
1,,,
1222
11211
?
??
?
?
ni
u
xub
x
u
uu
uuu
A
ii
n
ij
jiji
i
nn
n
n
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
③ ( n+ 1) n/2次运算
对方程组,作如下的变换,解不变
① 交换两个方程的次序
② 一个方程的两边同时乘以一个非 0的数
③ 一个方程的两边同时乘以一个非 0数,加到另一个方程
因此,对应的对增广矩阵 (A,b),作如下的变换,解不变
① 交换矩阵的两行
② 某一行乘以一个非 0的数
③ 某一个乘以一个非 0数,加到另一行
消元法 就是对增广矩阵作上述行的变换,变为我们已知的 3种类型之一,而后求根
? 高斯消元法:
思
路
首先将 A化为上三角阵,再回代求解 。
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nnnnn
n
n
baaa
baaa
baaa
?
?????
?
?
21
222221
111211
= ( 1 ) ( 1 ) ( 1 ) ( 1 ) ( 1 )
1 1 1 2 1 3 1 1
( 2 ) ( 2 ) ( 2 ) ( 2 )
2 2 2 3 2 2
( 3 ) ( 3 ) ( 3 )
3 3 3 3
( ) ( )
0
00
0 0 0
n
n
n
nn
n n n
a a a a b
a a a b
a a b
ab
??
??
??
??
??
步骤如下:
第一步:
nii
a
a i,,2,1
11
1 ????? 行第行第
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nnnnn
n
n
baaa
baaa
baaa
?
?????
?
?
21
222221
111211
?
?
?
?
?
?
?
?
?
?
?
?
?
?
)2()2()2(
2
)2(
2
)2(
2
)2(
22
111211
0
0
nnnn
n
n
baa
baa
baaa
?
?????
?
?
运算量,(n-1)*(1+n)
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
)3()3()3(
3
)3(
3
)3(
3
)3(
33
)2(
2
)2(
2
)2(
23
)2(
22
11131211
00
00
0
nnnn
n
n
n
baa
baa
baaa
baaaa
?
??????
?
?
?
运算量,(n-2)*(1+n-1)=(n-2)n
第二步:
niia a i,,3,2 )2(
22
)2(
2 ????? 行第行第
?
?
?
?
?
?
?
?
?
?
?
?
?
?
)2()2()2(
2
)2(
2
)2(
2
)2(
22
111211
0
0
nnnn
n
n
baa
baa
baaa
?
?????
?
?
第 k步:
nkii
a
a
k
kk
k
ik,,1,k
)(
)(
?????? 行第行第
类似的做下去,我们有:
运算量,(n- k)*(1+ n- k+ 1)=(n- k)(n- k+ 2)
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
)()(
)3(
3
)3(
3
)3(
33
)2(
2
)2(
2
)2(
23
)2(
22
11131211
000
00
0
n
n
n
nn
n
n
n
ba
baa
baaa
baaaa
?
??????
?
?
?n- 1步以后,我们可以得到变换后的矩阵为:
因此,总的运算量为:
??
?
???
1
1
)2)((
n
k
knkn
加上 解上述上三角阵的运算量 (n+1)n/2,总共为:
)(33 32
3
nOnnn ???
注意到,计算过程中 )(k
kka
处在被除的位置,因此整个计算过程要保证它不为 0
所以,Gauss消元法的可行条件为:
0)( ?kkka
就是要求 A的所有顺序主子式均不为 0,即
ni
aa
aa
iii
i
,,1,0d et
1
111
?
?
???
?
??
?
?
?
?
?
?
?
?
?
?
因此,有些有解的问题,不能用 Gauss消元求解
另外,如果某个 )(k
kka
很小的话,会引入大的误差
例,单精度解方程组
??
?
??
???
2
110
21
21
9
xx
xx
/* 精确解为 和 */.,,1000.,,00.1
101
1
91
???
??? ?x
8个
.,,8 9 99.,,99.02 12 ?????? xx
8个
用 Gaussian 消元法计算:
9112121 10/ ?? aam
? 999
2122 10101010.,,0.011 ???????? ?ma
8个
9212 1012 ????? ?mb
?????? ???
?
99
9
10100
1110
0,1 12 ??? xx
小主元可能导致计算
失败。
2,列主元消元法
在 Gauss消元第 k步之前,做如下的事情:
||||m a x )()( kjkkiknik aa ???若 交换 k行和 j行
行的交换,不改变方程组的解,同时又有效地克服了 Gauss消元地缺陷
例:
??????
?
211
1110 9
1,1 12 ??? xx
??????? 110
211
??????? ? 1110
211
9
?
3,Gauss-Jordan消元法
将在 Gauss消元第 k步,变为
nkkii
a
a
k
kk
k
ik,,1,1,,1,k
)(
)(
?? ?????? 行第行第
将该行上三角地部分也变为 0最后变为一个对角阵。
它地运算次数比 Gauss消元多。使用于计算多个系数一样地方程组,如
BAX ? X,B均为矩阵
5.2 直接分解法
Gauss消元法的第 k步:
nkiia a k
kk
k
ik,,1,k
)(
)(
?????? 行第行第
从矩阵理论来看,相当于左乘矩阵
nki
a
a
l
l
l
k
kk
k
ikk
ik
k
nk
k
kk
,,1,,
1
1
1
1
)(
)(
)(
)(
)(
1
?
???
????
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
因此,整个 Gauss消元法相当于左乘了一个单位下三角阵
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
1
1
1
1
)(
1
)(
1
)(
1
21
n
nn
k
nkn
k
kk
lll
l
l
L
??
??
?
?
所以有 ULAtsL ??,,L为单位下三角阵,U为上三角阵
LUAtsL ???,,
因此
??
?
?
?????
yUx
bLybL U xbAx
我们可以通过 2次反代过程求解方程组
注意:
分解的理论由 Gauss消元得出,因此分解能够进行的条件与 Gauss消元一样
1,Doolittle分解
L为下三角,U为单位上三角
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nn
n
n
nnnnnn
n
n
u
uu
uuu
ll
l
aaa
aaa
aaa
??
?
?
?
???
?
????
?
?
222
11211
21
21
21
22221
11211
1
1
1
比较第 1行:
jjjj aunjua 1111,,1 ???? ?
比较第 1列:
11
1
11111,,2 u u
alnila i
iii ???? ?
比较第 2行:
jjjjjj ulaunjuula 1211221212,,2 ?????? ?
比较第 2列:
22
1212
22221212
u,,3 uu
u
lalnilla ii
iiii
?????? ?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nn
n
n
nnnnnn
n
n
u
uu
uuu
ll
l
aaa
aaa
aaa
??
?
?
?
???
?
????
?
?
222
11211
21
21
21
22221
11211
1
1
1
比较第 k行,?? ?
?
?
?
??????
1
1
1
1
,,
k
r
rjkrkjkjkj
k
r
rjkrkj ulaunkjuula ?
比较第 k列:
kk
k
r
rkirik
ikkkik
k
r
rkirik u
ula
lnkiulula
?
?
?
?
?
?
?
??????
1
1
1
1
,,1 ?
K-1次
K-1+ 1次
分解过程完毕,加上两次反代过程
1,,,
,,1,
1
1
1
?
?
ni
u
xuy
x
niylby
ii
n
ij
jiji
i
i
j
jijii
?
?
?
???
?
?
??
?
?
总运算量为:
? ? 332 )1(2 )1()()1)(1( 231
1
nnnnnnnkknkknn
k
?????????????
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nnnn
n
n
ull
uul
uuu
?
????
?
?
21
22221
11211
存储在矩阵的原来位置,且不影响计算
2,Courant 分解
L为下三角,U为单位上三角
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
1
1
1
2
112
21
2221
11
21
22221
11211
??
?
?
?
???
?
????
?
?
n
n
nnnnnnnn
n
n
u
uu
lll
ll
l
aaa
aaa
aaa
两次反代过程
1,,,
,,1,
1
1
1
?
?
nixuyx
ni
l
ylb
y
n
ij
jijii
ii
i
j
jiji
i
???
?
?
?
?
?
??
?
?
下面,我们对一下 特殊 的矩阵,提出一些特定的分解法
比较第 k列,?? ?
?
?
?
??????
1
1
1
1
,,
k
r
rkirikikik
k
r
rkirik ulalnkilula ?
比较第 k行:
kk
k
r
rjkrkj
kjkjkk
k
r
rjkrkj l
ula
unkjulula
?
?
?
?
?
?
?
??????
1
1
1
1
,,1 ?
3,三对角阵的 追赶法
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
1
1
1
1
1
22
1
1
22
11
n
nnnn
n
ac
b
ac
ba
?
?
??
??
?
?
?
????
?
?
?
?
?
?
?
??
????
??
?
ni
b
cnica
nic
i
i
i
iiii
ii
,,1,
0,,,1,
,,2,
11
?
?
?
?
?
??
?
??
?
?
?
????
???
?
?
)0( 1,,,
,,1,)(
1
1
niiii
i
iii
i
nixyx
niycfy
??
?
?
?
所以,有计算过程如下:
1,,,
,,1
)(
1
1
1
?
?
nkxyx
ni
ycf
y
b
ca
kkkk
i
iii
i
i
i
i
iiii
???
?
?
?
?
?
?
?
?
?
?
?
??
?
?
?
?
?
?
?
??
3,对称正定阵的 LDLT分解
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nn
n
n
nnnnnnnn
n
n
T
l
ll
lll
lll
ll
l
aaa
aaa
aaa
LLA
??
?
?
?
???
?
????
?
?
222
12111
21
2221
11
21
22221
11211
若 A对称正定,则有下三角整 L,使得
所以有:
?
?
?
?
?
?
?
??
??
?
???
?
?
?
?
nki
l
lla
l
lal
kk
k
r
krikik
ik
k
r
krkkkk
,,1,
)(
)(
1
1
2
11
1
2
?
称为平方根法,
因为带了开方运算,
因此不常用
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nnnnnnnn l
l
l
ll
l
lll
ll
l
?
?
???
?
???
22
11
21
21
21
2221
11
1''
1'
1又
L~ D~
则有
TTT L D LLDDLA ?? ~~~~
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nnn d
d
d
D
ll
l
L
?
?
???
2
1
21
21
1
1
1
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
n
n
n
nnnnnn
n
n
d
ldd
ldldd
ll
l
aaa
aaa
aaa
??
?
?
?
???
?
????
?
?
222
112111
21
21
21
22221
11211
1
1
1
?
?
?
?
?
?
?
??
??
?
???
?
?
?
?
nki
d
dlla
l
dlad
k
k
r
rkrikik
ik
k
r
rkrkkk
,,1,
)(
1
1
1
1
2
?
比较等号两边后,有
1 21 31
21 2 32
31 32 3
6 3 2 1 1
3 5 1 1 1
2 1 6 1 1
d l l
A l d l
l l d
?? ? ? ? ? ? ? ?
? ? ? ? ? ? ? ???
? ? ? ? ? ? ? ?
? ? ? ? ? ? ? ?
? ? ? ? ? ? ? ?
1 11 6da? ? ? 2 1 2 1 1/ 1 / 2l a d? ? ? 31 31 1/ 1 / 3l a d? ? ?
2 2 2 2 1 1 2 1 1 3 / 2d a l d l? ? ?32 32 1 31 21 2( ) / 4 / 13l a d l l d? ? ?
3 3 3 3 1 1 3 1 3 2 2 3 2
26
39d a l d l l d l? ? ? ?
为了提高数值稳定性,可考虑列主元三角分解法,设已完成 A=LU的 k-1步分解计算,
矩阵分解成
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
????
nnnk2n1n
knkk2k1k
n1kk1k12k11k
n2k22221
n1k11211
aa
aa
uu
uuu
uuuu
??
????
??
??
?????
??
??
ll
ll
ll
l
设 )ua(m a xua 1k 1j jktjtkntk1k 1j jkijik ?? ?????? ??? ll令 rk?ri
相当于取 ????? 1k 1j jkijikkk uau l为第 k步分解的主元素,
但要注意方程组的常数项也要相应变换,
定义, 设 ‖ ?‖ 是一种向量范数
AxxAxA
xRxxRx nn 1,0,
s ups up
????
??
称之为由向量范数派生的 矩阵算子范数,
5.3 矩阵范数和条件数
定义 设 ‖ ?‖ 是以 n阶方阵为变量的实值函数,且满足条件,
(1)非负性, ‖ A‖ ?0,且 ‖ A‖ =0当且仅当 A=0
(2)齐次性, ‖ ?A‖ =|? |‖ A‖,??R
(3)三角不等式,‖ A+B‖ ?‖ A‖ +‖ B‖
(4)三角不等式,‖ AB‖ ?‖ A‖‖ B‖
则称 ‖ A‖ 为矩阵 A的 范数,
对应于 3种常见的向量范数,有 3种矩阵范数
?
???
?
n
i
ijnj aA
111
m ax 列和的最大值
?
????
?
n
j
ijni aA
11
m ax
行和的最大值
12 ??A
是 ATA的最大特征值,也称为谱范数
1?
矩阵范数的一些性质:
①
②
③
④
⑤
00,&,0 ???? AAA
RAA ???? ???,
nRBABABA ?????,,
nRBABAAB ????,,
nRxxAAx ????,
定理,若 ? 为 A 的特征值,则 A??
证:
xAx ????
xAx ????
xAxAxx ??????? ??
????
x
xAA
x为 A的特征值
#证毕
定义 5.2,谱半径
rnrA ?? ??? 1m a x)(
易知:
AA ?)(?
5.4 条件数和病态矩阵
定义 5.3:(条件数)
ppp AAAC on d
1)( ??? p?
表示某种范数
设 bAx?, A 引入误差 后A?,解引入误差 x?,则
bxxAA ??? ))(( ??
xAxAAxbxAA ????????? ???? )(
xAAAx ???? ? ??? 1)(
AAA
x
x ??? 1)( ???
AAAAI ?? 111 )( ?????
? ? 111 )()( ??? ??? AAAAIA ??
AAAAI ?? ???? ??? 111 )(
AAAAI ?? 111 )( ??? ???
注意到
BBI ???
?
1
1)( 1
因为:
BDDDBIIBID ???????? ? )(1 )( 1
)1( BDDBD ?????
AA
AA
x
x
?
?? 11
1
1
?
?
?
????
A
A
AA
A
A
AA
AA
AA
?
?
?
?
??
?
?
??
?
?
?
?
?
?
1
1
1
1
1
1
条件数
很小
A
A
AA
?
?? ? 1
条件数表示了对 误差的放大率
同样,类似有
b
b
AA
x
x
bxA
bbxxA
??
??
??
????
??
???
? 1
)(
注,一般判断矩阵是否病态,并不计算 A?1,而由经验
得出。
?行列式很大或很小(如某些行、列近似相关);
?元素间相差大数量级,且无规则;
?主元消去过程中出现小主元;
?特征值相差大数量级。
精确解 为,
1
1??
?
???
?
??x?
例
???
?
???
??
???
?
???
??
97.1
99.1,
98.099.0
99.01 bA ?
计算 cond (A)2 。
???????? ?? 100 00990 0 990 0980 0
A?1 =
解,考察 A 的特征根
??? 0)d e t ( AI?
0 0 0 0 5 0 5 0 4.0
9 8 0 0 5 0 5 0 4.1
2
1 ?????
??
2
1
2)( ?
?Ac o n d 39206 >> 1
?测试病态程度:
给 一个扰动b?
???
?
???
?
?
???
?
?
3
4
10106.0
1097.0b??,其相对误差为
%01.010513.0|||| |||| 4
2
2 ??? ?
b
b? ?? 此时 精确解 为
???
?
???
?
?? 020 3.1
3*x?
???
?
???
?
???? 0203.2
2* xxx ???? ??
2
2
||||
||||
x
x??? 2.0102 > 200%
为对称矩阵
实际中,存在大量的解线性方程组的问题。很多数值方
法到最后也会涉及到线性方程组的求解问题:如样条插值的
M和 m关系式,曲线拟合的法方程,方程组的 Newton迭代等
问题。
?
?
?
?
?
???
???
nnnnn
nn
bxaxa
bxaxa
?
?
?
11
11111
0)det( ?A
对线性方程组:
或者,bAx?
我们有 Gram法则:当且仅当 时,有唯一的解,而且解为:
?
?
?
?
?
?
?
?
?
?
???
??
??
nnninnin
nii
i
i
i
aabaa
aabaa
DAD
D
D
x
??
?????
??
111
11111111
d et),d et (,
但 Gram法则不能用于计算方程组的解,如 n= 100,1033次 /秒的计算机要算 10120年
解线性方程组的方法可以分为 2类:
① 直接法,准确,可靠,理论上得到的解是精确的
② 迭代法,速度快,但有误差
本章讲解直接法
5.1 消元法
我们知道,下面有 3种方程的解我们可以直接求出:
niabxaaad i a gA
ii
i
inn,,1,),,,( 2211 ?? ????
① n次运算
ni
l
xlb
x
lll
ll
l
A
ii
i
j
jiji
i
nnnn
,,1,
1
1
21
2221
11
?
?
???
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
② ( n+ 1) n/2次运算
1,,,
1222
11211
?
??
?
?
ni
u
xub
x
u
uu
uuu
A
ii
n
ij
jiji
i
nn
n
n
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
③ ( n+ 1) n/2次运算
对方程组,作如下的变换,解不变
① 交换两个方程的次序
② 一个方程的两边同时乘以一个非 0的数
③ 一个方程的两边同时乘以一个非 0数,加到另一个方程
因此,对应的对增广矩阵 (A,b),作如下的变换,解不变
① 交换矩阵的两行
② 某一行乘以一个非 0的数
③ 某一个乘以一个非 0数,加到另一行
消元法 就是对增广矩阵作上述行的变换,变为我们已知的 3种类型之一,而后求根
? 高斯消元法:
思
路
首先将 A化为上三角阵,再回代求解 。
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nnnnn
n
n
baaa
baaa
baaa
?
?????
?
?
21
222221
111211
= ( 1 ) ( 1 ) ( 1 ) ( 1 ) ( 1 )
1 1 1 2 1 3 1 1
( 2 ) ( 2 ) ( 2 ) ( 2 )
2 2 2 3 2 2
( 3 ) ( 3 ) ( 3 )
3 3 3 3
( ) ( )
0
00
0 0 0
n
n
n
nn
n n n
a a a a b
a a a b
a a b
ab
??
??
??
??
??
步骤如下:
第一步:
nii
a
a i,,2,1
11
1 ????? 行第行第
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nnnnn
n
n
baaa
baaa
baaa
?
?????
?
?
21
222221
111211
?
?
?
?
?
?
?
?
?
?
?
?
?
?
)2()2()2(
2
)2(
2
)2(
2
)2(
22
111211
0
0
nnnn
n
n
baa
baa
baaa
?
?????
?
?
运算量,(n-1)*(1+n)
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
)3()3()3(
3
)3(
3
)3(
3
)3(
33
)2(
2
)2(
2
)2(
23
)2(
22
11131211
00
00
0
nnnn
n
n
n
baa
baa
baaa
baaaa
?
??????
?
?
?
运算量,(n-2)*(1+n-1)=(n-2)n
第二步:
niia a i,,3,2 )2(
22
)2(
2 ????? 行第行第
?
?
?
?
?
?
?
?
?
?
?
?
?
?
)2()2()2(
2
)2(
2
)2(
2
)2(
22
111211
0
0
nnnn
n
n
baa
baa
baaa
?
?????
?
?
第 k步:
nkii
a
a
k
kk
k
ik,,1,k
)(
)(
?????? 行第行第
类似的做下去,我们有:
运算量,(n- k)*(1+ n- k+ 1)=(n- k)(n- k+ 2)
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
)()(
)3(
3
)3(
3
)3(
33
)2(
2
)2(
2
)2(
23
)2(
22
11131211
000
00
0
n
n
n
nn
n
n
n
ba
baa
baaa
baaaa
?
??????
?
?
?n- 1步以后,我们可以得到变换后的矩阵为:
因此,总的运算量为:
??
?
???
1
1
)2)((
n
k
knkn
加上 解上述上三角阵的运算量 (n+1)n/2,总共为:
)(33 32
3
nOnnn ???
注意到,计算过程中 )(k
kka
处在被除的位置,因此整个计算过程要保证它不为 0
所以,Gauss消元法的可行条件为:
0)( ?kkka
就是要求 A的所有顺序主子式均不为 0,即
ni
aa
aa
iii
i
,,1,0d et
1
111
?
?
???
?
??
?
?
?
?
?
?
?
?
?
?
因此,有些有解的问题,不能用 Gauss消元求解
另外,如果某个 )(k
kka
很小的话,会引入大的误差
例,单精度解方程组
??
?
??
???
2
110
21
21
9
xx
xx
/* 精确解为 和 */.,,1000.,,00.1
101
1
91
???
??? ?x
8个
.,,8 9 99.,,99.02 12 ?????? xx
8个
用 Gaussian 消元法计算:
9112121 10/ ?? aam
? 999
2122 10101010.,,0.011 ???????? ?ma
8个
9212 1012 ????? ?mb
?????? ???
?
99
9
10100
1110
0,1 12 ??? xx
小主元可能导致计算
失败。
2,列主元消元法
在 Gauss消元第 k步之前,做如下的事情:
||||m a x )()( kjkkiknik aa ???若 交换 k行和 j行
行的交换,不改变方程组的解,同时又有效地克服了 Gauss消元地缺陷
例:
??????
?
211
1110 9
1,1 12 ??? xx
??????? 110
211
??????? ? 1110
211
9
?
3,Gauss-Jordan消元法
将在 Gauss消元第 k步,变为
nkkii
a
a
k
kk
k
ik,,1,1,,1,k
)(
)(
?? ?????? 行第行第
将该行上三角地部分也变为 0最后变为一个对角阵。
它地运算次数比 Gauss消元多。使用于计算多个系数一样地方程组,如
BAX ? X,B均为矩阵
5.2 直接分解法
Gauss消元法的第 k步:
nkiia a k
kk
k
ik,,1,k
)(
)(
?????? 行第行第
从矩阵理论来看,相当于左乘矩阵
nki
a
a
l
l
l
k
kk
k
ikk
ik
k
nk
k
kk
,,1,,
1
1
1
1
)(
)(
)(
)(
)(
1
?
???
????
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
因此,整个 Gauss消元法相当于左乘了一个单位下三角阵
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
1
1
1
1
)(
1
)(
1
)(
1
21
n
nn
k
nkn
k
kk
lll
l
l
L
??
??
?
?
所以有 ULAtsL ??,,L为单位下三角阵,U为上三角阵
LUAtsL ???,,
因此
??
?
?
?????
yUx
bLybL U xbAx
我们可以通过 2次反代过程求解方程组
注意:
分解的理论由 Gauss消元得出,因此分解能够进行的条件与 Gauss消元一样
1,Doolittle分解
L为下三角,U为单位上三角
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nn
n
n
nnnnnn
n
n
u
uu
uuu
ll
l
aaa
aaa
aaa
??
?
?
?
???
?
????
?
?
222
11211
21
21
21
22221
11211
1
1
1
比较第 1行:
jjjj aunjua 1111,,1 ???? ?
比较第 1列:
11
1
11111,,2 u u
alnila i
iii ???? ?
比较第 2行:
jjjjjj ulaunjuula 1211221212,,2 ?????? ?
比较第 2列:
22
1212
22221212
u,,3 uu
u
lalnilla ii
iiii
?????? ?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nn
n
n
nnnnnn
n
n
u
uu
uuu
ll
l
aaa
aaa
aaa
??
?
?
?
???
?
????
?
?
222
11211
21
21
21
22221
11211
1
1
1
比较第 k行,?? ?
?
?
?
??????
1
1
1
1
,,
k
r
rjkrkjkjkj
k
r
rjkrkj ulaunkjuula ?
比较第 k列:
kk
k
r
rkirik
ikkkik
k
r
rkirik u
ula
lnkiulula
?
?
?
?
?
?
?
??????
1
1
1
1
,,1 ?
K-1次
K-1+ 1次
分解过程完毕,加上两次反代过程
1,,,
,,1,
1
1
1
?
?
ni
u
xuy
x
niylby
ii
n
ij
jiji
i
i
j
jijii
?
?
?
???
?
?
??
?
?
总运算量为:
? ? 332 )1(2 )1()()1)(1( 231
1
nnnnnnnkknkknn
k
?????????????
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nnnn
n
n
ull
uul
uuu
?
????
?
?
21
22221
11211
存储在矩阵的原来位置,且不影响计算
2,Courant 分解
L为下三角,U为单位上三角
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
1
1
1
2
112
21
2221
11
21
22221
11211
??
?
?
?
???
?
????
?
?
n
n
nnnnnnnn
n
n
u
uu
lll
ll
l
aaa
aaa
aaa
两次反代过程
1,,,
,,1,
1
1
1
?
?
nixuyx
ni
l
ylb
y
n
ij
jijii
ii
i
j
jiji
i
???
?
?
?
?
?
??
?
?
下面,我们对一下 特殊 的矩阵,提出一些特定的分解法
比较第 k列,?? ?
?
?
?
??????
1
1
1
1
,,
k
r
rkirikikik
k
r
rkirik ulalnkilula ?
比较第 k行:
kk
k
r
rjkrkj
kjkjkk
k
r
rjkrkj l
ula
unkjulula
?
?
?
?
?
?
?
??????
1
1
1
1
,,1 ?
3,三对角阵的 追赶法
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
1
1
1
1
1
22
1
1
22
11
n
nnnn
n
ac
b
ac
ba
?
?
??
??
?
?
?
????
?
?
?
?
?
?
?
??
????
??
?
ni
b
cnica
nic
i
i
i
iiii
ii
,,1,
0,,,1,
,,2,
11
?
?
?
?
?
??
?
??
?
?
?
????
???
?
?
)0( 1,,,
,,1,)(
1
1
niiii
i
iii
i
nixyx
niycfy
??
?
?
?
所以,有计算过程如下:
1,,,
,,1
)(
1
1
1
?
?
nkxyx
ni
ycf
y
b
ca
kkkk
i
iii
i
i
i
i
iiii
???
?
?
?
?
?
?
?
?
?
?
?
??
?
?
?
?
?
?
?
??
3,对称正定阵的 LDLT分解
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nn
n
n
nnnnnnnn
n
n
T
l
ll
lll
lll
ll
l
aaa
aaa
aaa
LLA
??
?
?
?
???
?
????
?
?
222
12111
21
2221
11
21
22221
11211
若 A对称正定,则有下三角整 L,使得
所以有:
?
?
?
?
?
?
?
??
??
?
???
?
?
?
?
nki
l
lla
l
lal
kk
k
r
krikik
ik
k
r
krkkkk
,,1,
)(
)(
1
1
2
11
1
2
?
称为平方根法,
因为带了开方运算,
因此不常用
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nnnnnnnn l
l
l
ll
l
lll
ll
l
?
?
???
?
???
22
11
21
21
21
2221
11
1''
1'
1又
L~ D~
则有
TTT L D LLDDLA ?? ~~~~
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
nnn d
d
d
D
ll
l
L
?
?
???
2
1
21
21
1
1
1
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
n
n
n
nnnnnn
n
n
d
ldd
ldldd
ll
l
aaa
aaa
aaa
??
?
?
?
???
?
????
?
?
222
112111
21
21
21
22221
11211
1
1
1
?
?
?
?
?
?
?
??
??
?
???
?
?
?
?
nki
d
dlla
l
dlad
k
k
r
rkrikik
ik
k
r
rkrkkk
,,1,
)(
1
1
1
1
2
?
比较等号两边后,有
1 21 31
21 2 32
31 32 3
6 3 2 1 1
3 5 1 1 1
2 1 6 1 1
d l l
A l d l
l l d
?? ? ? ? ? ? ? ?
? ? ? ? ? ? ? ???
? ? ? ? ? ? ? ?
? ? ? ? ? ? ? ?
? ? ? ? ? ? ? ?
1 11 6da? ? ? 2 1 2 1 1/ 1 / 2l a d? ? ? 31 31 1/ 1 / 3l a d? ? ?
2 2 2 2 1 1 2 1 1 3 / 2d a l d l? ? ?32 32 1 31 21 2( ) / 4 / 13l a d l l d? ? ?
3 3 3 3 1 1 3 1 3 2 2 3 2
26
39d a l d l l d l? ? ? ?
为了提高数值稳定性,可考虑列主元三角分解法,设已完成 A=LU的 k-1步分解计算,
矩阵分解成
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
????
nnnk2n1n
knkk2k1k
n1kk1k12k11k
n2k22221
n1k11211
aa
aa
uu
uuu
uuuu
??
????
??
??
?????
??
??
ll
ll
ll
l
设 )ua(m a xua 1k 1j jktjtkntk1k 1j jkijik ?? ?????? ??? ll令 rk?ri
相当于取 ????? 1k 1j jkijikkk uau l为第 k步分解的主元素,
但要注意方程组的常数项也要相应变换,
定义, 设 ‖ ?‖ 是一种向量范数
AxxAxA
xRxxRx nn 1,0,
s ups up
????
??
称之为由向量范数派生的 矩阵算子范数,
5.3 矩阵范数和条件数
定义 设 ‖ ?‖ 是以 n阶方阵为变量的实值函数,且满足条件,
(1)非负性, ‖ A‖ ?0,且 ‖ A‖ =0当且仅当 A=0
(2)齐次性, ‖ ?A‖ =|? |‖ A‖,??R
(3)三角不等式,‖ A+B‖ ?‖ A‖ +‖ B‖
(4)三角不等式,‖ AB‖ ?‖ A‖‖ B‖
则称 ‖ A‖ 为矩阵 A的 范数,
对应于 3种常见的向量范数,有 3种矩阵范数
?
???
?
n
i
ijnj aA
111
m ax 列和的最大值
?
????
?
n
j
ijni aA
11
m ax
行和的最大值
12 ??A
是 ATA的最大特征值,也称为谱范数
1?
矩阵范数的一些性质:
①
②
③
④
⑤
00,&,0 ???? AAA
RAA ???? ???,
nRBABABA ?????,,
nRBABAAB ????,,
nRxxAAx ????,
定理,若 ? 为 A 的特征值,则 A??
证:
xAx ????
xAx ????
xAxAxx ??????? ??
????
x
xAA
x为 A的特征值
#证毕
定义 5.2,谱半径
rnrA ?? ??? 1m a x)(
易知:
AA ?)(?
5.4 条件数和病态矩阵
定义 5.3:(条件数)
ppp AAAC on d
1)( ??? p?
表示某种范数
设 bAx?, A 引入误差 后A?,解引入误差 x?,则
bxxAA ??? ))(( ??
xAxAAxbxAA ????????? ???? )(
xAAAx ???? ? ??? 1)(
AAA
x
x ??? 1)( ???
AAAAI ?? 111 )( ?????
? ? 111 )()( ??? ??? AAAAIA ??
AAAAI ?? ???? ??? 111 )(
AAAAI ?? 111 )( ??? ???
注意到
BBI ???
?
1
1)( 1
因为:
BDDDBIIBID ???????? ? )(1 )( 1
)1( BDDBD ?????
AA
AA
x
x
?
?? 11
1
1
?
?
?
????
A
A
AA
A
A
AA
AA
AA
?
?
?
?
??
?
?
??
?
?
?
?
?
?
1
1
1
1
1
1
条件数
很小
A
A
AA
?
?? ? 1
条件数表示了对 误差的放大率
同样,类似有
b
b
AA
x
x
bxA
bbxxA
??
??
??
????
??
???
? 1
)(
注,一般判断矩阵是否病态,并不计算 A?1,而由经验
得出。
?行列式很大或很小(如某些行、列近似相关);
?元素间相差大数量级,且无规则;
?主元消去过程中出现小主元;
?特征值相差大数量级。
精确解 为,
1
1??
?
???
?
??x?
例
???
?
???
??
???
?
???
??
97.1
99.1,
98.099.0
99.01 bA ?
计算 cond (A)2 。
???????? ?? 100 00990 0 990 0980 0
A?1 =
解,考察 A 的特征根
??? 0)d e t ( AI?
0 0 0 0 5 0 5 0 4.0
9 8 0 0 5 0 5 0 4.1
2
1 ?????
??
2
1
2)( ?
?Ac o n d 39206 >> 1
?测试病态程度:
给 一个扰动b?
???
?
???
?
?
???
?
?
3
4
10106.0
1097.0b??,其相对误差为
%01.010513.0|||| |||| 4
2
2 ??? ?
b
b? ?? 此时 精确解 为
???
?
???
?
?? 020 3.1
3*x?
???
?
???
?
???? 0203.2
2* xxx ???? ??
2
2
||||
||||
x
x??? 2.0102 > 200%
为对称矩阵