漫步线性代数五——三角分解和行交换

来源：互联网发布：大数据在医疗领域应用编辑：程序博客网时间：2024/06/04 19:21

我们继续消去法，看看它对矩阵意味着什么。我们从Ax=b开始：

A x = ⎡ ⎣ ⎢ 24 - 2 1 - 6 7 102 ⎤ ⎦ ⎥ ⎡ ⎣ ⎢ u v w ⎤ ⎦ ⎥ = ⎡ ⎣ ⎢ 5 - 2 9 ⎤ ⎦ ⎥ = b (1)

然后有三个消去步骤，乘数分别为

2，−1，−1：
Step 1. 第二个方程减去第一个方程的2倍；
Step 2. 第三个方程减去第一个方程的-1倍；
Step 3. 第三个方程减去第二个方程的-1倍；
这个结果等价于

Ux=c，其中

U是一个新系数矩阵：

U p p e r t r i a n g u l a r U x = ⎡ ⎣ ⎢ 200 1 - 8 0 1 - 2 1 ⎤ ⎦ ⎥ ⎡ ⎣ ⎢ u v w ⎤ ⎦ ⎥ = ⎡ ⎣ ⎢ 5 - 12 2 ⎤ ⎦ ⎥ = c (2)

矩阵

U是上三角矩阵，主对角线下面的元素都为零。

对于原始向量b，用相同的方法(即将A变为U)推导出新值c，前向消元一共执行了三个行操作：

开始是A和b；
一次应用步骤1,2,3；
最后是U和c。

然后用回代求解Ux=c。这里我们先集中考虑A和U的关系。

第一步的矩阵E，第二步的矩阵F，第三步的矩阵G之前已经介绍了。他们都叫做初等矩阵，很容易看出他们是如何得出的。方程i减去方程j的ℓ 倍，也就是将数−ℓ 放到(i,j)的位置，其他地方跟单位矩阵保持一致，即对角线上为1其他地方为0，这样矩阵乘法就执行了行操作。

三次操作的结果是GFEA=U，注意E首先乘以A，然后是F，最后是G。我们可以将GFE乘在以前得到一个下三角矩阵(0省略掉了)：

G F E = ⎡ ⎣ ⎢ 1111 ⎤ ⎦ ⎥ ⎡ ⎣ ⎢ 1111 ⎤ ⎦ ⎥ ⎡ ⎣ ⎢ 1 - 2 11 ⎤ ⎦ ⎥ = ⎡ ⎣ ⎢ 1 - 2 - 1 111 ⎤ ⎦ ⎥ (3)

但是最重要的问题是反过来的：我们怎样从

U回到

A？如何撤销高斯消元法的步骤？

撤销第一步不难，不用减法而是用第二行加上第一行的二倍。做一次加法和一次减法最后回到了单位矩阵：

⎡ ⎣ ⎢ 120010001 ⎤ ⎦ ⎥ ⎡ ⎣ ⎢ 1 - 2 0 010001 ⎤ ⎦ ⎥ = ⎡ ⎣ ⎢ 100010001 ⎤ ⎦ ⎥ (4)

一个操作取消了另一个操作。用矩阵的语言来描述就是：一个矩阵是另一个矩阵的逆。如果初等矩阵

E在

(i,j)处有值

−ℓ，那么它的逆在同样的位置处为

ℓ。所以

E−1E=I，也就是等式(4)。

我们可以用E−1,F−1,G−1将消元的每步都取逆。最终的问题就是一次撤销整个操作，知道什么矩阵将U回到A。

因为从A到U的过程中步骤三是最后一步，所以矩阵G肯定是第一个取逆的。和正向的顺序正好相反！第二个取逆步骤是F−1，第三个是E−1：

E - 1 F - 1 G - 1 U = A i s L U = A (5)

现在我们将U变为A的矩阵称为L，因为它是下三角矩阵。通过将他们按顺序乘起来就能看出这个性质：

