26、MATLAB 线性代数：矩阵、向量空间与多项式求解

最新推荐文章于 2025-11-23 13:45:13 发布

postgres8guard

最新推荐文章于 2025-11-23 13:45:13 发布

阅读量11

点赞数

CC 4.0 BY-SA版权

分类专栏： MATLAB从入门到精通文章标签： MATLAB 线性代数矩阵分解

本文链接：https://blog.youkuaiyun.com/postgres8guard/article/details/155147465

MATLAB从入门到精通专栏收录该内容

29 篇文章 ¥499.90

订阅专栏¥69.90

会员秒杀 ¥9.9 重磅福利

超级会员免费看

MATLAB 线性代数：矩阵、向量空间与多项式求解

1. 标准矩阵与特殊矩阵

在数值模拟和分析问题中，标准矩阵和特殊矩阵具有重要作用。例如， ones() 、 eye() 、 zeros() 、 rand() 、 randn() 等标准矩阵常用于信号处理、数据分析和大规模计算中的内存分配。同时， gallery 函数可以生成一些特殊矩阵，如 Householder 矩阵、对称 4x4 Hankel 矩阵和烟雾矩阵。

>> B = gallery('house',rand(3,1))  % Householder matrix
B =
       2.1404
      0.28584
       0.7572
>> B = gallery('ris',4)  % Returns a symmetric 4-by-4 Hankel matrix
B =
    0.1429    0.2000    0.3333    1.0000
    0.2000    0.3333    1.0000   -1.0000
    0.3333    1.0000   -1.0000   -0.3333
    1.0000   -1.0000   -0.3333   -0.2000
>> gallery('smoke', 2) % Smoke matrix - complex, with "smoke ring" pseudo-spectrum.
ans =
  -1.0000 + 0.0000i   1.0000
   1.0000             1.0000

2. 向量空间

在信号处理、数值分析和计算机模拟模型构建中，向量空间非常重要。创建向量空间有几种简单的方法：
- 已知步长 ：如果已知步长 $\Delta w$，可以使用 W = w1: Δw: w2 来创建向量空间。例如， w = 1:0.1:13 定义了从 1 到 13，步长为 0.1 的向量空间。
- 已知点数 ：如果已知空间起始和结束边界之间的点数 $N$，可以使用 linspace() 函数。例如， w = linspace(1, 13, N) 创建一个包含 $N$ 个等间距元素的线性空间。如果未指定 $N$，默认值为 100。

>> fs = 3e4;          % Sampling frequency
>> f1 = 100; f2 = 200; f3 = 300; f4 = 400; f5 = 500; f6 = 600;
>> t = 0:1/fs:5;
>> x = sin(2*pi*t*f1)+sin(2*pi*t*f2)+sin(2*pi*t*f3)+sin(2*pi*t*f4)+sin(2*pi*t*f5);
>> [m, n] = size(x);  % Check the size of the created vector space
>> sound(x, fs)       % Play a created sound & hear from sound cards

此外，MATLAB 中还有 logspace() 函数，用于创建对数缩放的向量空间。例如， x = logspace(1,13,130) 创建一个包含 130 个对数间距元素的向量空间，范围从 1 到 13。如果未指定点数，默认值为 50。

3. 用向量表示多项式

在 MATLAB 中，多项式通过向量表示，向量元素为多项式系数，按降序排列。例如，一个五次多项式 $12x^5 + 13x^4 - 15x^2 + 17x - 13$ 可以表示为向量 f = [12 13 0 -15 17 -13] 。如果多项式缺少某个系数，需要用 0 填充。

以下是几种计算多项式根的方法：
- roots() 函数 ：使用 MATLAB 的 roots() 函数计算多项式的根。

>> f = [12  13 0 -15 17 -13];
>> x_sols = roots(f)
x_sols =
  -1.2403 + 0.9412i
  -1.2403 - 0.9412i
   0.7941
   0.3015 + 0.6869i
   0.3015 - 0.6869i

solve() 函数 ：使用 solve() 函数求解多项式的符号解，然后使用 double() 函数将解转换为数值数据。

