计算机科学计算答案 第二章 矩阵变换和计算(DOC)

第二章 矩阵变换和计算

一、内容提要

本章以矩阵的各种分解变换为主要内容，介绍数值线性代数中的两个基本问题：线性方程组的求解和特征系统的计算，属于算法中的直接法。基本思想为将计算复杂的一般矩阵分解为较容易计算的三角形矩阵. 要求掌握Gauss （列主元）消去法、矩阵的（带列主元的）LU 分解、平方根法、追赶法、条件数与误差分析、QR 分解、Shur 分解、Jordan 分解和奇异值分解.

(一) 矩阵的三角分解及其应用

1．矩阵的三角分解及其应用

考虑一个n 阶线性方程组b Ax =的求解，当系数矩阵具有如下三种特殊形状：对角矩阵D ，下三角矩阵L 和上三角矩阵U ，这时方程的求解将会变得简单.

???????

?

?=n d d d D 2

1, ??????? ??=nn n n l l l l l l L 21222111, ?????

?

?

??=nn n n u u u u u u U 22212111. 对于b Dx =，可得解为i i i d b x /=,n i ,,2,1 =. 对于b Lx =，可得解为1111/l b x =,ii i k k ik

i i l x l

b x /)(1

1∑-=-

=,n i ,,3,2 =.

对于b Ux =，可得解为nn n n l b x /=,ii n

i k k ik

i i l x l

b x /)(1

∑+=-

=,1,,2,1 --=n n i .

虽然对角矩阵的计算最为简单，但是过于特殊，任意非奇异矩阵并不都能对角化，因此较为普适的方法是对矩阵进行三角分解.

1）．Gauss 消去法

只通过一系列的初等行变换将增广矩阵)|(b A 化成上三角矩阵)|(c U ，然后通过回代求与b Ax =同解的上三角方程组c Ux =的解．其中第k 步消元过程中，在第1-k 步得到的矩阵)

1(-k A

的主对角元素)

1(-k kk

a

称为主元.从)

1(-k A

的第j 行减去第k 行的倍数)1()1(--=

k kk

k jk jk a a l （n j k ≤<）称为行乘数（子）.

2）．矩阵A 的LU 分解

对于n 阶方阵A ，如果存在n 阶单位下三角矩阵L 和n 阶上三角矩阵U ，使得LU A =, 则称其为矩阵A 的LU 分解，也称为Doolittle 分解．Gauss 消去法对应的矩阵形式即为LU 分解, 其中L 为所有行乘子组成的单位下三角矩阵, U 为Gauss 消去法结束后得到的上三角矩

阵. 原方程组b Ax =分解为两个三角形方程组??

?==y

Ux b

Ly .

3）．矩阵LU 分解的的存在和唯一性

如果n 阶矩阵A 的各阶顺序主子式),,2,1(n k k =D 均不为零, 则必有单位下三角矩阵L 和上三角矩阵U ，使得LU A =, 而且L 和U 是唯一存在的．

4）．Gauss 列主元消去法

矩阵每一列主对角元以下（含主对角元）的元素中, 绝对值最大的数称为列主元. 为避免小主元作除数、或0作分母，在消元过程中，每一步都按列选主元的Guass 消去法称为Gauss 列主元消去法．由于选取列主元使得每一个行乘子均为模不超过1的数，因此它避免了出现大的行乘子而引起的有效数字的损失.

5）．带列主元的LU 分解

Gauss 列主元消去法对应的矩阵形式即为带列主元的LU 分解,选主元的过程即为矩阵的行置换. 因此, 对任意n 阶矩阵A ，均存在置换矩阵P 、单位下三角矩阵L 和上三角矩阵U ，使得LU PA =．由于选列主元的方式不唯一, 因此置换矩阵P 也是不唯一的. 原方程组

b Ax =两边同时乘以矩阵P 得到Pb PAx =, 再分解为两个三角形方程组?

?

?==y Ux Pb

Ly . 5）．平方根法(对称矩阵的Cholesky 分解)

对任意n 阶对称正定矩阵A ，均存在下三角矩阵L 使T

LL A =，称其为对称正定矩阵

A 的Cholesky 分解. 进一步地, 如果规定L 的对角元为正数，则L 是唯一确定的．原方程

组b Ax =分解为两个三角形方程组?

??==y x L b Ly T .

利用矩阵乘法规则和L 的下三角结构可得

2

1

1

12?

??

? ??-=∑-=j k jk jj jj l a l , jj j k jk ik ij ij l l l a l /1

1???? ??-=∑-=, i=j +1, j +2,…,n , j =1,2,…,n . 计算次序为nn n n l l l l l l l ,,,,,,,,,2322212111 ．由于jj jk a l ≤，k =1,2,…,j ．因此在分解

过程中L 的元素的数量级不会增长，故平方根法通常是数值稳定的，不必选主元．

6）．求解三对角矩阵的追赶法

对于三对角矩阵???

??

??

?

??=---n n n n n b a c b a c b a c b 11122211

A , 它的LU 分解可以得到两个只有两条对角元素非零的三角形矩阵

???

????

?

??=???????? ??=--n n n n u d u d u d u l l l 112

21132,1111

U L . 其中???????=-====-==--n

i c l b u n i u a l b u n i c d i i i i i i i i i ,,3,2,,,3,2,/1,,2,1,111

1

计算次序是n n u l u l u l u →→→→→→→ 33221. 原方程组b Ax =分解为两个三角形方程组??

?==y

Ux b

Ly . 计算公式为

n i y l b y b y i i i i ,,3,2,,

11

1 =-==-,

.

1,,2,1,/)(,/1 --=-==+n n i u x c y x u y x i i i i i n

n n

该计算公式称为求解三对角形方程组的追赶法．当A 严格对角占优时，方程组b Ax =可用追赶法求解, 解存在唯一且数值稳定．

7）．矩阵的条件数

设A 为非奇异矩阵，?为矩阵的算子范数，称1)(cond -=A A A 为矩阵A 的条件数．矩阵的条件数是线性方程组b Ax =, 当A 或b 的元素发生微小变化，引起方程组解的变化的定量描述, 因此是刻画矩阵和方程组性态的量. 条件数越大, 矩阵和方程组越为病态, 反之越小为良态.常用的矩阵条件数为

∞-条件数: ∞

-∞∞=1

)(cond A A A ，

1-条件数: 11

11)(cond -=A A A ，

2-条件数: )

