4.1 因式分解
本节介绍线性代数的一些基本操作，包括行列式、逆和秩，LU分解和QR分解，以及范数等。其中LU分解和QR分解都是使用对角线上方或者下方的元素均为0的三角矩阵来进行计算。使用三角矩阵表示的线性方程组，可以通过向前或者向后置换很容易地得出结果。
4.1.1 行列式、逆和秩
在MATLAB中，用户可以通过以下命令来计算矩阵A的行列式、逆和矩阵的秩。
（1）det(A)：求方阵A的行列式。
（2）rank(A)：求矩阵A的秩，即A中线性无关的行数和列数。
（3）inv(A)：求方阵A的逆矩阵。如果A是奇异矩阵或者近似奇异矩阵，则会给出一个错误信息。
（4）pinv(A)：求矩阵A的伪逆。如果A是m×n的矩阵，则伪逆的尺寸为n×m。对于非奇矩阵A来说，有pinv(A)=inv(A)。
（5）trance(A)：求矩阵A的迹，也就是对角线元素之和。
下面用具体示例对矩阵行列式、逆和秩作简要的说明。
【例4-1】 矩阵的行列式计算示例。
det函数可以用来计算矩阵的行列式。
>>A1=[1 2;3 4] % 创建矩阵A1
A1 =
1 2
3 4
>> A2=[1 2;3,6] % 创建矩阵A2
A2 =
1 2
3 6
>> A3=[1 2;3 4;5 6] % 创建矩阵A3
A3 =
1 2
3 4
5 6
>>det1=det(A1) % 求方阵A1的行列式
det1 =
-2
>> det2=det(A2) %求方阵A2的行列式
det2 =
0
>> det3=det(A3) % 注意非方阵的行列式没有意义
??? Error using ==> det
Matrix must be square.
>> det_1=det(A1′) % 实数矩阵的行列式和它的转置的行列式相同
det_1 =
-2
>> det_2=det(A2′)
det_2 =
0
>> det_3=det(A3′)
??? Error using ==> det
Matrix must be square.
本例使用det函数计算A3的行列式时返回了错误信息，提醒用户A3必须是是方阵才可以调用det函数。
【例4-2】 矩阵的逆计算示例。
本例在上例创建的矩阵基础上进行演示。
>> inv1=inv(A1)
inv1 =
-2.0000 1.0000
1.5000 -0.5000
>>inv2=inv(A2) %A2行列式为0，不存在逆矩阵
Warning: Matrix is singular to working precision.
inv2 =
Inf Inf
Inf Inf
>> inv3=inv(A3) % 非方阵不存在逆矩阵
??? Error using ==> inv
Matrix must be square.
>> detinv1=det(inv(A1)) % A1的逆矩阵的行列式就等于1/det(A1)
detinv1 =
-0.5000
>> 1/det(A1)
ans =
-0.5000
【例4-3】 使用pinv函数计算矩阵的伪逆示例。
>> pinv1=pinv(A1) % A1的逆矩阵和它的伪逆是一样的
pinv1 =
-2.0000 1.0000
1.5000 -0.5000
>> pinv2=pinv(A2)
pinv2 =
0.0200 0.0600
0.0400 0.1200
>> pinv3=pinv(A3)
pinv3 =
-1.3333 -0.3333 0.6667
1.0833 0.3333 -0.4167
本例调用pinv函数计算了矩阵A1、A2、A3的伪逆。因为矩阵A2行列式为0，矩阵A3不是方阵正交矩阵的逆矩阵是它的转置，所以不能求矩阵A2和A3的逆，但是可以求这两个矩阵的伪逆。
【例4-4】 使用rank函数求解矩阵的秩示例。
>> rank1=rank(A1)
rank1 =
2
>> rank2=rank(A2)
rank2 =
1
>> rank3=rank(A3)
rank3 =
2
>> rank_1=rank(A1′)
rank_1 =
2
>> rank_2=rank(A2′)
rank_2 =
1
>> rank_3=rank(A3′)
rank_3 =
2
从本例中可以看出矩阵的秩和它的转置的秩相同。
通过上面的这4个例子，可以总结出以下规律。
（1）只有方阵的行列式才有意义。
（2）只有方阵的逆才有意义，但如果方阵的行列式为0，该方阵则不存在逆矩阵。
（3）如果方阵的逆矩阵存在，它的伪逆和逆相同。
（4）如果方阵的逆矩阵存在，它的逆矩阵的行列式det(A-1)等于1/det(A)。
（5）矩阵的秩和它的转置的秩相同。
（6）实数矩阵的行列式和它的转置矩阵的行列式相同。
4.1.2 Cholesky因式分解
分解法是将原方阵分解成一个上三角形矩阵（或是以不同次序排列的上三角阵）和一个下三角形矩阵，这样的分解法又称为三角分解法，它主要用于简化大矩阵行列式值的计算过程、求逆矩阵和求解联立方程组。需要注意的是：这种分解法所得到的上下三角阵并不是唯一的，可以找到多个不同的上下三角阵对，每对三角阵相乘都会得到原矩阵。
对线性系统的求解，MATLAB依据系数矩阵A的不同，而相应地使用不同的方法。如有可能，MATLAB会先分析矩阵的结构。例如A是对称且正定的，则使用Cholesky分解。
如果没有找到可以替代的方法正交矩阵的逆矩阵是它的转置，则采用高斯消元法和部分主元法，主要是对矩阵进行LU因式分解或LU分解。这种方法就是令A＝LU，其中A是方阵，U是一个上三角矩阵，L是一个带有单位对角线的下三角矩阵。
Cholesky因式分解是把一个对称正定矩阵A分解为一个上三角矩阵R和其转置矩阵的乘积，其对应的表达式为：A=R’R。从理论上说，并不是所有的对称矩阵都可以进行Cholesky因式分解，只有正定矩阵才可以。
Pascal矩阵就是一个有趣的例子。下面以Pascal矩阵为例，说明Cholesky因式分解的使用方法。
【例4-5】 Cholesky因式分解示例。
>> A = pascal(6) %创建Pascal矩阵
A =
1 11 1 1 1
1 23 4 56
1 36 10 1521
1 410 20 3556
1 515 35 70126
1 621 56 126252
矩阵A的元素是二项式系数，每一个元素都是上方和左方两个元素的和。在MATLAB中，进行Cholesky因式分解的是chol函数。矩阵A的Cholesky因式分解可以通过以下命令得到：
>> R = chol(A)
R =
1 11 1 11
0 12 3 45
0 01 3 610
0 00 1 410
0 00 0 15
0 00 0 01
得到的矩阵R的元素同样也是二项式系数。
Cholesky因式分解允许线性方程组A x = b被R’Rx=b代替。在MATLAB环境中，这个线性方程组可以通过以下命令来求解：
>> x=R(R’b)
4.1.3 LU因式分解
LU分解法主要用于简化大矩阵行列式值的计算过程、求逆矩阵和求解联立方程组。需要注意的是：这种分解法得到的上下三角阵对并不是唯一的，可以找到多个不同的上下三角阵对，每对三角阵相乘都会得到原矩阵。
在MATLAB中，求矩阵A的LU分解的调用函数是lu，调用格式如下：
[L,U]=lu(A)
另外，矩阵A的LU分解为线性系统A*x=b提供了以下表达式来快速求解：
x=U(Lb)
【例4-6】 矩阵A的LU分解示例。
>> A=[5 2 0;2 6 2;5 6 7]
A =
5 20
2 62
5 67
>> [L,U]=lu(A) % 分解所得L是带有单位对角线的下三角矩阵，U是上三角矩阵
L =
1.0000 0 0
0.40001.0000 0
1.0000 0.7692 1.0000
U =
5.0000 2.0000 0
0 5.2000 2.0000
0 0 5.4615
>> L*U % 验证结果
ans =
5 20
2 62
5 67
【例4-7】 矩阵A的LU分解实例。
>> A=[1 2 3;4 5 6;7 8 9];
>> [L,U]=lu(A);
>> B=[9 8 7;6 5 4; 3 2 1];
>> x=U(LB)
Warning: Matrix is close to singular or badlyscaled.
Results may be inaccurate. RCOND =1.586033e-017.
x =
-27 -26-17
42 4124
-16 -16-8
>> A*x %验证结果
ans =
9 87
6 54
3 21
4.1.4 QR因式分解
如果A是正交矩阵，那么它满足A’A=1。二维坐标旋转变换矩阵就是一个简单的正交矩阵：
矩阵的正交分解又称做QR分解，是将矩阵分解成一个单位正交矩阵和上三角形矩阵。假设A是m×n的矩阵，那么A就可以分解成：
A=QR
其中Q是一个正交矩阵，R是一个维数和A相同的上三角矩阵，因此Ax=B可以表示为QRx=B或者等同于Rx=QB。这个方程组的系数矩阵是上三角的，因此容易求解。
在MATLAB中，用户可以调用函数qr来求QR因式分解，这个命令可用于分解m×n的矩阵，假设A是m×n的矩阵。qr函数常用调用格式有以下几种。
（1）[Q,R]=qr(A)：求得m×m阶矩阵Q和m×n阶上三角矩阵R。Q和R满足A=QR。
（2）[Q,R,P]= qr(A)：求得矩阵Q，上三角矩阵R和置换矩阵P。R的对角线元素按大小降序排列，且满足AP=QR。
（3）[Q,R]= qr(A,0)：求矩阵A的QR因式分解。如果在m×n的矩阵A中行数小于列数，则给出Q的前n列。
（4）[Q1,R1]=gradelete(Q,R,j)：求去掉矩阵A中第j列之后形成的矩阵的QR因式分解，矩阵Q和R是A的QR因子。
（5）[Q1,R1]=qrinset(Q,R,b,j)：求在矩阵A的第j列前插入列向量b后形成的矩阵的QR因式分解，矩阵Q和R是A的QR因子。如果j=n+1，那么插入的一列放在最后。
【例4-8】 QR分解示例。已知魔方矩阵
对其进行QR分解。
用户只需要调用qr函数就可以实现对A进行QR分解。具体过程如下：
>> A=magic(3)
A =
8 16
3 57
4 92
>> [Q,R]=qr(A) % QR分解
Q =
-0.8480 0.52230.0901
-0.3180 -0.3655 -0.8748
-0.4240 -0.7705 0.4760
R =
-9.4340 -6.2540 -8.1620
0 -8.2394 -0.9655
0 0 -4.6314
>> Q*R % 验证结果
ans =
8.0000 1.0000 6.0000
3.0000 5.0000 7.0000
4.0000 9.0000 2.0000
【例4-9】 利用QR分解求线性方程组的解。
用户可以通过求A的QR分解并计算RQ’*B来求解X。具体过程如下：
>> A=[1 2 2;3 2 2;1 1 2]
A =
1 2 2
3 22
1 12
>> B=[7;9;5]
B =
7
9
5
>> [Q,R]=qr(A)
Q =
-0.3015 0.9239-0.2357
-0.9045 -0.3553-0.2357
-0.3015 0.14210.9428
R =
-3.3166 -2.7136-3.0151
0 1.27921.4213
0 00.9428
>> X=RQ’*B
X =
1.0000
2.0000
1.0000
>> AB
ans =
1.0000
2.0000
1.0000
4.1.5 范数
向量的范数是一个标量，用来衡量向量的长度。需注意不要把向量范数和向量中元素的个数相混淆。在MATLAB中，可以用命令norm得到不同的范数。
norm形式的表达式还有norm(x,-inf)，但它不是求向量的范数，而是求向量x的绝对值的最小值，即min(abs(x))。请读者注意区分。
【例4-10】 向量范数的求解示例。
>> x=[2 4 5]
x =
2 45
>> norm1=norm(x) % 欧几里德范数
norm1 =
6.7082
>> norm2=norm(x,1) % 1-范数
norm2 =
11
>> norm3=norm(x,inf) %¥-范数
norm3 =
5
>> norm4=norm(x,4) % p-范数
norm4 =
5.4727
>> norm5=norm(x,-inf) %向量绝对值的最小值
norm5 =
2
