第二章线性代数方程组的直接解法
教学目标:
1.了解线性代数方程组的结构、基本理论以及相关解法的发展历程;
2.掌握高斯消去法的原理和计算步骤,理解顺序消去法能够实现的条件,并在此基础上理解矩阵的三角分解(即LU分解),能应用高斯消去法熟练计算简单的线性代数方程组;
3.在理解高斯消去法的缺点的基础上,掌握有换行步骤的高斯消去法,从而理解和掌握选主元素的高斯消去法,尤其是列主元素消去法的理论和计算步骤,并能灵活的应用于实际中。
教学重点:
1. 高斯消去法的原理和计算步骤;
2. 顺序消去法能够实现的条件;
3. 矩阵的三角分解(即LU分解);
4. 列主元素消去法的理论和计算步骤。
教学难点:
1. 高斯消去法的原理和计算步骤;
2. 矩阵的三角分解(即LU分解);
3. 列主元素消去法的理论和计算步骤。
教学方法:
教具:
引言
在自然科学和工程技术中,许多问题的解决常常归结为线性方程组的求解,有的问题的数学模型中虽不直接表现为线性方程组,但它的数值解法中将问题“离散化”或“线性化”为线性方程组。例如,电学中的网络问题、船体数学放样中建立三次样条函数问题、最小二乘法用于求解实验数据的曲线拟合问题、求解非线性方程组问题、用差分法或有限元法求解常微分方程边值问题及偏微分方程的定解问题,都要导致求解一个或若干个线性方程组的问题。
目前,计算机上解线性方程组的数值方法尽管很多,但归纳起来,大致可以分为两大类:一类是直接法(也称精确解法);另一类是迭代法。例如线性代数中的Cramer法则就是一种直接法,但其对高阶方程组计算量太大,不是一种实用的算法。实用的直接法中具有代表性的算法是高斯(Gauss)消元法,其它算法都是它的变形和应用。
在数值计算历史上,直接法和迭代法交替生辉。一种解法的兴旺与计算机的硬件环境和问题规模是密切相关的。一般说来,对同等规模的线性方程组,直接法对计算机的要求高于迭代法。对于中、低阶(200
n )以及高阶带形的线性方程组,由于直接法的准确性和可靠性高,一般都用直接法求解。对于一般高阶方程组,特别是系数矩阵为大型稀疏矩阵的线性方程组用迭代法有效。
§2.1 基本定理和问题
设具有n 个未知数的n 个方程的线性方程组为
11112211
211222221122n n n n n n nn n n
a x a x a x
b a x a x a x b a x a x a x b +++=??+++=????++
+=? (2.1)
或
1
n
ij j i j a x b ==∑(1,2,
,i n =) (2.1)’
其矩阵形式可以表示为
Ax b = (2.2)
其中
1112121
22
212
n n n n nn a a a a a a A a a a ?? ? ?
= ?
???
12(,,,)T n x x x x = 12(,,
,)T n b b b b =
其增广矩阵为
1,1111212,121
222,112[|]n n n n
n n n n nn a a a a a a a a A A b a a a a +++?? ?
?== ? ? ???
(,1i n i a b +=,1,2,,i n =) (2.3) 则方程组(2.1)、(2.1)’或(2.2)和(2.3)是一一对应的。
若用()r A 表示矩阵A 的秩,则有关于线性方程组(2.2)的解的存在性的基本定理(在高等代数或线性代数中有证明):
定理2.1 (1)方程组(2.2)有解的充分必要条件是()()r A r A =; (2)若()()r A r A k n ==<,则方程组(2.2)具有一族解,其解可表示为
()1n k
i i i x x c x -*
==+∑
其中x *为Ax b =的任一特解,()i x 是齐次方程组0Ax =的解,且(1)(2)(),,,n k x x x -线
性无关,i c 为任意常数。
(3)若()()r A r A n ==,则方程组(2.2)有唯一解。
定义 2.1 如果两个方程组的解相同,则称这两个方程组为同解(等价)方程组。
不难证明:将方程组中任意两个方程交换次序,所得方程组和原方程组为同解方程组;将方程组中任一方程的一个倍数加到另一个方程上,所得方程组和原方程组为同解方程组。
用()i E 表示(2.1)的第i 个方程或(2.3)的第i 行,记(2.3)的第i 行与第j 行互换为()()i j E E ?,而(2.3)的第j 行乘以0α≠加到第i 行记为
()()j i i E E E α+→。这是矩阵的初等变换,相当于[|]A b 左乘一个初等矩阵。同样
的运算符号,我们也理解为方程组(2.1)作相应的变换。经过这些变换后得到的方程组与方程组(2.1)同解(或等价)。
对于线性方程组(2.2)的求解,在理论上并不存在困难。若()r A n =,即A 为非奇异(可逆)矩阵,它的行列式det 0D A =≠,则应用Cramer 法则可求得
i i D
x D =(1,2,,i n =)
其中i D 是用b 代替A 中第i 列而得到的相应的行列式。然而在实际中,当未知数的个数n 比较大时,按Cramer 法则进行计算,其工作量就会大得惊人,因而该方法在实际操作中并不可行。n 阶行列式共有!n 项,每项都有n 个因子,所以计算一个n 阶行列式需要做(1)!n n -?次乘法,我们共需要计算1n +个行列式,要计算出i x ,还要做n 次除法,因此用Cramer 法则求解线性方程组(2.2)就要做
2(1)(1)!(1)!N n n n n n n n =+?-?+=-?+
次乘除法(不计加减法)。如10n =时,359251210N =;当20n =时,
209.707310N ≈?,可见,在实际计算中Cramer 法则几乎没有什么用处。本章的
主要目的就是研究求解线性方程组(2.2)的有效算法。
某一算法的效率可以用下列两个主要的准则来判断:(1)该算法的计算速度如何?即计算中要设计多少次运算?(2)计算所得到解的精度如何?
这两个准则是针对在计算机上求解高阶方程组而提出的。由于线性方程组阶数很高时求解所需要的计算量极其巨大,因而很自然地提出准则(1)。准则(2)的提出,是由于在实际问题当中舍入误差的影响,可能使计算解产生偏离真实解的不可忽视的误差。特别地,在解高阶方程组时涉及大量的运算,舍入误差潜在地积累有可能造成计算解对真实解的严重偏离。后续章节还要详细地研究误差的影响。在研究(2.2)的数值方法之前,先考察一下求解中会遇到的一些困难,
这有助于理解后面将要提出的一些数值计算方法。
§2.2 Gauss 消去法
从这一节开始,我们来讨论线性方程组(2.1)或(2.2)的直接解法。所谓直接法,就是它只包含有限次的四则运算,若假定每一步运算过程中都不产生舍入误差,计算的结果就是方程组的精确解。这种方法中最基本和最简单的就是Gauss 消去法及其变形。
2.2.1 Gauss 消去法的计算过程
设方程组(2.1)或(2.2)的系数矩阵A 非奇异,并记
(1)ij ij a a =,1,2,
,i n =,1,2,,,1j n n =+
(1)(1)[|][|]A b A b =
这样方程组(2.2)又记为(1)(1)A x b =。
要完成Gauss 消去法的第1步须先假定(1)11
0a ≠,令(1)1
1(1)11
i i a l a =(2,
,i n =)
,则运算11()()i i i l E E E -+→(2,
,i n =)将(1)(1)[|]A b 变换为
(1)
(1)(1)(1)1,111
121(2)(2)(2)
2,1(2)(2)22
2(2)(2)(2)
,12
[|]n n
n n n n n nn a a a a a a a A b a a a +++??
? ?
=
? ? ??
?
(2.4) 其中
(2)
(1)(1)11ij
ij
i j
a
a l a =-,(1)1
1(1)11
i i a l a =,2,,i n =,2,,1j n =+
对应(2.4)的方程组是(2)(2)A x b =,它与(2.2)等价,而其第2至第n 个方程中的1x 项已经消去。
一般地,设消去法已进行了1k -步,得到方程组()()k k A x b =。此时对应的增广矩阵为
(1)
(1)(1)(1)1,111121(2)(2)(2)2,122
2()()
()()()
,1()()()
,1[|]n n
n n
k k k k k k n kk
kn k k k n n nk
nn a a a a a a a A b a a a a a a ++++?? ? ?
?
?= ?
? ? ??
?
(2.5) 假设()0k kk
a
≠,令()()k ik
ik k kk
a l a =(1,
,i k n =+),则运算()()
ik k i i l E E E -+→(1,,i k n =+)的结果是方程组(1)(1)k k A x b ++=,对应增广矩阵为(1)(1)[|]k k A b ++,
其中的元素为
()
(1)()(),1,,;1,,1
,1,,;1,,10,1,,k ij k k k ij
ij ik kj a i k j n a a l a i k n j k n i k n +?==+?=-=+=++??=+?
(2.6)
()()
/,1,
,k k ik ik kk l a a i k n ==+
如果依次有()
0k kk
a ≠(1,,1k n =-)
,则可进行第1n -步,得到与( 2.2)等价的方程组()()n n A x b =,其中()n A 是一个上三角阵,且
(1)
(1)(1)(1)1,111
121(2)(2)(2)2,1()()22
2()()
,1[|]n n
n n n n n n n n nn a a a a a a a A b a a +++??
? ?
=
? ? ??
?
这样就完成了消去过程。因为A 非奇异,故有()
0n nn
a ≠。接下来解()()n n A x
b =,因为()n A 是上三角阵,这只要用逐次向后代入的方法即可,这个过程称为回代过程,其计算公式为
()
,1
()
()()
,11
()
,1,2,,1n n n n n nn
n i i i n ij j
j i i i ii a x a a a x x i n n a ++=+?=????-?
?==--??
∑ (2.7) 以上有消去过程和回代过程合起来求出(2.2)的解的过程就称为Gauss 消去法,或称为顺序Gauss 消去法。
从(2.6)可以看出,消去过程的第k 步共含有除法运算n k -次,乘法和减
法运算各()(1)n k n k -+-次,所以消去过程共含有乘除法次数为
321
1
11
5()()(1)326n n k k n n n
n k n k n k --==-+-+-=+-∑∑ 含加减法次数为
31
1
()(1)33n k n n
n k n k -=-+-=-∑
而回代过程含乘除法次数为
(1)2n n +,加减法次数为(1)2
n n -,所以Gauss 消去法总的乘除法次数为332
333n n n n +-≈,加减法次数为32353263
n n n n +-
≈。 当10n =时,用Cramer 法则需要8359251210 3.610≈?乘除法,而用Gauss 消去法仅需430次乘除法运算。
例2.1 用Gauss 消去法解方程组:
(1)12312314012824113x x x ?????? ??? ?= ??? ? ??? ???????
(2)1234123341518311151111631112x x x x -?????? ? ? ?---- ? ? ?= ? ? ? ? ? ?
-????
?? 解:(1)对增广矩阵进行初等变换
133(2)1231412314012801282411300515E E E -+→????
? ??????→ ? ? ? ?--????
得等价的方程组1232332314
28515
x x x x x x ++=??
+=??-=-?
,解得33x =,22x =,11x =。
(2)对增广矩阵进行初等变换
122133144
3
()21
()121
()4
12334153715123341505
22218311155321911116044343111277
70
0444E E E E E E E E E +→-+→-+→-??
??-? ?- ? ?---- ? ???????→ ? ? ? ?
?-?? ?-- ?
?
?
233244344
5
()6
77()()6
11
123341233
4151537370505151522222211
29
11
29
00
0111136367
3579100
00
003633E E E E E E E E E +→+→-+→--????
?
? ? ?-- ? ? ? ??????→??????→ ? ?
?
? ?
? ? ????
?
得等价方程组
123423434412334153715
52221129
1136
910
33x x x x x x x x x x -++=??
?-++=??
?+=??
?=?? 回代得40x =,33x =,22x =,11x =。
2.2.2 消去法的进一步讨论,矩阵的LU 分解
从上面的消去过程可以看出,Gauss 消去步骤能顺序进行的条件是
(1)(2)(1)
11221,1,,
,n n n a a a ---全不为零。设矩阵A 的顺序主子式为i ?,即
11
11
i
i i ii
a a a a ?=
,1,2,
,i n = (2.8)
则有下面的定理:
定理 2.2 ()
i ii a (1,2,
,i k =)全不为零的充分必要条件是A 的顺序主子式
0i ?≠(1,2,
,i k =),其中k n ≤。
证明:设()
0i ii
a ≠(1,2,,i k =),则可以进行消去法的1k -步,每步()m A 由A
逐次实行()()ij j i i l E E E -+→的运算得到,这些运算不改变相应顺序主子式之值,
所以有
(1)(1)(1)11
121(2)(2)(1)(2)()
22
21122()m m m
m mm
m mm
a a a a a a a a a ?=
= 这样便有0m ?≠(1,2,
,m k =)
,必要性得证。 用归纳法证明充分性。1k =时显然成立。设命题对1k -成立。现设
110,
,0,0k k -?≠?≠?≠。由归纳假设有(1)
(1)
111,10,
,0k k k a a ---≠≠,
Gauss 消去法就可以进行1k -步,A 约化为
()
()
()
1112
()22k k k k A A A O A ??= ???
其中()11k A 是对角元为(1)(2)
(1)11221,1
,,
,k k k a a a ---的上三角阵。因为()
k A 是通过消去法由A 逐步得到的,A 的k 阶顺序主子式等于()k A 的k 阶顺序主子式,即
()
(1)
(1)()11
111,1,()*k k k k k k k k k kk
A a a a O
a ---?=
=
由0k ?≠可推出()
,0k k k a ≠。
定理 2.3 对方程组Ax b =,其中A 非奇异,若A 的顺序主子式均不为零,则可以Gauss 消去法求出方程组的解。
定义 2.2 设矩阵()ij n n A a ?=每一行对角元素的绝对值都大于同行其它元素绝对值之和,即1||||n
ii ij j j i a a =≠>∑,1,2,
,i n =,则称A 为严格(行)对角占优矩阵。
定理 2.4 设线性方程组Ax b =的系数矩阵为严格对角占优矩阵,则用顺序
Gauss 消去法求解时()
i ii a (1,2,
,i k =)全不为零。
证明:(略)
下面我们来讨论Gauss 消去过程用矩阵运算表示的形式。 第1步令(1)A A =,作一次运算11()()i i i l E E E -+→,2,,i n =,这相当于A
左乘矩阵
1
11
1i i M l ?? ? ?
?=- ? ? ??
?
,2,,i n =
第1步的全过程相当于(1)(1)(2)(2)1[|][|]L A b A b =,其中
21
11
21111n n n l L M M M l -?? ?- ?== ? ?-??
设1k -步后系数矩阵化为()k A ,其分块形式写成 ()()
(1)()
1112
12
1()22k k k k k k A A L L L A A
O A --??== ???
其中()11k A 为上三角的1k -阶方阵,()
22k A 为1n k -+阶方阵,设其左上角元素()
0k kk a ≠,则下一步的乘数为()()/k k ik ik kk l a a =,1,
,i k n =+。若记
(0,
,0,1,0,
,0)T k e =是第k 个分量为
1的单位向量,记
()1,(0,
,0,,
,)k T k k nk l l l +=,其前k 个分量为零,从而有()
0T k k
e l =。第k 步中系数矩阵的约化可用矩阵运算描述为 (1)
(1)
()
(1)
1112
(1)22k k k k k k A A L A
A
O A ++++??== ???
其中(1)11k A +是上三角的k 阶矩阵,(1)
22k A +是n k -阶方阵,而
()1,1
11
1k T
k k k k nk
L I l e l l +?? ? ? ?
=-=
?- ? ? ? ?-?
?
(2.9) 可以验证
()()()()()()()()k T k T k T k T
k
k k k I l e I l e I l e I l e -+=+- ()()()k T k T
k k I l e l e I =-=
即有1()1()()k T k T
k k k L I l
e I l e --=-=+。这样,经过1n -步得到(1)()n 121n n L L L A A --=,
这里的()n A 是上三角阵,记()n U A =,又记
1(1)(1)1
111()()
()T
n T
n n L L L I l e I l e ----==++
21313212
,11
11
11n n n n l l l l
l l -??
? ?
?= ? ? ???
(2.10) L 是一个对角线元素全为1的下三角阵,这种矩阵称为单位下三角阵。L 的对角线以下元素就是各步消去的乘数。最后我们可以得到A LU =,其中L 是一个单位下三角阵,U 是一个上三角阵。
定义 2.3 将矩阵A 分解为一个下三角矩阵L 和一个上三角矩阵U 的乘积A LU =,称为对矩阵A 的LU 分解或三角分解。当L 是单位下三角矩阵时称为Doolittle (杜里克尔)分解。当U 是单位上三角矩阵时称为Crout (克洛特)分解。
由上面的分析过程知,Gauss 消去法的实质是将系数矩阵分解为一个下三角矩阵和一个上三角矩阵相乘,即将系数矩阵进行LU 分解。
在矩阵A 的LU 分解A LU =中,我们将U 写成U DU =,其中D 是对角矩阵,
U 是单位上三角阵,进一步记L LD =,它是一个下三角阵,这样有
()A LU LDU LD U LU ====
其中L 是一个下三角阵,U 是单位上三角阵,此即A 的Crout 分解。
在矩阵A 的Doolittle 分解A LU =中,我们将上三角阵U 写成DU 的形式,这里的D 为对角阵,U 为单位上三角阵,这样得到A LDU =,其中L 为单位下三角阵,D 为对角阵,U 为单位上三角阵,称其为A 的LDU 分解。
定理2.5 设非奇异矩阵n n A R ?∈,若其顺序主子式(1,
,1)i i n ?=-都不等于
零,则存在唯一的单位下三角阵L 和上三角阵U ,使A LU =。
证明:上面的分析过程已经说明了非奇异矩阵A 可作LU 分解,下面只需证明分解的唯一性。
设A 有两个分解式11A LU =和22A L U =,
其中12,L L 都是单位下三角阵,12,U U 都是上三角阵,则有1122LU L U =。因为A 非奇异,从而12,L L ,12,U U 都可逆。故
在1122
LU L U =两边同时左乘1
1L -和右乘12U -,这样得到111212U U L L --=。因为12U -仍为上三角阵,故112U U -也是上三角阵,同理可得1
1
2L L -是单位下三角阵,结合