()()(cond min max 2

1

2

2A A A A A

A

A H H λλ==-．

矩阵的条件数具有如下的性质： (1) 1)(cond ≥A ;

(2) )(cond )(cond 1

-=A A ;

(3) )(cond )(cond A A =α，0≠α，R ∈α;

(4) 如果U 为正交矩阵，则1)(cond 2=U ，)(cond )(cond )(cond 222A AU UA ==.

一般情况下，系数矩阵和右端项的扰动对解的影响为

定理2.5 设b Ax =，A 为非奇异矩阵，b 为非零向量且A 和b 均有扰动．若A 的扰动δA 非常小,使得11<-A A δ，则

)(

)

(cond 1)(cond b

δb

A δA A

A

A A x

δx

+-≤δ. 关于近似解的余量与它的相对误差间的关系有

定理2.6 设b Ax =，A 为非奇异矩阵，b 为非零向量，则方程组近似解x ~的事后估计式为

b

x A b A x x x b x A b A ~)

cond(~~)cond(1

-≤-≤-. 其中称x A b ~-为近似解x ~

的余量，简称余量。 8）．矩阵的QR 分解

利用正交变换保条件数的性质, 将满秩矩阵化为主对角元都大于零的上三角矩阵, 保持矩阵条件数不变.

设A 是n 阶可逆实矩阵, 则存在正交阵Q 和对角元都大于零的上三角阵R ，使得

QR A =, 称其为矩阵A 的QR 分解, 并且)(cond )(cond 22R A =.

为实现矩阵一般的QR 分解，我们引入Householder 矩阵T

T

-

=ωωω

ωI ωH 2)(, 其中0,≠∈ωωn R . 该矩阵具有如下性质：

（1） 特征值为：)(2

1))((T T H ωωλω

ωωλ-

= 即，12

1-=-

ωωω

ωT T ，

个

11,,1-n ； （2） )()(ωH ωH =T

, 即H 阵为对称阵；

（3） n I ωH ωH =T

)()(，即H 阵为正交阵；

（4） 如果y x ωH =)(，则22

x y

= (不变长度，镜面反射)；

（5） 设n

n x x x R ∈=T ),,,(21 x 且0x ≠，取12e x x ω-=，则

（6） .00)()(12212e x x x e x x H ωH =????

??

? ??=-= x

提示：Householder 变换并不是直接变换n 阶矩阵A , 而是通过重复变换矩阵的下三角部分

的列向量得到上三角矩阵, 因此, 每次变换的Householder 矩阵

)(,),(),(1-n 21ωH ωH ωH 在逐渐降阶, 然后将它们分别“嵌入”n 阶单位矩阵得到相应的n 阶正交阵1-n 21Q Q Q ,,, , 最后得到正交阵1-n 21Q Q Q Q ,,, =.具体变换过程见例子.

(二) 特殊矩阵的特征系统

特征系统即为矩阵的特征值和特征向量, 本节主要介绍与其计算相关的Schur 分解. 矩阵变换的思想主要为两点: 一是三角矩阵的主对角元素即为其所有特征值, 二是矩阵的特征多项式和特征值在相似变换下是不变的. 因此, 理论上获得矩阵特征值的方法就是通过相似变换将其变为一个三角矩阵.

Schur

定理

:

设

n

n A ?∈C ，则存在酉阵

n

n U ?∈C 使得

H URU A =, 其中n n R ?∈C 为上三角矩阵．

由于实矩阵的特征值可能是复数, 因此通常在复数域中考虑Schur 分解. 复数域中相应的矩阵名称及记号为:

U 的共轭转置: T H U U =, 它在实数域即为转置矩阵. U 为酉阵: 若I UU U U H H ==, 它在实数域即为正交阵.

A 为正规矩阵: 若H H AA A A =.常见的Hermite 阵（A A =H ）

、实对称矩阵（A A =T

）、斜Hermite 阵（A A -=H

）、实反对称矩阵（A A -=T

）、酉阵（I AA A A ==H

H

）和正交矩阵（I AA A A ==T

T

）等均为正规矩阵． Schur 分解的一些特殊情况如下:

● 上三角矩阵R 为正规矩阵当且仅当R 为对角矩阵. ● n 阶方阵A 为正规矩阵当且仅当存在酉阵U 使得H UDU A =,D 为n 阶对角阵. ● n 阶方阵A 为Hermite 阵当且仅当存在酉阵U 使得H UDU A =,D 为n 阶实对角阵. ●

n 阶方阵

A 为酉阵当且仅当存在酉阵U 使得H UDU A =,D 为n 阶对角阵,且对角元的

模均为1.

(三) 矩阵的Jordan 分解介绍

矩阵的每一个特征值有两个重要的指标: 代数重数和几何重数. 一个特征值作为矩阵多项式的根个重数称为代数重数; 它对应的特征子空间的维数称为几何重数. 它们分别刻画了特征值在矩阵特征系统中的代数和几何的性质. 一般有, 代数重数≥几何重数. 当一个特征值的代数重数=几何重数, 称它为半单的; 而当代数重数>几何重数时称它为亏损的.

n 阶方阵

A 可对角化当且仅当它的所有特征值都是半单的, 此时称A 为单纯矩阵; 否

则, A 不可对角化当且仅当它有亏损的特征值, 此时称A 为亏损矩阵.

对于亏损矩阵, 只能将其经过相似变换为一个三角矩阵, 即为其Jordan 标准型. Jordan 标准型是一个块对角矩阵，每一个块称为Jordan 块, 其对角元便为矩阵的特征值．所谓矩阵A 的Jordan 分解即为通过可逆变换矩阵T 化为与之相似的Jordan 标准型J , 使得

1-=TJT A .

1. 关于Jordan 标准型J ．对于特征值i λ, 它的代数重复度就是Jordan 标准型中以i λ为特征值的Jordan 块阶数的和，而其几何重复度（即与i λ相对应的线性无关的特征向量的个数）恰为以i λ为特征值的Jordan 块的个数．J 中以i λ为特征值、阶数为l 的Jordan 块的个数为l l l r r r 211

-+-+,其中l i l I r )(rank A -=λ, n I I r i ==-=)(rank )(rank 00A λ．

2. 关于变换矩阵T 可以通过Jordan 链得到． 将T 按J 的对角线上的Jordan 块相

应地分块为()k T T T T ,,,21 =, 其中T i 为n ×n i 型矩阵．记()

i n i i i i t t t T ,,,21 =, 则

???????+=+==-i n i n i i n i i i i i

i i i i i

112211t t At t t At t At λλλ n

i j C t ∈, k i ,,2,1 =, i n j ,,2,1 =

我们称向量i

n i

i

i t t t ,,,21 为关于特征值i λ的长度为i n 的Jordan 链．显然该Jordan 链的第一个向量就是矩阵A 的关于特征值i λ的特征向量，称其为链首．而链中的第j 个向量则可由等价的方程

()i i

j i

j n i n j ,,3,2,

1 ==--t t I A λ （2-45）

求出．

但是应当注意:

1) Jordan 链的链首i

1t 不仅要求是一个特征向量，而且还要求利用（2-45）可以求出Jordan 链中的其它向量i n i

i t t ,,2 (即不是任何一个特征向量都可作为Jordan 链的链首).

2) 对应于某个特征值i λ 的Jordan 链虽然一定存在，但当与i λ 相对应的线性无关的特征向量的个数大于或等于2时，关于特征值i λ的特征向量中的任何一个有可能都不能作为链首.