>> syms x
>> Sol = solve(12*x^5+13*x^4-15*x^2+17*x-13)
Sol =
 root(z^5 + (13*z^4)/12 - (5*z^2)/4 + (17*z)/12 - 13/12, z, 1)
 root(z^5 + (13*z^4)/12 - (5*z^2)/4 + (17*z)/12 - 13/12, z, 2)
 root(z^5 + (13*z^4)/12 - (5*z^2)/4 + (17*z)/12 - 13/12, z, 3)
 root(z^5 + (13*z^4)/12 - (5*z^2)/4 + (17*z)/12 - 13/12, z, 4)
 root(z^5 + (13*z^4)/12 - (5*z^2)/4 + (17*z)/12 - 13/12, z, 5)
>> double(Sol)
ans =
   7.9411e-01 + 0.0000e+00i
  -1.2403e+00 + 9.4121e-01i
  -1.2403e+00 - 9.4121e-01i
   3.0155e-01 + 6.8690e-01i
   3.0155e-01 - 6.8690e-01i

zero() 函数 ：使用 Control Toolbox 的 zero() 函数计算多项式的根。

>> f_tf = tf(f,1, 's')
Transfer function:
12 s^5 + 13 s^4 - 15 s^2 + 17 s - 13
>> x_sols = zero(f_tf)
x_sols =
  -1.2403 + 0.9412i
  -1.2403 - 0.9412i
   0.7941
   0.3015 + 0.6869i
   0.3015 - 0.6869i

4. 使用 Simulink 模型求解多项式

可以使用 Simulink 建模来求解多项式。具体步骤如下：
1. 使用 MATLAB Fcn 块、Constant 块输入多项式系数，Display 块查看计算结果。
2. 在 MATLAB Function 块中嵌入以下命令：

function y = fcn(u1, u2, u3, u4, u5, u6)
y = roots([u1, u2, u3, u4, u5, u6]);

调整求解器类型为固定步长类型：通过 Simulation ➤ Model Configuration Parameters ➤ Solver Selection ➤ Fixed Step Solver 进行设置。
更改输出变量 y 的大小：点击 Model Explorer ➤ [Model Hierarchy] ➤ Polynomial_Solver.slx ➤ MATLAB Function ➤ y Output ➤ Size ，将大小设置为 5。
运行模型，查看结果。

5. 特征值和特征向量

特征值和特征向量不仅在线性代数中，而且在许多工程问题中都有广泛应用，如振动、模态分析、控制应用、机器人等。

定义：对于方阵 $A$，特征值 $\lambda$ 和特征向量 $\nu$ 满足 $A\nu = \lambda\nu$。

计算特征值和特征向量的方法：
- 手动计算 ：通过求解特征多项式 $\det(A - \lambda I) = 0$ 得到特征值，再将特征值代入方程求解特征向量。
- eig() 函数 ：使用 MATLAB 的 eig() 函数计算特征值和特征向量。

>> A = ([2.3, 3.4, 5; 3, 2.4, -1.5; 2, -0.4, -7.2])
A =
    2.3000    3.4000    5.0000
    3.0000    2.4000   -1.5000
    2.0000   -0.4000   -7.2000
>> [v, lambda] = eig(A)
v =
    0.4725   -0.7649    0.6510
   -0.2479   -0.6366   -0.7276
   -0.8458   -0.0983    0.2164
lambda =
   -8.4345         0         0
         0    5.7726         0
         0         0    0.1619
>> A*v - v*lambda % Verify: eigen-vectors and eigen-values;
ans =
  1.0e-014 *
    0.9326    0.4441    0.1346
   -0.2220   -0.0888    0.1402
   -0.3553   -0.0222    0.0520

eig() 函数有多种语法形式，还有 eigs() 函数用于计算最大特征值和特征向量。

以下是 eig() 和 eigs() 函数的常见语法：
| 函数 | 语法 | 说明 |
| ---- | ---- | ---- |
| eig() | d = eig(A) | 计算方阵 $A$ 的特征值 |
| eig() | [V,D] = eig(A) | 计算方阵 $A$ 的特征值和特征向量 |
| eig() | [V,D] = eig(A,'nobalance') | 不进行平衡处理计算特征值和特征向量 |
| eig() | [V,D] = eig(A,B) | 计算广义特征值和特征向量 |
| eigs() | d = eigs(A) | 计算方阵 $A$ 的最大特征值 |
| eigs() | [V,D] = eigs(A) | 计算方阵 $A$ 的最大特征值和特征向量 |

