第3章 线性方程组的直接解法

第三章 线性方程组的数值方法

在自然科学和工程技术中很多问题的解决常常归结为求解线性代数方程组.

?a11x1?a12x2??a1nxn??a21x1?a22x2??a2nxn???ax?ax??ax?n11n22nnn?b1?b2??bn

当它的系数行列式不为零时，由克莱姆法则可以给出方程组的唯一解，但是这一理论上完 善

的结果，在实际计算中可以说没有什么用处。因此如何建立在计算机上可以实现的有效而实用的解法，具有极其重要的意义。这些方法大致可分为两类：一类是直接法，就是经过有限步算术运算,可求得方程组精确解的方法(如果每步计算都是精确进行的话)；另一类是迭代法，就是用某种极限过程去逐步逼近其精确解的方法.

本章将阐述这两类算法中最基本的高斯消元法及其变形、矩阵分解法、雅可比迭代法、高斯 -塞德尔迭代法等.

第一节 高斯消元法

一 回代过程

 设系数矩阵为n阶上三角矩阵的线性方程组

?a11x1?a12x2??a1nxn?a22x2??a2nxn????annxn?如果a11,a22,...ann都

?b1?b2??bn (1)

(1)自下而上可以逐次求出xn,xn?1,...x1为

?bnxn??ann?n?bk??akjxj?j?k?1?,k?n?1,n?2,...1xk?akk?按上述公式求方程组(1) 解的过程称为回代过程.

?2?

11n?n?1?n?n?1?不难看出，解方程组(1)共需2次加法和2次乘法.这恰好是用一个n阶三

2角方阵乘n维向量所需的运算次数。当n较大时,n??n，同时加法运算速度远快于乘法的运

12n2算速度，所以，可用次乘法来近似表示回代过程的运算量.

二 消元过程

设有线性方程组

?a11x1?a12x2??a1nxn?b1??a21x1?a22x2??a2nxn?b2?3?????ax?ax??ax?b?n11n22nnnn ?0??0??aij,ain?1?bi?i?1,2,...n;j?1,2...n?, 则原方程组改写成

为了符号统一，记aij?0??0?0?x1?a12x2??a1?nxn?a11??0??0??0?)?a21x1?a22x2??a2nxn???a?0?x?a?0?x??a?0?x?111n22nnn0??a1?n?1?0??a2n?1??0??ann?1?3?'如果a11量.令

?0?

?0，那么就可以保留其中第一个方程并利用它分别与其余方程消去第一个未知

ai1?,li1?0?a11则以

?0?i?2,3,..,n

?4?

?li1乘第一个方程加到第

i个方程中，就把方程组

?0??a1n?1?3?'化为

?0??0?0?x1?a12x2??a1?nxn?a11??1??1?a22x2??a2x?nn???1??1??ax??ax?n22nnn?a?21n??1?1??a?nn?1?5?

其中

i?2,3,...n,j?2,3,...n?1 ?6?

?0??0?aa1111由方程组(3′)化为(5)的过程中，元素起着特殊的作用，特把元素称为主元素.

?1??1?如果方程组(5)中a22?0, 则以a22为主元素，并利用类似的方法消去第3,4,...n个方程

?1??0??0?aij?aij?li1a1j,中的第二个未知量，即令

ai2i?3,4,..,nli2??1?,a22 ?7?

则以?li2乘以第二个方程加到第i个方程中,于是得到新的方程组

?0??0??0??0??0?x1?a12x2?a13?a11x3???a1nxn?a1n?1??1??1??1??1??ax???ax?axa2222n?12332nn???2??2??2?a33x3?....?a3nxn?a3n?1?????2??2??2??ax??ax?ann?1n33nnn??1??8?

其中

?2??1??1?aij?aij?li2a2j,i?3,...n,j?3,...n?1 ?9?

重复上述过程

n?1步后，我们得到原方程组等价的系数矩阵为三角形方阵的方程组

?0??0??0??0??0?x1?a12x2?a13?a11x3???a1nxn?a1n?1??1??1??1??1??ax???ax?axa2222n?12332nn???2??2??2?a33x3?....?a3nxn?a3n?1?????n?1??n?1??ax?ann?1nnn??10?

其中

把方程组(3)逐步化为方程组(10)的过程称为消元过程.最后，由回代过程可求得原方程组的

解为

?n?1??ann?1?xn??n?1?ann?n??k?1??k?1?axjkn?1??akj?j?k?1,?xk??k?1?akk?aiklik??k?1?akk ?11?

?k??k?1??aij?lika?kjk?1??aij??k?1,2,...n?1?i?k?1,k?2,...n;j?k?1,k?2,...n?1? ?12?

?k?1?k?n?1.n?2,...2,1这种通过消元、再回代的求解方法称为高斯(Gauss)消元法(其特点是始终消去主对角线下方

的元素).

注意到，上标k仅仅用来识别一次消元前后系数矩阵的变化，而aij再使用，所以在计算机存贮中只要用aij?k?1? ?13?

?k?1??k?变为aij后,aij?k?不

即可；另一方面，主元素所在列中主元素下

面的各元素在消元过程中必然是零，而且在后面将要列出的回代过程中也不用它们，所以没有必要通过计算得到它们，从而在消元过程中

例1 用高斯消元法求解方程组

冲掉aij?k?j就可从k?1开始，这样做还可以节约计算时间.

?2x1?8x2?2x3?14??x1?6x2?x3?13?2x?x?2x?53?12