因此我们必须从i λ的特征子空间中选取适当的向量作为Jordan 链的链首．

(四) 矩阵的奇异值分解

对于方阵,利用其特征值和特征向量可以刻画矩阵的结构.对非方阵情形,这些方法已经

不适用.而推广的特征值--矩阵的奇异值分解理论能改善这种情况. 利用奇异值和奇异向量不仅可以刻画矩阵的本身结构,而且还可以进一步刻画线代数方程组的解的结构,是构造性的研究线代数问题的有利的工具.

设n

m ?∈C A , Hermite 半正定矩阵A A H 的特征值为021≥≥≥≥n λλλ , 称非负实

数i i λσ=

)(A (n i ,,2,1 =)为矩阵A 的奇异值．

奇异值分解: 设A n

m ?∈C

, 且其秩rank(A )=r , 则存在m 阶、n 阶酉阵U 、V 使得

H V ΣU A ???

? ??=00

0, 其中),,,(diag 21r σσσ =Σ，),,2,1(r i i =σ为矩阵A 的非零奇异值．U 与V 的列向量m u u u ,,,21 和n v v v ,,,21 分别称为矩阵A 的与奇异值i σ对应的左奇异向量和右奇异向量．

利用矩阵的奇异值讨论矩阵的性质:

(1) 矩阵A 的非零奇异值的个数恰为矩阵A 的秩．

(2) },,,{span )(21r u u u A =R , )(A N },,,{span 21n r r v v v ++=，其中n

m R A ?∈?,

},|{)(n m R x y Ax R y ∈?=∈=A R 为由A 的列向量生成的子空间，称为A 的值域或

像空间，即},,,{span )(21n a a a A =R 。)(A N }{0x |=∈=A R x n

称为A 的零空间或核，即}{)(0Ax x A ==N 。

(3) 设021>≥≥≥r σσσ ，则12σ=A , F A =

22

221r σσσ+++ ．

(4) 如果A 为Hermite 矩阵，则A 的奇异值即为A 的特征值的绝对值． (5) 如果A 为n 阶方阵，则 ∏==n

i i

1

)det(σ

A ．

(6) 秩为

r 的

m×n

矩阵A 可以表示为

r 个秩为

1

的矩阵的和

H r r r H

H v u v u v u A σσσ+++= 222111.

(7) n

n C

A ?∈?为正规阵,

λ是A 的特征值, x 是相应于λ的特征向量, 则λ是H A 的特征

值, 相应于λ的特征向量仍为x . (8) n

n C

A ?∈?为正规阵, μλ,是A 的特征值, y x ,是相应的特征向量, 如果μλ≠, 则 x

与y 正交.

2.2典型例题分析

例 1 证明在对矩阵n n j i a ?=)(,A 进行Gauss 消去法的过程中, 主元)