E - 1 F - 1 G - 1 = ⎡ ⎣ ⎢ 1211 ⎤ ⎦ ⎥ ⎡ ⎣ ⎢ 1 - 1 11 ⎤ ⎦ ⎥ ⎡ ⎣ ⎢ 1 1 - 1 1 ⎤ ⎦ ⎥ = ⎡ ⎣ ⎢ 12 - 1 1 - 1 1 ⎤ ⎦ ⎥ = L (6)

对角线下面的元素分别是

ℓ=2,−1,−1。当矩阵相乘时，通常无法直接给出答案。但这里比较特殊直接按从左到右的顺序立马写出来。如果计算机存储每个乘数

ℓij(也就是第

i 行减去第

j行的倍数)，并且在第

i,j的位置得到零，那么这些乘数就记录了消元法。

8、没有行交换的三角分解A=LU。L是下三角矩阵，其中对角线元素为1，下面的元素为ℓij。U是上三角矩阵，它的对角线元素是主元。

例1：

A = ⎡ ⎣ ⎢ 111122123 ⎤ ⎦ ⎥ = ⎡ ⎣ ⎢ 111011001 ⎤ ⎦ ⎥ ⎡ ⎣ ⎢ 100110111 ⎤ ⎦ ⎥ = L U

从

A到

U执行行的减法操作，从

U到

A执行加法操作。

例2：当U是单位矩阵时，L和A 相等

A = ⎡ ⎣ ⎢ 1 ℓ 21 ℓ 31 01 ℓ 32 001 ⎤ ⎦ ⎥

在

A上的消元步骤很容易：

(i)第二行减去第一行的

ℓ21倍，

(ii) 第三行减去第一行的

ℓ31 倍，

(iii)第三行减去第二行的

ℓ32 倍。结果是单位矩阵

U=I。将他们的逆相乘即可得到

⎡ ⎣ ⎢ 1 ℓ 21 11 ⎤ ⎦ ⎥ \times ⎡ ⎣ ⎢ 1 ℓ 31 11 ⎤ ⎦ ⎥ \times ⎡ ⎣ ⎢ 1 1 ℓ 32 1 ⎤ ⎦ ⎥ = ⎡ ⎣ ⎢ 1 ℓ 21 ℓ 31 01 ℓ 32 001 ⎤ ⎦ ⎥

三角分解

A=LU非常重要，所以我在多说一点。以前在线性代数中不讲这些知识，或许是因为它太难(但是相信大家已经掌握了)。如果例2中的

U是任何一个矩阵而不是单位矩阵

U=I，那么我们就是会看到一般的处理过程：

A = L U ⎡ ⎣ ⎢ 1 ℓ 21 ℓ 31 01 ℓ 32 001 ⎤ ⎦ ⎥ ⎡ ⎣ ⎢ r o w 1 o f U r o w 2 o f U r o w 3 o f U ⎤ ⎦ ⎥ = A (7)

它的证明就是利用消元步骤，右边将

A变为

U，左边将

L 变为

I，和例2一样。最终方程两边都等于

U，这些步骤都是双效的，所有(7)是正确的即

A=LU

A=LU非常重要，形式也很美。我们目前用的是3×3 矩阵，但是对于更大的矩阵它也适用。这里我们给出一个例子

例3：

A = ⎡ ⎣ ⎢ ⎢ ⎢ 1 - 1 - 1 2 - 1 - 1 2 - 1 - 1 2 ⎤ ⎦ ⎥ ⎥ ⎥ = ⎡ ⎣ ⎢ ⎢ ⎢ 1 - 1 1 - 1 1 - 1 1 ⎤ ⎦ ⎥ ⎥ ⎥ ⎡ ⎣ ⎢ ⎢ ⎢ 1 - 1 1 - 1 1 - 1 1 ⎤ ⎦ ⎥ ⎥ ⎥

A有三条对角线，

L,U有两条。

一个线性系统=两个三角系统

A=LU有一个很实际的应用，它不仅记录了消元法的步骤；L,U也是求解Ax=b的矩阵。事实上A可以被抛弃了！通过前向消元(利用L)可以将b变成c，然后利用回代(利用U)可以将c变为x：

L c = b t h e n U x = c (8)

用

L乘以第二个方程得

LUx=Lc，其中