解 用第一个方程消去后两个方程中的x1，得 ?2x1?8x2?2x3?14?2x2?2x3?6???9x2??9 ?再用第二个方程消去第三个方程中的x2，得

?2x1?8x2?2x3?14?2x2?2x3?6???9x3?18 ?最后，经过回代求得方程组的解为

三 高斯消元法的条件与运算量

18??2?96?2x3x2??1214?8x2?2x3x1??52 x3??k?1?从消元过程可以看出，对于n阶线性方程组，只要各步主元素不为零，即akk求得原方程组的解.因此,有下面结论.

?0，经过

n?1步消元，就可以得到一个等价的系数矩阵为上三角形阵的方程组，然后再利用回代过程可

?k?1?a定理1 如果在消元过程中A的主元素不为零，即kk?0(k=1,2,…,n),则可通过高斯

消元法求出Ax?b的解.

?k?1?aA矩阵在什么条件下才能保证kk?0,下面的定理给出了这一条件.

?k?1?aA引理 在高斯消元过程中系数矩阵的主元素不为零，即kk?0(k=1,2,…，n)的充要

条件是矩阵A的各阶顺序主子式不为零，即

D1?a11?0a11...a1kDk?.........?0,k?2,3...nak1...akk证明 首先利用归纳法证明引理的充分性，显然 ，当n?1时，引理的充分性是成立的，现假设引理对n?1时也成立，求证引理对n也成立，由归纳法假设有

?k?1?akk?0, k?1,2...n?1

于是可用高斯消元法将

A化为

(0)a12?(1)a22?且

(0)?a11??A??????a1(n0?)1(1)a2n?1??(n?2)an?1n?1a1(n0)??(1)a2n????(n?2)an?1,n?(n?1)?ann?

?0??0?D1?a11?a11

a11?D20?0?a12?0??1??0??1??a11a22a22

(0)(0)a12?a1(n0?)1a1(n0)??a11??(1)(1)(1)a22?a2n?1a2n???0??1??n?1??...?????a11aaDn??22nn????(n?1)?a?nn?

?n?1?由假设Dk?0,k?1,2,...n,所以有ann?0.

反过来,由上式可知必要性是显然的.

定理2 如果n阶矩阵A的所有顺序主子式均不为零，即Dk高斯消元法求出的解.

下面考虑求解(3)的高斯消元法的运算量.消元过程需要除法

?0,k?1,2,...n,则可通过

?n?1???n?2??...?2?1?1n?n?1?2

次，而需要的乘法和加法的次数都是

1n??n?1???n?1???n?2??...?2?1?n?n2?1?3

加上回代过程的运算次数，共需乘、除法的次数为

1111n?n?1??n?n2?1??n?n?1??n?n2?3n?1?2323

加、减法的次数为

111n?n2?1??n?n?1??n?2n2?3n?5?326 32当n较大时,n??n,消元过程的运算量远大于回代过程，从而，高斯消元法中乘除法的次数

13n与加减法的次数近似为3.

第二节 第二节 高斯主元素消元法

一 问题的提出

由高斯消元法可知，在消元过程中如果出现akk舍入误差的扩散，最后也使得计算结果很不可靠.

?k?1??0的情况，这时消元法将无法进行；另

?k?1?a一方面，即使主元素kk?0，但很小时，用其作除数，会导致其它元素数量级的严重增长和

例1 例1 解方程组

?0.0003x1?3.0000x2?2.0001??1.0000x1?1.0000x2?1.0000