6. 矩阵分解

矩阵分解在线性代数和工程问题求解中有广泛应用，如求解线性方程组、线性最小二乘、非线性优化等。常见的矩阵分解类型包括 QR、LU、Cholesky 等。

6.1 QR 分解

QR 分解也称为正交三角分解，将矩阵 $A$ 分解为正交矩阵 $Q$ 和上三角矩阵 $R$ 的乘积，即 $A = QR$，其中 $Q^TQ = I$。

在 MATLAB 中，使用 qr() 函数进行 QR 分解，有多种语法形式：
- [Q,R] = qr(A) ：产生上三角矩阵 $R$ 和单位矩阵 $Q$。
- [Q,R] = qr(A,0) ：产生经济规模的分解。
- [Q,R,E] = qr(A) ：产生 $Q$、$R$ 和置换矩阵 $E$，使得 $A E = Q R$。

graph LR
    A[矩阵 A] -->|qr(A)| Q[正交矩阵 Q]
    A -->|qr(A)| R[上三角矩阵 R]
    A -->|qr(A,0)| Q1[经济规模 Q]
    A -->|qr(A,0)| R1[经济规模 R]
    A -->|qr(A,E)| Q2[带置换 Q]
    A -->|qr(A,E)| R2[带置换 R]
    A -->|qr(A,E)| E[置换矩阵 E]

示例：计算 5x5 矩阵的 QR 分解

>> format short
>> A = randn(5,5)
A =
    1.2280   -0.5261    0.6770    0.9856    0.8005
    1.8766    1.3087    0.1504    2.0377   -0.0944
   -0.8958    0.6341    0.5682   -0.9155   -1.2224
    0.2791    1.0490   -0.1173    0.7488    1.2225
   -0.7456    0.3522    0.1206   -0.4902    0.1131
>> [Q, R] = qr(A)
Q =
   -0.4829    0.4184   -0.6628   -0.3755   -0.1068
   -0.7380   -0.5121   -0.0038    0.2917    0.3286
    0.3523   -0.4417   -0.7048    0.3635   -0.2279
   -0.1098   -0.5418    0.1633   -0.5900   -0.5654
    0.2932   -0.2720   -0.1932   -0.5419    0.7135
R =
   -2.5428   -0.5003   -0.1895   -2.5283   -0.8486
         0   -1.8346   -0.0140   -0.4991    0.2301
         0         0   -0.8922    0.2012    0.5091
         0         0         0   -0.2847   -1.5550
         0         0         0         0   -0.4484

6.2 LU 分解

LU 分解也称为高斯消元法的改进形式，将矩阵 $A$ 分解为下三角矩阵 $L$ 和上三角矩阵 $U$ 的乘积，即 $A = LU$。

在 MATLAB 中，使用 lu() 函数进行 LU 分解，有多种语法形式：
- Y = lu(A) ：对于稀疏矩阵 $A$，产生矩阵 $Y$，$Y$ 仅包含 $L$。
- [L,U] = lu(A) ：产生 $U$ 和 $L$。
- [L,U,P] = lu(A) ：产生 $U$ 和 $L$，并带有单位对角和置换矩阵 $P$。

graph LR
    A[矩阵 A] -->|lu(A)| Y[稀疏矩阵 L]
    A -->|lu(A)| L[下三角矩阵 L]
    A -->|lu(A)| U[上三角矩阵 U]
    A -->|lu(A,P)| L1[带置换 L]
    A -->|lu(A,P)| U1[带置换 U]
    A -->|lu(A,P)| P[置换矩阵 P]

示例：计算 3x3 Pascal 矩阵的 LU 分解

% LU_decomposition.m
A = input('Enter rectangular matrix: ');
if  issparse(A)
    Y = lu(A)            %#ok
    [L,U,P,Q] = lu(A)    %#ok
    disp('  oops more  ')
    [L,U,P,Q,R] = lu(A)  %#ok
    [L, U, P, Q, R] = lu(A,'vector') %#ok
else
    [L,U] = lu(A)        %#ok
    [L,U,P] = lu(A)      %#ok
    % Check evaluation results:
    ERROR = P*A - L*U    %#ok
    [L,U,P] = lu(A, 'vector')        %#ok