Ax=b。每一个三角系统可以很快的求解出来，这样我们的消元法就变为：

分解(从A中求出因子L,U)
求解(从L,U,b中求出解x)

求解子过程满足方程(8)：两个三角方程每个需要n2/2步，所以对于右边任何一个b，我们只需要n2步就能找出答案。这远小于n3/3。

例4：考虑上个例子中的矩阵，右边的b=(1,1,1,1)

A x = b x 1 - x 1 - + - x 2 2 x 2 x 2 - + - x 3 2 x 3 x 3 - + x 4 2 x 4 = = = = 1111 分 解 成 L c = b, U x = c L c = b c 1 - c 1 + - c 2 c 2 + - c 3 c 3 + c 4 = = = = 1111 求 出 c = ⎡ ⎣ ⎢ ⎢ ⎢ 1234 ⎤ ⎦ ⎥ ⎥ ⎥ U x = c x 1 - x 2 x 2 - x 3 x 3 - x 4 x 4 = = = = 1234 求 出 x = ⎡ ⎣ ⎢ ⎢ ⎢ 10974 ⎤ ⎦ ⎥ ⎥ ⎥

对于这些特殊的三角矩阵，操作步骤由

n2变为

2n。可以清楚的看到

Lc=b 通过前向步骤求解出来

c1在

c2 之前)，而

Ux=c 通过回代过程求解出来

x4在

x3 之前)。

注解1：LU关于在对角线上是不对称的：L是1而U 是主元。我们可以让它变得一样，提出一个主元矩阵D即可：

U = ⎡ ⎣ ⎢ ⎢ ⎢ ⎢ ⎢ d 1 d 2 ⋱ d n ⎤ ⎦ ⎥ ⎥ ⎥ ⎥ ⎥ ⎡ ⎣ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ 1 u 12 / d 1 1 u 13 / d 1 u 23 / d 2 ⋱ ⋮ ⋮ ⋮ 1 ⎤ ⎦ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ (9)

上面的例子中主元为1，所以

D=I。但是这只是个特例，一般情况下

LU和

LDU(也写成

LDV)是不相同的。

当看到LDU,LDV时，可以将U,V理解为每行是除以主元得到的。举例说明

A = [1324] = [131] [1 2 - 2] = [131] [1 - 2] [121] = L D U

L,U的对角线上都是1，

D的对角线上是主元1和-2。

注解2：我们已经给出了每个消元步骤的表达式，其中计算必须按顺序进行。对于这一点不全对，有一种克劳特算法计算方法跟他有点不同，在顺序上稍微自由一下。但是最终结果L,D,U不存在自由：

9、如果A=L1D1U1,A=L2D2U2，其中L是下三角矩阵，U是上三角矩阵，D是对角矩阵且对角线上没有零元素，那么L1=L2,D1=D2,U1=U2。LDU分解和LU 分解由A唯一确定。

这个证明需要用到逆矩阵，等学到时再给出详细证明。

行交换和置换矩阵

我们现在必须面对无法避免的问题：主元可能是零。这种情况可能发生在中间的计算过程中，如果a11=0那就是在开始就发生了。如简单的一个例子

[0324] [u v] = [b 1 b 2]

明显这种情况比较麻烦，我们无法利用第一个方程消去系数3。

但是修正也比较明显，交换两个方程，将元素3上移变为主元。对于这个例子矩阵将变成上三角矩阵：

3 u + 4 v 2 v = = b 2 b 1

为了用矩阵形式表示，我们需要置换矩阵

P来产生行交换，通过对单位矩阵

I 进行交换即可得到：

P = [0110] P A = [0110] [0324] = [3042]

对于

b，

P可以产生相同的效果，即交换

b1,b2。新的系统就变为

PAx=Pb，行交换中未知量

u,v没有改变顺序。

置换矩阵P和单位矩阵有相同的行，每行每列都有1。最普通的置换矩阵是单位矩阵P=I(没有发生任何交换)。两个置换矩阵的乘积是另一个置换矩阵。