12x1?,x2?33) (它的精确解为

11?8????n?1??n?k?1???n3?nk?122

1312??nnn3,加减法的次数在第二节中我们曾指出，高斯消元法的运算量，乘除法次数为313125n?n?n26，从而，高斯—若当消元法比高斯消元法的运算量乘除法多为31312213121??nnnn?n?n623次，加减法多623次。当n值较大时，高斯消元法比高斯

13n6—若当消元法节省次乘除法和加减法，这个运算量是十分可观的.

n二 逆矩阵

高斯—若当消元法对求 一个矩阵的逆矩阵，或对求解仅常数项不同的很多方程组及矩阵方程是非常有用的. 求矩阵矩阵分块

A的逆矩阵A?1，即求n阶矩阵X，使AX?In，其中In为n阶单位矩阵.将

X??x?1?,x?2?,...x?n?? In??e1,e2,...en?

于是，求解AX?In等价于求解n个方程组

Ax?j??ej,j?1,2,...,n

AX?In的增广矩阵为C??A?In?，如果对C应用

?1???B?B. InA高斯-若当方法化为，则

定理 设A为非奇异矩阵，方程组例2 例2 用高斯—若当消元法求

由线性代数理论，我们有下面结论.

?132?A??254?????365??

的逆矩阵

解

A?1.

?132100?C??A|I3???254010?????365001?? ?132100???0?10?210?????0?3?1?301???102?530???0102?10?????00?13?31??

?9??10?

所以

?1001?32???0102?10?????001?33?1???1?32?A?1??2?10??????33?1??

?11?

为了节省存储单元，可不必将单位矩阵存放起来.作为第一步结果的式(9)中第1列已无用

处，而第4列又相当于逆矩阵所求第1列的中间结果，把它移到第1列不影响简化过程的实质.而且第5、6两列的常数项可取消，它们对简化也无实质影响，所以，最终按原位记法(9)式的结果可存放为

同理,式(10)中的n?则按原位记法为

?2n?阶矩阵将第4列移到第1列，第5列移到第2列，取消第6列，

2???53?2?10?????3?3?1?? ?1?32??2?10??????33?1??

32??1??2?10??????3?3?1??

式(11)的原位结果为

即为我们所要求的逆矩阵A,所以在计算中求逆矩阵的过程可简记为

32??132??1A??254????2?10???????365?????3?3?1??

2???53?1?32??2?10??2?10??????1?3?3?1??????33?1??=A ? ?一般地，逆矩阵的计算公式为

?1a?'kjakk?'akjakk1,j?1,...,k?1,k?1,n?12?

akk?13?

aij?aij?akjaik,''i?1,...,k?1,k?1,nj?1,...,k?1,k?1,ni?1,...,k?1,k?1,nk?1,2,...n

?14??15?

a??aika,'ik'kk式(12)—(15)就是求逆矩阵的基本计算公式，对应于每一个k值就是完成了一个消元过程，在消元过程进行中i和

j变量在规定的范围内进行循环.

我们知道,只要矩阵A的行列式det?0,则A总是可逆的,然而,当主元素为零或绝对值

太小时,按上述方法计算机可能要溢出，因此，在约化的过程中也应采用选主元素的方法,如果互换矩阵的两行,对方程组的解来说,这样的对换对结果没有影响;而对求逆矩阵来说,这样的对换改变了所要求的逆矩阵.事实上,是逆矩阵也作了相应两列的互换,所以,计算逆矩阵也可以通过列主元素消元法,只要记住行的交换,然后在结果中施行相应的列交换即可.

第四节 第四节 矩阵分解 一 矩阵的LU分解

高斯消元过程实际上是对方程组的增广矩阵施行初等行变换，也就相当于用相应的初等矩阵

?1??0??1??0?x?x?bbAA左乘增广矩阵.如果对施行第一次消元后化为, 则存在L1, 使得

?0??1??0??1???bbLAAL11 ,

其中

?1????1l21???1L1???l31???????1???ln1?

?k??k?x?kbA一般地，施行第次消元后化为，则有 ?k?1??k??k?1??k???bbLAALkk ，

其中

重复这一过程，最后得到

?1????????1Lk????1lk?1k?????????1lnk??

Ln?1?L2L1A?0??A?n?1?

Ln?1?L2L1b?0??b?n?1?

?n?1?将上三角矩阵A记为U,则

A?LU

其中

L?L1L2...Ln?1为单位下三角矩阵。

?1?1?1?1??l?1?21?????1??ll?1??n1n2?

这就是说，高斯消元法实质上产生了一个将作用.

A分解为两个三角形矩阵相乘的因式分解,称为

D?0?i?1,2,..n?则A可

A的三角分解或LU分解.于是我们有如下的重要定理,它在解方程组的直接法中起着重要的

定理1(矩阵的LU分解) 设A为n阶矩阵，如果A的顺序主子式i分解为一个单位下三角矩阵L和一个上三角矩阵U的乘积,且这种分解是唯一的. 证明 根据以上高斯消元法的矩阵分析, A?一性，设

LU的存在性已经得到证明,下面证明分解的唯

A?LU?L1U1

L,L1为单位下三角矩阵;U,U1为上三角矩阵，由于A可逆，从而L1与U可逆,故 其中

?1LL?UU1 1

?1上式右边为上三角矩阵,左边为单位下三角矩阵,因此,上式两边都必须等于单位矩阵,于是

L1?L,U1?U例1 例1 求矩阵

的LU分解.

解 由高斯消元法

?111?A??04?1?????2?2?1??

l21?a21a11?0,l31?a31a11?2

且

?111?A?A?0???04?4??A?1?????0?4?1??

?1??4al32?32???1?1?4a22进一步，有，且

?111?A?1???04?1??A?2?????0?4?1??

所以, A?LU,其中

?100??100?L??l2110???010???????2?11?? ?l31l321???U?A?2?

如果把A分解为乘积LU后，求解Ax?b的问题可以看作是相继求解具三角状系数的方程组

的问题,也就是用

?Ly?b??Ux?y

求y，再用

?y1?b1?k?1?lksyj?yk?bk??j?1????k?2,3?,n?

求x,因此,利用LU分解在求解具相同系数矩阵而有不同常数列的方程组时,只要保留L与U的记录就不必要做A的分解及约化的重复工作(这是比较费时的),只要做求解两个三角状系数

方程组的工作就行了,而这两件工作是比较容易和省时的.

二 LU分解的计算公式

定理1告诉我们，如果A的各阶主子式不为0，则存在唯一的LU分解，所以，矩阵分解不一定采用高斯消元法，下面给出一种直接计算方法，设 A?其中

yn?xn??unn??n?/u?x??y?uy??kkksj?kk?j?k?1????k?n?1,n?2?,1?LU ?1?

?1??u11u12u13?u1n??l???1uu?u22232n?21????U???L????????????????????unn??ln1ln2??1????

利用矩阵乘法及矩阵相等则对应元素相等的事实,可以逐一求出L与U的各个元素.首先，从第1行得出U的第1行元素

?2? u1j?a1j,j?1,2,?,n再从第1列算出L的第1列元素

ai1?3?li1?,i?2,3,?,nu11

其次，从第2行算出U的第2行元素

?4? u2j?a2j?l21u1j,j?1,2,?,n再从第2列算出L的第2列元素

li2?ai2?li1u12u22k?1i?1,i?3,?,n?5?

一般地,设已经给出U的第1行到第k行元素为

?1行元素与L的第1列到第k?1列元素，则U的第kukj?akj??lkiuij,j?k,k?1,?,n?6?L的第k列元素为

lik?aik??lijujkj?1k?1

ukk,i?k?1,k?2,?,n?7?

总结上述讨论,可得用直接三角分解求L,U的计算公式

j?1,2,...n?u1j?a1j,??l?ai1,i?2,3,...ni1?u11?k?1?ukj?akj??lkiuij,j?k,k?1,...,n,k?2,3...,ni?1?k?1?aik??lijujk?j?1l?,i?k?1,...,n?ikukk?

uaula由于在电算时，当kj计算好后,kj就不用了,而ik算好后ik也不再使用,因此,计算好的kj与lik可以存放在A的相应位置,例如

?a11a12a13a14??u11u12?aa22a23a24??l21u2221???A???a31a32a33a34??l31l32???aaaa?41424344??l41l42最后,在存放A的数组中,得到分解矩阵L,U的元素. 例2 例2 将矩阵A进行LU分解

?4215??87210??A???4836???1261120??

解

u13u23u33l43u14?u24??u34??u44?

?4?8A???4??12215??47210??8???836??4??61120??12215?7210??836??61120?

?4?2???1??3215??427210??23???836??18??61120??36?4215??4?2300??2??????1236??1???301120???3?4215??4?2300??2??????1221??1???30420???3?1?2L???1??3所以

15?00??36??1120? 215?300??221??01120? 215?300??221??041?

000??4215??0300?100??,?U???0021?210????041?0001??

下面考虑求解矩阵的LU分解及由此求解方程组Ax?b的运算量.第1步要进行除法n?1次，第2步要进行除法n?2次，乘法与加法均为?n?1???n?2?次；一般地.第k步要进行除法n?k次,乘法与加法?k?1??n?k???k?1??n?k?1?次，所以LU分解共

需加减法

k?1次，乘除法

1????????k?1n?k?k?1n?k?1?n?n3337?n2?n26

3k?1?n?k?1??n?k???k?1??n?k?1???n?k???1n3次.而由此求解

Ax?b还需进行加减法

k?12???k?1???n?k???n?nn2?n2?n3

次;乘除法

k?12???k?1???n?k??1??nn

13213121n?nn?n?n263；次.故由LU分解求解线性方程组的运算量乘除法为3加减法3次.与利用高斯消元法的计算量基本相同. 从直接三角分解公式可看出，当

ukk?0时，计算将中断或者当ukk绝对值很小时,按分解公

式计算可能引起舍入误差的累积.因此，对非奇异矩阵

A可采用与列主元素消元法类似的方法，

将直接三角分解法修改为列主元三角分解法.

?1步分解已完成，这时有

?u11u12?lu22?21???A??uk?1k?1uk?1k?lkk?1akk?????ank?ln1ln2?lnk?1u为了避免用小的数kk作除数,第k步分解需引入

设第ku1n?u2n??????uk?1n??akn?????ann??

si?aik??lijujk,j?1k?1i?k,k?1,?,nsi,skk?i?n?9??10?于是有

ukk?sk,令

lik?i?k?1,?,n

sik?maxsi则交换

A的第k行与ik行元素,然后再进行第k步分解计算,于是有

lik?1,i?k?1,...,nL1P1iA?0??A?1?1

下面用矩阵运算来描述列主元素法为

?11? ?12?

?13?

LkPkiA?k?1??A?k?,kk?1,2,...n?1k其中第

Lk的元素满足lik?1,i?k?1,...n,Pki是初等排列矩阵(由交换单位矩阵的第k行与

ik行得到，从而

Ln?1Pn?1i...L2P2iL1P1iA?0??A?n?1??Un?121简记为L~?14?

?U,其中

L?Ln?1Pn?1i...L2P2iL1P1in?12~1令

?15?

Ln?1?Ln?1

Ln?2?Pn?1iLn?2Pn?1in?1~~n?1

n?2n?1Ln?3?Pn?1iPn?2iLn?3Pn?2iPn?1in?1n?2~ ……

?16?

n?1L1?Pn?1iPn?2i...P2iL1P2i...Pn?2iPn?1in?1n?222n?2~

则Lk是单位下三角矩阵，且

~Ln?1Ln?2...L1Pn?1iPn?2i...P1in?1n?2~~~1

~1记

?Ln?1Pn?1iPn?2i...L2P2iL1P1i?L

n?1n?22P?Pn?1iPn?2i...P1in?1n?21?17?

L?Ln?1Ln?2...L1则

?1~~~?18?

?19?

PA?LU其中P为排列矩阵,L为单位下三角矩阵,U为上三角矩阵,总结以上的讨论有 定理2(列主元素的三角分解定理) 如果A为非奇异矩阵，则存在排列矩阵P,使得

PA?LU其中L为单位下三角阵, U为上三角阵.

在列主元素的三角分解中，L的元素存放在数组A的下三角部分；U的元素存放在A的上三

角部分，而P可通过整型数组P(n)表示. 三 平方根法

应用有限元法解结构力学问题时,最后归结为求解线性方程组,系数矩阵大多具有对称正定性，所谓平方根法就是利用对称正定矩阵的三角分解而得到的求解对称正定方程组的一种有效方法.目前在计算机上广泛应用平方根法解此类方程组.

不难证明,在满足定理 1的条件下,有下面结果.

定理3(矩阵的LDU分解) 设A为n阶矩阵，如果

A的各阶顺序主子式，

Di?0,i?1,2,...n，则A可唯一地分解为

A?LDU其中L为单位下三角阵,U为单位上三角阵,D为对角阵.进一步,如果A是对称矩阵,则

?20?

U?LT,即

A?LDLT

定理4(对称正定矩阵的三 角分解) 如果A为n阶对称正定矩阵,则存在一个实的非奇异下三

角阵L，使得

A?LLT当限定L的对角元素为正时,这种分解是唯一的.

由矩阵乘法,可直接获得L的计算公式

2?lkk???akk??lkj?j?1??k?112?21?

?22?i?k?1,...n

lik?(aik??lijlkj)j?1k?1lkkk?1,2,...n按此方法进行的矩阵分解称为平方根法.由于

2akk??lkj,j?1k?23?

k?1,2,...n所以

2lkj?akk?max?akk?1?k?n于是

lkj??max?akk?max?2k,j1?k?n因此,分解过程中元素kj的数量级不会增长,且对角元素

换,不选主元的平方根法就是一个数值稳定的方法. 例3 例3 用平方根法分解对称正定矩阵

l

lkk恒为正数,于是无需进行行的交

?11??4A???14.252.75?????12.753.5??

l11?a11?4?2l21?解

a21?1???0.5l112 a311??0.5l112

l31?2l22?a22?l21?4.25?0.25?2

a?ll2.75?0.5??0.5?l32?323121??1.5l222

22l33?a33?l31?l32?3.5?0.25?2.25?1

T于是A?LL，其中

00??2L???0.500?????0.51.51??

1n?n?1?由于A为对称矩阵,因此,在电算时只要存储A的下三角部分,其需要存储2个元素,

可用一维数组存放,即

?1?A?n?n?1????a11,a21,...,an1,an2,...ann??2? ?1?1A?n?n?1??i?i?1??jaij?的第2矩阵元素存放在?2个位置,L的元素存放在A的相应位

置上.另外,平方根法的运算量是

开平方 n次；

13321n?n?n23次； 乘除法 6137n?n2?n6次. 加减法 6

当n比较大时,平方根法的运算量和存贮量约为高斯消元法的二分之一，因此它是求解对称正定矩阵比较好的方法.

为了避免开方运算,我们可以采用下面的分解式

A?LDLT对于i其中L是单位下三角阵,D是对角阵,由矩阵乘法,可得L与D的计算公式.

?24? ?25??1,2,...,n,有

k?1j?1lik?(aik??lijdjlkj)dki?1j?1k?1,2,...,i?1

di?aii??lij2dj,为了避免重复计算,我们引入

?26?

tij?lijdj于是上述公式可改写成

对于i?27?

k?1,2,...,i?1?1,2,...,n,有

tik?aik??tijlkj,j?1k?1?28??29?

tiklik?,dki?1j?1k?1,2,...,n

di?aii??tijlij,计算出T?30??LD的第i行元素tik,k?1,2,...i?1后,存放在A的第i行相应位置，然后再

ltalt计算L的第i行元素ik仍然存放在A的第i行，即用ik冲掉ik,再用ik冲掉ik,D的对角线

元素存放在A的相应位置上.

对称正定矩阵A按LDL的分解和按LL分解其计算量差不多,但LDL分解不需要开方

计算，它称为改进的平方根法. 四 追赶法

在计算样条函数,解常微分方程边值问题,解热传导方程等都会要求解系数矩阵呈三对角线形的线性方程组,这时

TTT

的LU分解中,矩阵L和U分别取下二对角线和上二对角线形式,设

?a11?a?21A??????a12a22?a23?an?1n?2?an?1n?1ann?1?????an?1n?ann??

?l11??1u12??l?????l2221??L??U???????1un?1n?????ll1?nn?1nn?, ??

由A?LU得计算公式

a11?l11 aii?1?lii?1,i?2,3,...,n aii?lii?1ui?1i?lii,,i?2,3,...,n aii?1?li1uii?1,i?1,2,...,n?1即

l11?a11

ua12?12l11

lii?1?aii?1

lii?aii?lii?1ui?1i

uaii?1ii?1?liii?2,3,...,

n

此时,求解Ax?b等价于解两个二对角线方程组

??Ly?b?Ux?y自上而下解方程组

Ly?b形象地称为“追”.

yb11?l11

ybi?lii?1yi?1i?l,i?2,3,...nii自下而上解方程组

Ux?y称为“赶”.

xn?ynxi?yi?uii?1xi?1,i?n?1,...,2,1习惯,上述求解方法称为“追赶法”. 例4 用追赶法解三对角线方程组

??2x1?x2?1???x1?2x2?x3?0??x2?2x3?x4?0???x3?2x4?1解 由三对角分解公式有

l11?a11?2

u12?a12l11??12 l21?a21??1

l22?a22?l21u12?2?12?32u23?a23l22??23 l32?a32??1

?31??32??33?

而由“追”公式有

l33?a33?l32u23?43 u34?a34l33??34 l43?a43??1

l44?a44?l43u34?54

最后，由“赶”公式得原方程组的解

b11?l112

b?ly1y2?2211?l223 b3?l32y21y3??l334 b?lyy4?4433?1l44 y1?x4?y4?1 x3?y3?u34x4?1 x2?y2?u23x3?1 x1?y1?u12x2?1

追赶法公式实际上就是把高斯消元法用到求解三对角线方程组上去的结果，这时由于A特别简单，因此使得求解的计算公式非常简单，而且计算量仅有5n仅占5n?4次乘除法,3n?3次加减法,

?2个存贮单元,所以可以在小机器上解高阶三对角线形的线性代数方程组.

第五节 向量和矩阵的范数

为了研究线性方程组近似解的误差估计和迭代法的收敛性，我们需要对R（n维向量空间）

n中向量及R（n?n维矩阵空间）中矩阵的“大小”引进某种度量——向量与矩阵范数的概念.向量范数概念是三维欧氏空间中向量长度概念的推广,在数值分析中起着重要作用. 一 向量的范数

n?nn????x?x,x,...,x,y?y,y,...,y?R12n12n 定义1 设，记

TT

称为向量x与

?x,y??xT?y?y?x??xiyii?1Tn?1?y的内积.记非负实数

x2?

称为向量x的欧氏范数.

?x,x????x?ni?1122i?2?

nx,y?R下面结果可在线性代数书中找到，设,k为实数.

① ① 非负性

x2?0且

x2?02当且仅当x?0时成立.

2② ② 齐次性

③ ③ 柯西—许瓦兹不等式

k?x?k?x

?x,y??④ ④ 三角不等式

xy等式当且仅当x与线性相关时成立.

2y2

x?y2?x2?yn2

我们还可以用其它办法来度量R中向量的“大小”,例如,用一个x的函数N的“大小”.在应用中对N的一般定义.

定义2 设

?x?常要求是非负的,齐次的且满足三角不等式,下面我们给出向量范数

n?x?来度量xx① 非负性

n是R上定义的一个实值函数，如果对任意的x,y?R，k?R,满足 x?0x?0?3?x?0且当且仅当时成立.

② 齐次性

③ 三角不等式

k?x?k?x

?4?

?5?

?6?

x?y?x?y则称

是向量x的一个范数(或模). 由定义可推出不等式

xx?yx② 1-范数,也称为绝对和范数

?x?y

实践中最常用的范数有以下几种：

① ∞-范数,也称为最大范数或切彼晓夫范数

??maxxi1?i?n

?7? ?8?

x1??xii?1n③ 2-范数,也称为欧几里德范数

④ p-范数

x2??xi?1??n122i

p

其中关系

xp??xii?1?n?1pp??1,???，可以证明,上述三种范数都是p-范数的特殊情况,而且它们还满足下列

?9?

x??x1?nx?x??x2?nx?1x1?x2?x1n一般地,有下面结论

引理(向量范数的连续性) 设非负函数

?10??11?

?12?

xx为R上的任一向量范数,则

nx是x的分量

x1,x2,...,xn的连续函数.

定理1(向量范数的等价性) 设

xs，

t为R上的任意两个向量范数,则存在常数

nC1,C2?0,使得对一切x?Rn恒有

c1xs?xt?c2xs （13）

证明 实际上,只要能证明一切范数对某一个固定的范数等价,那么任意两种范数都必然等价,因此,可取

记

S?xx?Rn,x??1??，则S是一个有界闭集,由于函数x1x?c2x?t?xs?x??maxxi1?i?n

t在S上连续,所以在S上必

达到其最大值

C2与最小值C1.设x?Rn，且x?0,则

1x?Sx?,从而有

c1?由范数的齐次性,有

?c1xn?xt?c2x

对一切x?R,x?0成立.而对x?0上式显然成立.

注意 定理1不能推广到无穷维空间.由定理还可以看到,对某个向量x来说如果它的某一种范数小(或大),那么它的任 一种范数也不会很大(或很小).

有了范数的概念,我们就可以来讨论收敛的问题.

????为Rx定义3 设

kn中一向量序列,

k??x??Rn,如果

（14）

limx?k??x??0则称x?k??????收敛x义下该向量序 列也收敛.因此,一般按计算的需要采用不同的范数,而且把向量序列

k依范数‖·‖收敛于x.

从上面范数的等价性可以推出,如果在某种范数意义下向量序列收敛,则在任何一种范数意

于向量x记为

?limx?k??x?k??而不强调是在那种范数意义下收敛.

m定理2 设A为m×n阶矩阵,其列向量为线性无关的,如果 ‖·‖是R中范数,则

N?x??Ax,便是R中的一种范数.

nx?Rn

n证明 因为A的列向量线性无关，所以对任意的非零向量x?R有Ax从而

?0

?x?具有非负性;又

N?kx??A?kx??k?Ax?kN?x?

即N?x?具有齐次性；进一步

N?x?y??A?x?y??Ax?Ay?N?x??N?y?

所以, N?x?具有范数的全部基本性质,是一种 范数.

即NN?x??Ax?0

推论 设A为一个n×n阶正定矩阵,则

N?x???xAx?T12, x?R

n是R中的一种范 数.

证明 因为A为正定矩阵,所以有非奇异矩阵L,使得A?nLLT,于是

2

由定理2可知,N二 矩阵的范数

mn

N?x???xTAx???LTx???LTx??LTxT12?x?是R????12n中的一种范数.

m?n一个m×n阶的矩阵也可以看作是m×n维的向量，用R表示m×n阶矩阵的集合,本质上

是和R一样的向量空间，因此可以按向量的办法来定义其上的范数,但是,矩阵还有矩阵间的乘法运算.所以,对于n×n阶的方阵我们定义范数如下 定义4(矩阵的范数) 如果矩阵

① ① 非负性 ② ② 齐次性 ③ 三角不等式 ④相容性

nA?Rn?n的某个非负的实值函数AA?0A?0A?0且

当且仅当

,满足条件

(15)

k?A?k?A (16)

(17)

A?B?A?BA?B?A?Bn?n (18)

则称‖A‖是R的一个矩阵范数.

由向量的2-范数可以得到R中矩阵的一种范数.

由于在大多数与估计有关的问题中,矩阵和向量会同时参与讨论,所以希望引进一种矩阵的范数,它是和向量范数相联系而且和向量范数相容的,即

Ax?A对任意向量x?R及A?Rnn?nx

n?n都成立.为此我们再引进一种矩阵的范数.

n定义5(矩阵的算子范数) 设x?R,A?R,‖·‖是R上的向量范数,记

nA??maxx?0Ax?x??maxAx?x??1 (19)

则

A?是矩阵范数,称为

A的算子范数.进一步它还满足相容性条件

Ax??A?x? （20）

从而

A?也称为从属于向量范数

n?nx?n的矩阵范数.

定理3 设A?R,则

1?i?nj?1① 称为A的行范数. ② 称为A的列范数.

A??max?aijA1?max?aij1?j?ni?1n (21)

(22)

T??AA?表示ATA的最大特征值. 称为A的2-范数,其中max③

A2??max?ATA? (23)

证明 只就①，③给出证明，②同理.

① ① 设x??x1,x2,...,xn??0,不妨设A?0记

Tt?maxxi1?i?n,

n??max?aij1?i?nj?1n

n则

Ax??max?aijxj?t?max?aij1?i?n1?i?nj?1j?1这说明对任意非零x?R，有

n

Axx另一方面，设

0??? ,取向量

?????aijj?10n?,j?1,2,nxj?si?agijn...x0,显然

0x0??x1,x2,...,xn?T,其中

?1n0,且

Ax0的第i0的i0个分量为

?aijxj??aij??j?1j?1n这说明

Ax0??2,即

A????max?aij1?i?nj?1Tn③由于设ATA2??Ax??Ax??xT?ATA?x?0.

nTx?R,对一切,从而AA是半正定的.

A的特征值为

?1??2?...??n?0

?,?,...,?n为ATA的分别对应于?1,?2,...,?n的正交规范的特征向量，则对任一向量再设12x?Rn，有

x??ki?ii?1n,

nTki为组合系数

2A2x另一方面，取

222?x?AA?xi?1?n??1T2xx?kiTi?1?ki?ix??1则

A2Axx2?1T?ATA??1???12T?1?1x2

2A2?maxx?02??1??max?ATA?.

故

?,1还是比较容易的，由定理3看出,计算一个矩阵的而矩阵的2-范数

上不方便,但由于它有许多好的性质.所以,它在理论上是有用的.

AAA2在计算

定义6设A?Rn?n的特征值为

?i?i?1,2,...n?称

（24）

??A??max?i1?i?n为A的谱半径.

定理4(特征值上界) 设A?Rn?n,则

??A??A (25)

即A的谱半径不超过A的任何一种算子范数.

证明 设?是A任一特征值,x为相应的特征向量，则

Ax??x，由定义得

即

??x??x?Ax?A??A，所以??A??A.

进一步，有下面结果 定理5 设

x

A?为与某种向量范数

???相容的矩阵范数，则

??A??infA? （26）

证明 对任意矩阵A，总存在相似变换化A为上三角形,即存在可逆矩阵P，使

PAP?1???U

??diag??1,?2,...,?n?，各?i为A的特征值，U是主对角元为零的上三角矩阵

其中,

?u?ijn?n.

做矩阵

D?diag?1,??1,...,?1?n?,于是

其中矩阵B的元素为

??0

D???U?D?1???B?C

令

?0,bij???uij?j?ij?i????uij?j?i?1?,j?iQ?DP,则有

QAQ?1?C

矩阵A和C相似,所以有相同的特征值.

今规定向量范数令

??如下：

y?Qx，则有

A??maxx?0Ax?x?22?maxx?0QAxQx22

?maxy?0Cyy?C2??2?B2

因为?为对角矩阵，所以

?2????????A?其次,存在常数k>0,使

B2???K从而有

?A????A???K

?????0K, 上式就化为 今对给定的，取

??A??A????A???从而便证明了（26）.

定理6 如果A?Rn?n

为对称矩阵，则

A2???A?事实上,

2 (27)

2A2???ATA????A2?????A??由于

A2???ATA?2

,所以

A2也常记为

ASP并称为谱范数.

定理7 如果

B?1则I?B为非奇异矩阵,且

?I?B??1?11?B第六节 误差分析

(28)

其中‖·‖是指矩阵的算子范数.

前几节讨论了求解线性代数方程组的直接法.给出系数矩阵A和自由项b,求未知向量x.实践中,A和b往往是实验观测数据或是计算所得结果.因此我们处理的线性方程组成了

Ax?b实际上变

?A?A???x?x????b??b? （1）

?b，其中

?121314?A??131415?????141516???A或?b与?x的关系怎样,是人们十分关心的问题.

例1 例1 解方程组Ax现用绝对精确的计算(即不带任何舍入误差的计算)求解,可以看出

?x1?72b1?240b2?180b3??x??240b?900b?720b?123??2??x2?180b1?720b2?600b3??

此时,我们发现对于两组不同的自由项.

它的差只有?b?b?b????,?,???,而所得解x与x之差却是

T?x?~x?x???492?,1860?,?1500??T~T~Tb??b1,b2,b3?,b??b1??,b2??,b3???

~

换句话说,两组不同的右端其分量之差不过是,可是解的差却高达之1860倍.

对于这样的方程组,不管用什么样的数值方法，我们总很难(甚至不可能)算出合理的(与真正精确解相差不大的)解,像这样的方程组或矩阵A就叫做病态的.

定义1 如果矩阵A或自由项b的微小变化,引起方程组Ax?b解的巨大变化,则称此方 程组为病态方程组,矩阵A称为病态矩阵,否则称方程组为良态方程组,A为良态矩阵.

应该注意,矩阵的病态性质是矩阵本身的特性.下面我们研究方程组(1),希望能找出刻划矩阵病态性质的量,为了简单,先假设?A?0,讨论自由项对x的影响，再假设?b?0,讨论系数矩阵与解x的关系.

如果方程组(1)中系数矩阵A是精确的，自由项b有误差?b，相应的解为x???x?，则

A?x?x???b??b (2)

利用关系式

Ax?b

A?x??b即?x?A?1?b

两边取范数,有

?xx??bAx??bAA?1b

?xx故

?A?1x1A?b?AA?1Ax?b?AA?1b?b

A?1??bb??xx?AA?1?bb (3)

A是微小误差?A,相应的解为x?x?

?A??A???x?x???b (4)

这种情况比前一种情况复杂,A??A可能是奇异矩阵.为了简单,我们要求?A能够保证使A??A非奇异,因为

A??A?A?I?A?1?A? A?1?A?1I?A?1?AA??A如果b是精确的,由第六节定理7，当

时,

非奇异,从而

非奇异,且

?A??A??1?A?1?1?I?A??A?1?1?A??1?A?11?A?1?A

下面的讨论使用了更强的条件,即

A,于是

?x?即

A?1??A?x1?A?A?1?A?1??A?x1?A?1?A

?xx另一方面,(4)利用两边取范数,有

A??1?A??AA1?A?1?A??AA (5)

Ax?b有

?x??A?1?A?x?x??

?x?A?1?1?A??A?A?A?x??xA由此可以看出,量

（6）

A?1?A越小由A或b的相对误差引起的解的相对误差就越小,量

越大，解的相对误差就可能越大.所以量

化的灵敏程度,刻划了方程组的病态程度.于是引进下述定义

定义2 设A为非奇异矩阵,称数

A?1?AA?1?A实际上刻划了解对原始数据变

cond?A??A?1cond??A??A?1?A (7)

为矩阵A的条件数.

矩阵的条件数与范数有关,常使用的条件数有

??A?1??8??9?

cond2?A??A2?A当A为对称矩阵时

2?max?ATA???min?ATA??1cond2?A???2其中

不难证明,条件数具有下列性质: ① ① 对任何非奇异矩阵 A，都有cond (10)

?1,?n分别为A的绝对值最大和最小的特征值.

② ②

③ ③ 如果A为正交矩阵，则

?A??1; 11)

设A为非奇异矩阵,C≠0为常数，则cond?CA??cond?A? (12)

cond2?AB??cond2?AB??cond2?B? (14)

cond2?A??1 (13)

其中B为非奇异矩阵.

Ax?b的解对初始数据误差的灵敏度，其

值越大，这种灵敏度越高，即对很小的初始误差?b或?A，解x的相对误差就有可能很大，从

据以上的讨论可以看出，cond(A)反映线性方程组

而大大破坏了解的精确度.当cond(A)接近于1时，矩阵是良态的，否则是病态的.当A是正交矩阵时，cond2A?1,故正交矩阵是最稳定的一类矩阵。对病态系数矩阵解方程组或求逆矩阵，舍入误差的影响十分明显，因此应特别注意. 例2 例2 求矩阵A的条件数，其中

?? 解 因为

?121314?A??131415?????141516?? ?121314?A??131415?????141516??

所以

13?1?i?3j?112

?240180??72A?1???240900?720?????180?720600?? A?max?aij?3于是

A?1??1860,从而

cond??A??A所以A是病态的.

例3 例3 研究方程组

?1?A???2015

系数矩阵为

?10?4x1?x2?1?x1?x2?2 ??10?41?A???11??

1?51?5?1?,?2?,22有特征从而

如果用列主元素消元法,则交换方程次序

?1cond2?A???2.62?2

x1?x2?2???1?10x1?x2?1

有系数矩阵

1??1A1???1??101?

?1??1??2?2??1?10,??1?10;所以 2及相应的特征值1?1??1cond2?A1???1??1.022?2

cond2?A1?比cond2?A设x是方程组

~?小，并且cond?A?接近于1，这说明选主元的必要性.

21Ax?b的近似解，令

~ e?x?x (15)

称为解的误差向量,而令

r?b?A~x?Ae (16)

称为解的剩余向量或残余向量.

下面讨论剩余向量r和误差向量e的关系.由(3)式,有

rer1??cond?A?cond?A?bxbr当cond (17)

?A?接近于1时,

exb与ex相差不大.然而,当矩阵A的条件数cond?A?比较

rbex比较小时,

rb也可能很大.所以对病

大时，尽管比较小,但却可能很大；反之，当

态的方程组,不能单从r或e来刻划解的近似程度.

~对给定的近似解x，用(16)的左等式计算r,然后再用右等式计算e，如果r和e都是精确

~值, 则x?e是方程组Ax?b的真解.实践中由于舍入误差的影响,只能得到e的近似e,从而

~~~xx?e仍是x的近似值,但它的近似程度显然比

要高.

~~如果用高斯消元法或其它方法得到解方程

Ae?m?Ax?b的近似解x,记为x?1??~x,则m步的剩余向量

r?m??b?Ax?m? (18)

?r?m? (19)

?m??m?xe得近似解和的修正向量

x?m?1??x?m??e?m? (20)

e?m?公式(18),(1 9),(20)给出改善近似解x的迭代过程,当过程停止，这就是所谓迭代校正法.

~e?m?或

x?m?小于误差限?时，迭代

本文来源：https://www.bwwdw.com/article/6iz6.html