end

运行脚本，输入 pascal(3) 作为矩阵：

Enter rectangular matrix: pascal(3)
L =
    1.0000         0         0
    1.0000    0.5000    1.0000
    1.0000    1.0000         0
U =
    1.0000    1.0000    1.0000
         0    2.0000    5.0000
         0         0   -0.5000

6.3 Cholesky 分解

Cholesky 分解对于蒙特卡罗模拟和卡尔曼滤波器设计特别重要，仅适用于对称正定矩阵。将矩阵 $A$ 分解为下三角矩阵 $L$ 与其转置 $L^T$ 的乘积，即 $A = LL^T$。

在 MATLAB 中，使用 chol() 函数进行 Cholesky 分解，有多种语法形式：
- R = chol(A) ：产生上三角矩阵 $R$，满足 $R’ R = A$。
- L = chol(A,'lower') ：产生下三角矩阵 $L$，满足 $L L’ = A$。

graph LR
    A[矩阵 A] -->|chol(A)| R[上三角矩阵 R]
    A -->|chol(A,lower)| L[下三角矩阵 L]
    A -->|chol(A,P)| R1[带标志 R]
    A -->|chol(A,P)| P[标志 P]
    A -->|chol(A,S)| R2[带置换 R]
    A -->|chol(A,S)| S[置换矩阵 S]

示例：计算 4x4 Pascal 矩阵的 Cholesky 分解

% Chol_decomposition.m
clearvars; clc
disp('Note your matrix must be square & positive definite!!!')
disp('NB: Positive means all determinants must be positive.')
disp('You can enter as matrix elements ')
disp('or define your matrix 1st, ')
disp('and then just enter your matrix name')
disp('      ')
A = input('Enter a given Matrix: ');
[rows, cols] = size(A);
for k = 1:rows
    % Determinants are computed
    Det_A(k) = det(A(1:k, 1:k));
end
if rows == cols
  if true(Det_A > 0) == 1
     if issparse(A)
            [R,p,S] = chol(A) %#ok
            [R,p,s] = chol(A,'vector');
            [L,p,s] = chol(A,'lower','vector');
else
            R = chol(A) %#ok % Upper triangular matrix R: R'*R=A
            L = chol(A,'lower') %#ok % Lower triangular matrix R.
            [R,p] = chol(A);
% Verify:
Error_up  = A - R'*R;
Error_low = A - L*L';
disp('Error is with upper triangular matrix: ')
disp(Error_up)
disp('Error is with lower triangular matrix:')
disp(Error_low)
       end
else
warndlg('Sorry your matrix is not positive and definite!')
warndlg('Try again!!!')
run('Chol_decomposition')
    end
end

运行脚本，输入 pascal(4) 作为矩阵：

Enter a given Matrix: pascal(4)
R =
     1     1     1     1
     0     1     2     3
     0     0     1     3
     0     0     0     1
L =
     1     0     0     0
     1     1     0     0
     1     2     1     0
     1     3     3     1
Error is with upper triangular matrix:
     0     0     0     0
     0     0     0     0
     0     0     0     0
     0     0     0     0
Error is with lower triangular matrix:
     0     0     0     0
     0     0     0     0
     0     0     0     0
     0     0     0     0

综上所述，本文介绍了 MATLAB 中矩阵操作、向量空间创建、多项式求解、特征值和特征向量计算以及矩阵分解的相关知识和方法，并通过具体示例进行了说明。这些方法在数值模拟、工程计算等领域具有广泛的应用。

MATLAB 线性代数：矩阵、向量空间与多项式求解

7. 不同矩阵分解的对比与应用场景

不同的矩阵分解方法在不同的场景下具有各自的优势。下面通过表格对比 QR、LU 和 Cholesky 分解的特点和应用场景：
| 分解类型 | 特点 | 应用场景 |
| ---- | ---- | ---- |
| QR 分解 | 将矩阵分解为正交矩阵和上三角矩阵的乘积，稳定性好 | 求解线性方程组、最小二乘问题、特征值计算的预处理等 |
| LU 分解 | 将矩阵分解为下三角矩阵和上三角矩阵的乘积，计算效率高 | 求解线性方程组、矩阵求逆等 |
| Cholesky 分解 | 仅适用于对称正定矩阵，分解为下三角矩阵与其转置的乘积 | 蒙特卡罗模拟、卡尔曼滤波器设计等 |