1(,-k k k a (n k ,,2,1 =)

均不为零的充要条件是A 的各阶顺序主子式k D (n k ,,2,1 =) 均不为零.

证明 利用归纳法, 当1=k 时, )0(1,11,11a a D ==, 结论显然成立. 假设结论直到1-k 成立,

则Gauss 消去法可以进行到1-k 步, 即存在1-k 个Gauss 变换11,,-k L L , 使得

???

?

??==-----)1(2,2)1(2,1)

1(1

,111)

1(k k k

k k A 0

A A A L L A

, 其中)

1(1

,1-k A 是对角元为)

1(,-i i

i a (1,,2,1-=k i )的上三角阵, 于是)1(-k A 的k 阶顺序主子阵为

k

k k k k k

k k

a ?---???? ??=)1(,)

1(1

,1)

1(*0A A

. 另一方面, 将)

1(111)(---=k k A L L A 的两端在第k 行k 列处分块有

???

?

?????? ??==----****)()

1(*2*

1

)

1(111k k k k A L 0L A

L L A , 其中*

1L 为k 阶单位下三角阵. 因此A 的k 阶顺序主子式

)

1(,)0(1,1)1()1(*1)det()det()det(---===k k k k k k k k a a D A A L ,

由归纳假设知, 主元0)

1(,≠-k k

k a 当且仅当0≠k D , 即结论对k 成立. 故由归纳法, 0)

1(,≠-k k

k a (n k ,,2,1 =)当且仅当0≠k D (n k ,,2,1 =) .

例2 证明: 若n n j i a ?=)(,A 为可逆矩阵, 则A 可进行LU 分解的充要条件是A 的各阶顺序主子式k D (n k ,,2,1 =) 均不为零.

证明 充分性. 由例1结论知如果k D (n k ,,2,1 =) 均不为零, 则主元0)

1(,≠-k k k a , 于是可

对A 进行Gauss 消去法, 从而得到A 的LU 分解.

必要性. 若存在单位下三角阵L 和上三角阵U 使得LU A =, 则

)

1(,)0(1,1)det()det()det()det(-===n n n a a U U L A ,

由A 可逆知主元0)

1(,≠-k k

k a (n k ,,2,1 =), 再由例1可得A 的各阶顺序主子式

0≠k D (n k ,,2,1 =).