如果P=I，那么最简单的置换就是只交换两行，其他的置换交换更多。对于大小为n的矩阵存在n!=(n)(n−1)⋯(1)中置换结果，第一行有n种选择，那么第二行有n−1种选择，最后一行只有1种选择。下面我们给出所有的3×3置换矩阵(有3!=(3)(2)(1)=6个矩阵)：

I = ⎡ ⎣ ⎢ 111 ⎤ ⎦ ⎥ P 31 = ⎡ ⎣ ⎢ 111 ⎤ ⎦ ⎥ P 21 = ⎡ ⎣ ⎢ 111 ⎤ ⎦ ⎥ P 32 = ⎡ ⎣ ⎢ 111 ⎤ ⎦ ⎥ P 32 P 21 = ⎡ ⎣ ⎢ 111 ⎤ ⎦ ⎥ P 21 P 32 = ⎡ ⎣ ⎢ 111 ⎤ ⎦ ⎥

如果

n=4，将会有24种置换矩阵，如果

n=2那么只有两种置换矩阵，即

[1001] [0110]

如果我们知道逆和转置(之后将他们定义为A−1和AT)，那么我们将发现一个重要的事实：P−1总是等于PT。

主元位置上是零会产生两种可能：这个问题可能很容易修改，或者比较麻烦。这完全取决于零元素，如果一列下方有非零元素，那么就能执行行交换，将非零元素变成主元然后就能继续消元过程：

A = ⎡ ⎣ ⎢ 00 d a 0 e b c f ⎤ ⎦ ⎥ d = 0 \Rightarrow 没 有 第 一 主 元 a = 0 \Rightarrow 没 有 第 二 主 元 c = 0 \Rightarrow 没 有 第 三 主 元

如果

d=0，那么这个问题将无法解决，矩阵就是奇异的，对

Ax=b不存在唯一解。如果

d不是零，那么通过交换1,3行可以将

d变成主元。然而下一个主元是零，它下面是

a，如果

a不是零，那么用置换矩阵

P23进行行交换：

P 13 = ⎡ ⎣ ⎢ 001010100 ⎤ ⎦ ⎥ P 23 = ⎡ ⎣ ⎢ 100001010 ⎤ ⎦ ⎥ P 23 P 13 A = ⎡ ⎣ ⎢ d 00 e a 0 f b c ⎤ ⎦ ⎥

还有一点是：置换矩阵P23P13可以执行一次行交换：

P 23 P 13 A = ⎡ ⎣ ⎢ 100001010 ⎤ ⎦ ⎥ ⎡ ⎣ ⎢ 001010100 ⎤ ⎦ ⎥ = ⎡ ⎣ ⎢ 010001100 ⎤ ⎦ ⎥ = P

这样的话我们就直接用

P乘以

A，执行了正确的行交换后，那么就可以进行消元了。

消元法：PA=LU

主要的观点是这样的：如果在行交换的帮助下可以完成消元，那么我们可以首先用P来完成这个步骤。矩阵PA不需要行交换。换句话说，PA可以分解成标准的L,U。高斯消元理论可以概括为下面几行：

10、对于奇异的情况，存在一个置换矩阵，它重排矩阵A的行来避免主元位置出现零，然后Ax=b就存在唯一解：

提前进行行排序，PA可以分解成LU。

如果不存在P使得产生所有主元集合，那么消元法将失效。

实际中，当原始主元接近零时，我们也考虑行交换，选择一个更大的主元是为了减小舍入误差。

仔细留意一下L，如果消元过程是先从第二行减去第一行，也就是ℓ21=1，然后2,3行进行交换。或者如果提前进行交换，那么乘因子就变为ℓ31=1。

例5：

A = ⎡ ⎣ ⎢ 112115138 ⎤ ⎦ ⎥ \to ⎡ ⎣ ⎢ 100103136 ⎤ ⎦ ⎥ \to ⎡ ⎣ ⎢ 100130162 ⎤ ⎦ ⎥ = U (10)

然后利用行交换来恢复

LU，但是此时

ℓ31=1,ℓ21=2：

P = ⎡ ⎣ ⎢ 100001010 ⎤ ⎦ ⎥ ⎡ ⎣ ⎢ 121010001 ⎤ ⎦ ⎥ P A = L U (11)

0 0