在实际应用中，需要根据矩阵的性质和具体问题选择合适的分解方法。例如，对于一般矩阵的线性方程组求解，QR 分解和 LU 分解都可以使用，但如果矩阵是对称正定的，Cholesky 分解会更加高效。

8. 向量空间创建方法的选择

在创建向量空间时， linspace() 和 logspace() 函数各有特点。 linspace() 用于创建线性间距的向量空间，适用于需要等间距采样的场景，如信号处理中的时间序列生成。而 logspace() 用于创建对数间距的向量空间，适用于需要对数尺度采样的场景，如频率响应分析。

以下是一个对比示例：

% 线性间距向量空间
linear_space = linspace(1, 10, 10);
% 对数间距向量空间
log_space = logspace(0, 1, 10);

disp('线性间距向量空间:');
disp(linear_space);
disp('对数间距向量空间:');
disp(log_space);

运行上述代码，可以看到线性间距向量空间的元素是等间距的，而对数间距向量空间的元素是按对数尺度分布的。

9. 多项式求解方法的性能比较

在求解多项式的根时， roots() 、 solve() 和 zero() 函数各有优缺点。 roots() 函数是 MATLAB 中最常用的求解多项式根的函数，它直接计算多项式的数值根，计算速度快。 solve() 函数用于求解多项式的符号解，适用于需要精确解或需要进行符号运算的场景，但计算速度相对较慢。 zero() 函数是 Control Toolbox 中的函数，主要用于控制系统中的传递函数求解，在处理控制系统相关问题时更加方便。

以下是一个性能比较示例：

% 定义多项式系数
f = [12  13 0 -15 17 -13];

% 使用 roots() 函数求解
tic;
x_sols_roots = roots(f);
time_roots = toc;

% 使用 solve() 函数求解
syms x;
tic;
Sol = solve(12*x^5+13*x^4-15*x^2+17*x-13);
x_sols_solve = double(Sol);
time_solve = toc;

disp(['roots() 函数求解时间: ', num2str(time_roots)]);
disp(['solve() 函数求解时间: ', num2str(time_solve)]);

运行上述代码，可以看到 roots() 函数的求解时间通常比 solve() 函数短。

10. 特征值和特征向量计算的注意事项

在计算特征值和特征向量时，需要注意以下几点：
- 矩阵的性质 ：不同性质的矩阵（如对称矩阵、正定矩阵等）在计算特征值和特征向量时可能有更高效的方法。例如，对称矩阵的特征值都是实数，且可以使用更稳定的算法进行计算。
- 数值稳定性 ：在使用 eig() 和 eigs() 函数时，可能会遇到数值不稳定的问题，特别是对于病态矩阵。可以使用 eig(A,'nobalance') 语法来避免平衡处理，提高数值稳定性。
- 计算效率 ： eigs() 函数用于计算矩阵的部分特征值和特征向量，在处理大规模矩阵时可以提高计算效率。

以下是一个数值稳定性示例：

% 病态矩阵
A = [1 1e-10; 0 1];

% 不进行平衡处理计算特征值和特征向量
[V,D] = eig(A,'nobalance');

disp('特征向量:');
disp(V);
disp('特征值:');
disp(D);

运行上述代码，可以看到不进行平衡处理可以得到更准确的特征值和特征向量。

11. 总结与展望

本文详细介绍了 MATLAB 中线性代数的相关知识，包括标准矩阵与特殊矩阵、向量空间的创建、多项式的求解、特征值和特征向量的计算以及矩阵分解等内容。通过具体的代码示例和解释，展示了这些方法的使用和应用场景。

在实际应用中，需要根据具体问题选择合适的方法和函数，同时注意数值稳定性和计算效率。未来，可以进一步探索 MATLAB 在更复杂的线性代数问题中的应用，如大规模矩阵的处理、高维向量空间的分析等。此外，结合其他工具箱和技术，如深度学习工具箱、并行计算技术等，可以进一步提高线性代数计算的效率和精度。

希望本文对读者在 MATLAB 线性代数的学习和应用中有所帮助。通过不断实践和探索，读者可以更好地掌握这些知识和技能，解决实际问题。