25、MATLAB 线性代数操作与矩阵运算详解

MATLAB 线性代数操作与矩阵运算详解

1. 嵌入 MATLAB 函数块计算行列式和求解线性方程

在 MATLAB 中,用于计算矩阵行列式、矩阵逆或线性系统解的函数/命令,都可以通过 MATLAB 函数块嵌入到 Simulink 中。下面以计算任意大小矩阵行列式的 det() 函数和求解线性方程的 linsolve() 函数为例,将它们嵌入到 Simulink 模型中。

步骤如下
1. 定义矩阵 A1 B1

A1=[16 2 -3 13 ; -5 11 10 -8; 9 7 -6 12; -4 14 15 1 ];
B1 = [3; 2; 4; 5];
  1. 构建 Simulink 模型:使用三个 Content 块、两个 MATLAB 函数块和两个 Display 块构建模型。
  2. 编辑 MATLAB 函数块:
    • 双击 MATLAB 函数块打开 MATLAB 编辑器窗口。
    • 对于上侧的 MATLAB 函数块(有一个输入),输入以下脚本:
function y = fcn(u)
y = det(u);
- 对于下侧的 MATLAB 函数块(有两个输入 `A1` 和 `B1`),输入以下脚本:
function y = fcn(A1, B1)
y = linsolve(A1, B1);
  1. 保存并执行模型:编辑完块的代码后保存,然后执行最终的模型。上侧的 Display 块显示行列式,下侧的 Display 块显示给定系统的解。

验证结果 :将 Simulink 模型的计算结果与 MATLAB 计算结果进行比较:

>> det(A1)
ans =
      -18812
>> linsolve(A1, B1)
ans =
    -3.614714012332536e-03
     2.740803742292154e-01
     6.272591962577077e-02
     2.075271103550925e-01

结果显示,行列式计算和线性 MATLAB 求解器的计算结果与 Simulink 模型的结果在 13 位小数内匹配。

2. 线性方程求解函数的准确性

为了找出计算解时哪个函数/工具(方法)更准确,我们使用 MATLAB 的 magic() randi() 函数生成一个 13×13 的矩阵 A 和一个 13×1 的矩阵 B ,并使用 norm() 函数计算给定线性系统及其计算解的范数。以下是完整的解决方案脚本 LA_ex2.m

%% Given 13-by-13 system of linear equations
clearvars; clc
A = magic(13);
B = randi([-169,169], 13,1); % Elements of B vary within [-169, 169]
%% 1-Way: inv()  or pinv()     %% INVERSE matrix method
x1a = inv(A)*B;
Err_INV = norm(A*x1a-B)/norm(B)     %#ok: ERROR checking
x1b = pinv(A)*B;
Err_PINV = norm(A*x1b-B)/norm(B)     %#ok: ERROR checking
%% 2-Way: \     %% backslash
x2  = A\B;
Err_BACKSLASH = norm(A*x2-B)/norm(B)    %#ok: ERROR checking
%% 3-Way: mldivide()          %% Left divide function
x3 = mldivide(A, B);
Err_MLDIVIDE = norm(A*x3-B)/norm(B)   %#ok: ERROR checking
%% 4-Way: Using linsolve()
x4 = linsolve(A, B);
Err_LINSOLVE = norm(A*x4-B)/norm(B)   %#ok: ERROR checking
%% 5-Way: Using lsqr()
x5 = lsqr(A, B);
Err_LSQR = norm(A*x5-B)/norm(B)   %#ok: ERROR checking
%% 6-Way: Using lu()
[L, U, P] = lu(A);
y = L\(P*B);
x6 = U\y;
Err_LU = norm(A*x6-B)/norm(B)   %#ok: ERROR checking
%% 7 - Way: Using rref()
MA = [A, B];
xyz = rref(MA);
x7= xyz(:,end);
Err_RREF = norm(A*x7-B)/norm(B)   %#ok: ERROR checking
%%  8 - Way: Using svd()
[U, S, V]= svd(A);
x8 = V*inv(S)*U'*B;
Err_SVD = norm(A*x8-B)/norm(B)   %#ok: ERROR checking
%% 9 - Way: Using chol()
[U, L] = chol(A);
x9 = U\(U'\B);
Err_CHOL = norm(A*x9-B)/norm(B)   %#ok: ERROR checking
%% 10 - Way: Using qr()
[Q, R] = qr(A);
x10  = R\Q.'*B;
Err_QR = norm(A*x10-B)/norm(B)   %#ok: ERROR checking
%% 11 - Way: Using decomposition()
x11 = decomposition(A)\B;
Err_DECOMPOSITION = norm(A*x11-B)/norm(B) %#ok: ERROR checking
%% 12 - Way: Using bicg()
x12 = bicg(A,B);
Err_BICG = norm(A*x12-B)/norm(B)   %#ok: ERROR checking
%% 13-Way: solve()            %% SOLVE() symbolic math method
x = sym('x', [1, 13]); x=x';
Eqn = A*(x); Eqn = Eqn - B;
Solution = solve(Eqn); SOLs = struct2array(Solution);
SOLs = double(SOLs); x13 = SOLs';
Err_SOLVE = norm(A*x13-B)/norm(B)   %#ok: ERROR checking

计算结果

Err_INV =
   5.8087e-16
Err_PINV =
   3.7982e-15
Err_BACKSLASH =
   3.0569e-16
Err_MLDIVIDE =
   3.0569e-16
Err_LINSOLVE =
   3.0569e-16
lsqr converged at iteration 7 to a solution with relative residual 3.5e-07.
Err_LSQR =
   3.4959e-07
Err_LU =
   3.0569e-16
Err_RREF =
   1.1576e-05
Err_SVD =
   3.7982e-15
Err_CHOL =
    2.2400
Err_QR =
   7.0615e-16
Err_DECOMPOSITION =
   3.0569e-16
bicg stopped at iteration 13 without converging to the desired tolerance 
1e-06
because the maximum number of iterations was reached.
The iterate returned (number 13) has relative residual 9.6e-06.
Err_BICG =
   9.5856e-06
Err_SOLVE =
   1.5109e-16

从计算的误差可以看出, RREF() BICG() LSQR() 函数在计算给定系统的解时误差在 $10^{-5}…10^{-7}$ 范围内,而其他方法的误差在 $10^{-15}…10^{-16}$ 范围内。

3. 线性方程求解函数的效率

为了演示哪种方法在计算时间上更高效,我们使用随机整数生成函数 randi() 生成两个大矩阵:一个 1000×1000 的矩阵 A 和一个 1000×1 的矩阵 B 。同时,使用 tic toc 函数记录每个计算方法的耗时。以下是完整的解决方案脚本 LA_ex3.m

clearvars
A=randi([-100,100],1000); B=randi([-100, 100], 1000, 1);
%% 1) inv() or pinv()
tic; Ai = inv(A); xyz1=Ai*B; T_inv=toc
%% 2) bacslash oerator: \
clearvars
A=randi([-100,100],1000); B=randi([-100, 100], 1000, 1);
tic; xyz2 = A\B; T_backslash = toc
%% 3) mldivide()
clearvars
A=randi([-100,100],1000); B=randi([-100, 100], 1000, 1);
tic; xyz3= mldivide(A, B); T_mld = toc
%% 4) linsolve()
clearvars
A=randi([-100,100],1000); B=randi([-100, 100], 1000, 1);
tic; xyz4 = linsolve(A, B); T_linsolve = toc
%% 5) lsqr()
clearvars
A=randi([-100,100],1000); B=randi([-100, 100], 1000, 1);
tic; xyz5 = lsqr(A, B); T_lsqr = toc
%% 6) lu()
clearvars
A=randi([-100,100],1000); B=randi([-100, 100], 1000, 1);
tic; [L, U, P]=lu(A); y=L\(P*B); xys6=U\y; T_lu=toc
%% 7) rref()
clearvars
A=randi([-100,100],1000); B=randi([-100, 100], 1000, 1);
tic; MA = [A, B];xyz7 = rref(MA); XYZ7=xyz7(:, end); T_rref=toc
%% 8) svd()
clearvars
A=randi([-100,100],1000); B=randi([-100, 100], 1000, 1);
tic; [U S V] = svd(A); xyz8 = V*inv(S)*U'*B; T_svd=toc
%% 9) chol()
clearvars
A=randi([-100,100],1000); B=randi([-100, 100], 1000, 1);
tic; [U L]= chol(A); xyz9 = U\(U'\B); T_chol=toc
%% 10) qr()
clearvars
A=randi([-100,100],1000); B=randi([-100, 100], 1000, 1);
tic; [Q R] = qr(A); xyz10 = R\Q.'*B ; T_qr=toc
%% 11) decomposition()
clearvars
A=randi([-100,100],1000); B=randi([-100, 100], 1000, 1);
tic; xyz11 = decomposition(A)\B; T_decom = toc
%% 12) bicg()  Gradient methods
clearvars
A=randi([-100,100], 1000); B=randi([-100, 100], 1000, 1);
tic; xyz12 = bicg(A, B); T_bicg=toc
%% 13) solve()
A=randi([-100,100],100); B=randi([-100, 100], 100, 1);
tic;
x = sym('x', [1, 100]); x=x';
Eqn = A*(x); Eqn = Eqn - B;
Solution = solve(Eqn); SOLs = struct2array(Solution);
SOLs = double(SOLs); x13 = SOLs';
T_solve=toc

计算结果

T_inv =
    0.0390
T_backslash =
    0.0173
T_mld =
    0.0171
T_linsolve =
    0.0171
lsqr stopped at iteration 20 without converging to the desired tolerance 
1e-06
because the maximum number of iterations was reached.
The iterate returned (number 20) has relative residual 0.24.
T_lsqr =
    0.0236
T_lu =
    0.0235
T_rref =
   10.1406
T_svd =
    0.4263
T_chol =
    0.0330
T_qr =
    0.1045
T_decom =
    0.0459
bicg stopped at iteration 20 without converging to the desired tolerance 1e-06
because the maximum number of iterations was reached.
The iterate returned (number 0) has relative residual 1.
T_bicg =
    0.0195
T_solve =
   14.6306

从这些计算结果可以看出, linsolve() mldivide \ (反斜杠运算符)(高斯消元法)在所有测试方法中速度最快。最慢且计算成本最高的是符号数学工具箱的 solve() 运算符,即使系统规模小了 10 倍也是如此。值得注意的是, rref() 函数(简化行阶梯形法)是仅次于 solve() 运算符的次慢方法。

4. 通过改变 [b] 的值求解线性方程

本练习由两部分组成:
1. 求解给定线性系统中的未知数 a b c

0.072*a - c = -12
0.12*b - c = -9
a + b = 50
  1. 求解给定系统中的未知数 a b c ,其中第三个方程的值在 50 到 250 之间变化。

以下是解决方案脚本 LA_Ex4.m

% PART I.
% The given system is written from the Ax=B as [A]*[abc]=[B]
A=[.072, 0, -1; 0, .12, -1; 1 1 0];
B=[-12, -9, 50];
abc1=A\B'                %#ok  % BACKSLASH \
abc2 = linsolve(A,B')    %#ok  % LINSOLVE()
abc3 = inv(A)*B'         %#ok  % INV
% SOLVE() in symbolic MATH
syms a b c
abc4=solve(0.072*a-c+12, 0.12*b-c+9, a+b-50);
abc4=double([abc4.a; abc4.b; abc4.c]) %#ok  % SOLVE()
%% Part II. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% BACKSLASH \ ; LINSOLVE(); INV
tic
Bk=50:250;
a=zeros(numel(Bk),1); b=zeros(numel(Bk),1); c=zeros(numel(Bk),1);
A=[.072, 0, -1; 0, .12, -1; 1 1 0];
for ii=1:numel(Bk)
    B=[-12; -9; Bk(ii)];
    abc=A\B;
    a(ii)=abc(1,:);
    b(ii)=abc(2,:);
    c(ii)=abc(3,:);
end
Time1=toc;
fprintf('Computation time with BACKSLASH: %3.3f  \n', Time1)
clearvars
tic
Bk=50:250;
a=zeros(numel(Bk),1); b=zeros(numel(Bk),1); c=zeros(numel(Bk),1);
A=[.072, 0, -1; 0, .12, -1; 1 1 0];
for ii=1:numel(Bk)
    B=[-12; -9; Bk(ii)];
    abc=linsolve(A,B);
    a(ii)=abc(1,:);
    b(ii)=abc(2,:);
    c(ii)=abc(3,:);
end
Time2=toc;
fprintf('Computation time with LINSOLVE:  %3.3f  \n', Time2)
clearvars
tic
Bk=50:250;
a=zeros(numel(Bk),1); b=zeros(numel(Bk),1); c=zeros(numel(Bk),1);
A=[.072, 0, -1; 0, .12, -1; 1 1 0];
for ii=1:numel(Bk)
    B=[-12; -9; Bk(ii)];
    abc=inv(A)*B;
    a(ii)=abc(1,:); b(ii)=abc(2,:); c(ii)=abc(3,:);
end
Time3=toc;
fprintf('Computation time with INV:       %3.3f  \n', Time3)
%% SOLVE() from symbolic math
clearvars
tic;
Bk=50:250;
a1=zeros(numel(Bk),1);b1=zeros(numel(Bk),1); c1=zeros(numel(Bk),1);
syms a b c
for ii=1:numel(Bk)
    abc=solve(0.072*a-c+12,0.12*b-c+9,a+b-Bk(ii));
    a1(ii)=double(abc.a);
    b1(ii)=double(abc.b);
    c1(ii)=double(abc.c);
end
Time4=toc;
fprintf('Computation time with SOLVE:     %3.3f  \n', Time4)

计算结果

abc1 =
   15.6250
   34.3750
   13.1250
abc2 =
   15.6250
   34.3750
   13.1250
abc3 =
   15.6250
   34.3750
   13.1250
abc4 =
   15.6250
   34.3750
   13.1250
Computation time with BACKSLASH: 0.002
Computation time with LINSOLVE:  0.002
Computation time with INV:       0.002
Computation time with SOLVE:    22.066

从使用四种方法计算具有三个变量和 201 种可能情况的给定线性系统的解所花费时间来看,最不高效的求解线性方程的方法是使用符号数学工具箱的 solve() 函数。反斜杠运算符 \ linsolve() inv() 方法的性能相似,它们比 solve() 函数高效且快 11033 倍以上。

5. 矩阵运算

矩阵运算涵盖了矩阵、向量和特征向量的一般数学运算和计算。以下是一些常见矩阵运算及其 MATLAB 命令语法:
| 运算名称 | MATLAB 第一种方式 | MATLAB 第二种方式 |
| — | — | — |
| 矩阵乘法 | A B | mtimes(A,B) |
| 数组乘法 | A.
B | times(A,B) |
| 矩阵右除 | A/B | mrdivide(A,B) |
| 数组右除 | A./B | rdivide(A,B) |
| 矩阵左除 | A\B | mldivide(A,B) |
| 数组左除 | A.\B | ldivide(A,B) |
| 矩阵幂 | A^B | mpower(A,B) |
| 数组幂 | A.^B | power(A,B) |
| 复转置 | A’ | ctranspose(A) |
| 矩阵转置 | A.’ | transpose(A) |
| 二进制加法 | A+B | plus(A,B) |
| 一元加 | +A | uplus(A) |
| 二进制减法 | A-B | minus(A,B) |
| 一元减 | -A | uminus(A) |
| 行列式 | det(A) | det(A) |
| 旋转 90 度 | rot90(A) | rot90(A) |
| 数组复制和拼接 n 次 | repmat(A, n) | repmat(A, n) |
| 矩阵左右翻转 | fliplr(A) | fliplr(A) |
| 矩阵上下翻转 | flipud(A) | flipud(A) |

以下是一些数值示例,展示了如何在命令窗口中执行矩阵与标量的运算,包括加法、减法、幂、乘法和除法,以及数组运算:

>> A=[8,1,6; 3,5,7; 4,9,2]  % Matrix 3-by-3
A =
         8     1     6
         3     5     7
         4     9     2
>> a = 2; b = 2+3i; c = 5j;
>> B = A^a % NB:  Difference between .^ and ^
B =
    91    67    67
    67    91    67
    67    67    91
>> C = A.^a % Elementwise op. NB:  Difference between .^ and ^
C =
    64     1    36
     9    25    49
    16    81     4
>> D = A*a+B/b
D =
  30.0000 -21.0000i  12.3077 -15.4615i  22.3077 -15.4615i
  16.3077 -15.4615i  24.0000 -21.0000i  24.3077 -15.4615i
  18.3077 -15.4615i  28.3077 -15.4615i  18.0000 -21.0000i
>> E = C./c
E =
   0.0000 -12.8000i   0.0000 - 0.2000i   0.0000 - 7.2000i
   0.0000 - 1.8000i   0.0000 - 5.0000i   0.0000 - 9.8000i
   0.0000 - 3.2000i   0.0000 -16.2000i   0.0000 - 0.8000i
>> F = C/c
E =
   0.0000 -12.8000i   0.0000 - 0.2000i   0.0000 - 7.2000i
   0.0000 - 1.8000i   0.0000 - 5.0000i   0.0000 - 9.8000i
   0.0000 - 3.2000i   0.0000 -16.2000i   0.0000 - 0.8000i
6. 执行矩阵运算示例

给定六个数组: A (4×3)、 B (3×4)、 C (4×4)、 D (4×3)、 E (3×3) 和 F (3×3),我们可以执行一些矩阵运算,如求和、减法、乘法、幂、标量乘法、平方根、均值、四舍五入、标准差以及矩阵的复制、旋转和翻转等操作:

>> A=[2 -3 1; 3 2 5; 1 3 4; -3 -2 3] ;
>> B=[3,4,-2 1;2,5,4,-6;4,-3, 1,2] ;
>> C=[16,2,3,13;5,11,10,8;9 4 7 14;6 15 12 1] ;
>> D=[1 2 3; 2 3 4; 4 3 1; -2 -3 1] ;
>> E=[8, 1, 6; 3, 5, 7; 4, 9, 2];
>> F=[3 7 3; 3 2 8; 9 2 1];
>> M_AB = A*B
M_AB =
     4   -10   -15    22
    33     7     7     1
    25     7    14    -9
    -1   -31     1    15
>> M_BA = B*A
M_BA =
    13    -9    18
    41    28    25
    -6   -19    -1
>> M_S = M_AB-C
M_S =
   -12   -12   -18     9
    28    -4    -3    -7
    16     3     7   -23
    -7   -46   -11    14
>> M_S = M_BA-C
Matrix dimensions must agree.
>> CM=C*M_S     % Not equivalent to M_S*C
CM =
  -179  -789  -416   243
   352  -442  -141  -150
    18  -747  -279    88
   533  -142   -80  -313
>> CM1=M_S*C    % Not equivalent to C*M_S
CM1 =
  -360   -93  -174  -495
   359  -105   -61   283
   196  -252  -149   307
  -357  -354  -390  -599
>> CM2=M_S.*C  % Elementwise operation: NOT equivalent to M_S*C
CM2 =
  -192   -24   -54   117
   140   -44   -30   -56
   144    12    49  -322
   -42  -690  -132    14
>> MDE=M_S./C  % Elementwise operation: NOT equivalent to M_S/C
MDE =
   -0.7500   -6.0000   -6.0000    0.6923
    5.6000   -0.3636   -0.3000   -0.8750
    1.7778    0.7500    1.0000   -1.6429
   -1.1667   -3.0667   -0.9167   14.0000
>> MD=M_S/C   % Not equivalent to M_S./C
CM =
    1.9275    8.5704   -5.6271   -5.8414
    1.4496   -6.4076    1.5420    3.8277
   -1.0389  -10.8246    5.0116    6.9401
   -4.9118  -12.9118   12.6765    3.6765
>> M_AD = A.*D % Elementwise operation: matrix multiplication
M_AD =
     2    -6     3
     6     6    20
     4     9     4
     6     6     3
>> MM_AD = A*D % Error due to size mismatch of [A] and [D]
Error using  *
Incorrect dimensions for matrix multiplication. Check that the
number of columns in the first matrix matches the number of rows
in the second matrix. To perform elementwise multiplication, use
'.*'.
>> M_EF=E.*F     % Elementwise multiplication of square matrices
M_EF =
    24     7    18
     9    10    56
    36    18     2
>> MM_EF=E*F     % Square matrices can be multiplied matrix-wise
MM_EF =
    81    70    38
    87    45    56
    57    50    86
>> Csqrt=sqrt(C) % Not equivalent to sqrtm(C)
Csqrt =
    4.0000    1.4142    1.7321    3.6056
    2.2361    3.3166    3.1623    2.8284
    3.0000    2.0000    2.6458    3.7417
    2.4495    3.8730    3.4641    1.0000
>> Csqrt=sqrtm(C) % Not equivalent to sqrt(C)
Csqrt =
  Columns 1 through 3
   3.8335 - 0.0167i   0.0738 + 0.7839i   0.1262 + 0.3666i
   0.3251 + 0.0011i   2.6850 - 0.0526i   1.6850 - 0.0246i
   1.3123 - 0.0237i   0.7322 + 1.1107i   1.9687 + 0.5194i
   0.5925 + 0.0373i   2.0922 - 1.7477i   1.7997 - 0.8172i
  Column 4
   1.7975 - 1.1337i
   1.1359 + 0.0761i
   1.8178 - 1.6064i
   1.3466 + 2.5276i
>> C_E1 = expm(C)   % Matrix exponential not equal to exp(C)
C_E1 =
   1.0e+14 *
    1.5718    1.3711    1.3622    1.5295
    1.5718    1.3711    1.3622    1.5295
    1.5718    1.3711    1.3622    1.5295
    1.5718    1.3711    1.3622    1.5295
>> C_E2 = exp(C) % Exponential of a matrix: not equal to expm(C)
C_E2 =
   1.0e+06 *
    8.8861    0.0000    0.0000    0.4424
    0.0001    0.0599    0.0220    0.0030
    0.0081    0.0001    0.0011    1.2026
    0.0004    3.2690    0.1628    0.0000
>> S=[A(1,1:3); B(2,1:3);C(3,2:4)]; % Created from the existed
>> Y=[A(1), 1.3];
>> Arot90=rot90(A)      % Matrix rotate
Arot90 =
     1     5     4     3
    -3     2     3    -2
     2     3     1    -3
>>Crep=repmat(C, 2,1)   % Matrix replication/copy
Crep =
     5    11    10     8
     9     4     7    14
     6    15    12     1
    16     2     3    13
     5    11    10     8
     9     4     7    14
     6    15    12     1
>>Bflip=fliplr(B)       % Matrix flip
Bflip =
     1    -2     4     3
    -6     4     5     2
     2     1    -3     4
>> Cud=flipud(Crep)     % Matrix flip up or down
Cud =
     6    15    12     1
     9     4     7    14
     5    11    10     8
    16     2     3    13
     6    15    12     1
     9     4     7    14
     5    11    10     8
    16     2     3    13

许多这些矩阵运算也可以在 Simulink 环境中执行。Simulink 库中包含用于求和、乘除、幂、指数和拼接的块。例如,我们可以在命令窗口中定义矩阵 A D ,然后在 Simulink 中计算它们的和与差:

>> A=[2 -3 1; 3 2 5; 1 3 4; -3 -2 3] ;
>> D=[1 2 3; 2 3 4; 4 3 1; -2 -3 1] ;

需要注意的是,矩阵 A D 是通过命令窗口和工作区定义的,Simulink 中计算的和与使用 MATLAB 命令窗口计算的结果相匹配。

7. 标准矩阵生成器

MATLAB 有许多标准的数组和矩阵生成器,可用于生成各种矩阵。以下是一些常见的矩阵生成函数及其示例:

>> eye(3)
ans =
     1     0     0
     0     1     0
     0     0     1
>> magic(5)         % Magic matrix in a size of 5 by 5
ans =
    17    24     1     8    15
    23     5     7    14    16
     4     6    13    20    22
    10    12    19    21     3
    11    18    25     2     9
>> A=pascal(4)       % Pascal matrix in a size of 4 by 4
A =
     1     1     1     1
     1     2     3     4
     1     3     6    10
     1     4    10    20
>> A=pascal(4,2)       % Pascal matrix in a size of 4 by 4
A =
    -1    -1    -1    -1
     3     2     1     0
    -3    -1     0     0
     1     0     0     0
>> zeros(3)            % Zero matrix in a size of 3 by 3
ans =
     0     0     0
     0     0     0
     0     0     0
>> zeros(2,3)          % Zero matrix in a size of 2 by 3
ans =
     0     0     0
     0     0     0
>> ones(3)             % Ones matrix in a size of 3 by 3
ans =
     1     1     1
     1     1     1
     1     1     1
>> ones(3,2)           % Ones matrix in a size of 3 by 2
ans =
     1     1
     1     1
     1     1
>> eye(4,3)            % Unit diagonal matrix of size 4 by 3
ans =
     1     0     0
     0     1     0
     0     0     1
     0     0     0
>> eye(4,5)
ans =
     1     0     0     0     0
     0     1     0     0     0
     0     0     1     0     0
     0     0     0     1     0
>> rand(2)
ans =
   0.814723686393179   0.126986816293506
   0.905791937075619   0.913375856139019
>> rand(2,3)
ans =
   0.957166948242946   0.800280468888800   0.421761282626275
   0.485375648722841   0.141886338627215   0.915735525189067
>>A=round(rand(2))
A =
     0     1
     0     0
>>A_rep=repmat(A, 2, 3) % replicating the matrix A by making its 
          % replication 2 times of rows and 3 times of columns
A_rep =
     0     1     0     1     0     1
     0     0     0     0     0     0
     0     1     0     1     0     1
     0     0     0     0     0     0
>> randi([-13, 13], 5) % Random integers within [-13, 13]
ans =
    -7    10    12     0    10
     6    -8     5     1    -2
   -11     5    -9    -6     4
    -9    11    -6     8     3
   -10     0     3    -1    -2
>> C=eye(2); B=magic(3); A=ones(4);
>> D=blkdiag(A,B,C)    % combine matrices in diagonal directions to 
     % create a block diagonal matrix.
D =
     1     1      1     1     0     0     0     0     0
     1     1      1     1     0     0     0     0     0
     1     1      1     1     0     0     0     0     0
     1     1      1     1     0     0     0     0     0
     0     0      0     0     8     1     6     0     0
     0     0      0     0     3     5     7     0     0
     0     0      0     0     4     9     2     0     0
     0     0      0     0     0     0     0     1     0
     0     0      0     0     0     0     0     0     1
>> K=reshape(randperm(9), 3,3) % Change the size (reshape) of array % to 
make 3 by 3 matrix by random permutation
K =
     4     9     5
     3     8     6
     7     1     2

此外,还有几十个矩阵生成函数,如 gallery 函数可以生成各种测试矩阵。例如:

>> L=gallery('leslie', s, x)  % This is the 3-by-3 Leslie pop.         
% matrix taken from the model with average birth numbers s(1:n) 
and          % survival rates x(1:n-1)
L =
            0            2            0
          0.5            0            0
            0         0.78            0
>> C = gallery('chebspec',3,1) % Chebyshev spectral differentiation   
% matrix of order 3
C =
     -0.33333           -1      0.33333
            1      0.33333       -1
      -1.3333            4      -3.1667
>> L=gallery('cauchy', s, y) % Returns 3-by-3 matrix, C(i,j) =  
%1/(s(i)+y(j)). Arguments s and y are vectors of length 3. If you 
         % pass in scalars for s and y, they are interpreted as vectors 1:s 
         % and 1:y.
L =
            2       1.2821      0.76923
          0.4      0.35971      0.30303
            2       1.2821      0.76923
>> B = gallery('krylov', rand(3)) % Returns 3 by 3 Krylov matrix
B =
            1       1.3462       1.6804
            1       1.2107       1.578
            1       1.3387       1.7209

要获取有关 gallery 矩阵的更多信息,可以在命令窗口中输入 help gallery doc gallery

通过以上内容,我们详细介绍了 MATLAB 中线性方程求解、矩阵运算以及矩阵生成的相关知识和操作方法,希望对大家有所帮助。在实际应用中,我们可以根据具体需求选择合适的方法和函数来完成相应的计算任务。

MATLAB 线性代数操作与矩阵运算详解

8. 矩阵运算在 Simulink 中的应用

在 Simulink 环境中,许多矩阵运算也能够顺利执行。Simulink 库为我们提供了一系列用于完成矩阵运算的块,涵盖了求和、乘除、幂、指数以及拼接等操作。下面,我们将深入探讨如何在 Simulink 中进行矩阵运算。

操作步骤
1. 定义矩阵 :首先,在命令窗口中定义矩阵 A D

>> A=[2 -3 1; 3 2 5; 1 3 4; -3 -2 3] ;
>> D=[1 2 3; 2 3 4; 4 3 1; -2 -3 1] ;
  1. 计算矩阵的和与差 :在 Simulink 中,利用相应的块来计算矩阵 A D 的和与差。
graph LR
    A[矩阵 A] -->|输入| SumSubtract[求和/差块]
    D[矩阵 D] -->|输入| SumSubtract
    SumSubtract -->|输出| Result[计算结果]

需要注意的是,矩阵 A D 是通过命令窗口和工作区进行定义的,Simulink 中计算得到的和与使用 MATLAB 命令窗口计算的结果是一致的。

  1. 矩阵乘法 :对于矩阵乘法操作,以图 7 - 18 所示为例,需要将 Multiply 块从元素级乘法(. )切换为矩阵乘法( ),如下面的流程图所示:
graph LR
    A[矩阵 A] -->|输入| Multiply[乘法块]
    D[矩阵 D] -->|输入| Multiply
    Multiply -->|输出| Result[计算结果]
    style Multiply fill:#E5F6FF,stroke:#73A6FF,stroke-width:2px
    click Multiply "将 Multiply 块切换为矩阵乘法(*)"

否则,由于 A D 的尺寸不匹配,乘法操作将无法执行。而且,计算结果与使用 MATLAB 计算的结果是相符的。

  1. 指数和平方运算 :在图 7 - 20 所示的操作中,指数和幂运算可以通过一个 Math Function 块来完成,只需选择其 Function type 为 [pow](用于幂运算)和 [square](用于平方运算),具体流程如下:
graph LR
    A[矩阵 A] -->|输入| MathFunction[Math Function 块]
    MathFunction -->|输出| Result[计算结果]
    style MathFunction fill:#E5F6FF,stroke:#73A6FF,stroke-width:2px
    click MathFunction "选择 Function type 为 [pow] 和 [square]"
  1. 矩阵拼接 :使用 Matrix Concatenate 块,我们可以从计算得到的矩阵和(4×3)、平方(4×4)和矩阵除法(4×3)中创建一个新的 4×10 矩阵,其流程如下:
graph LR
    Sum[矩阵和] -->|输入| Concatenate[矩阵拼接块]
    Square[矩阵平方] -->|输入| Concatenate
    Division[矩阵除法] -->|输入| Concatenate
    Concatenate -->|输出| NewMatrix[新矩阵]
9. 其他矩阵操作

除了上述常见的矩阵运算,还有一些其他的矩阵操作可以帮助我们创建新的矩阵。例如,我们可以使用 diag() 函数提取现有矩阵的对角线元素,或者提取矩阵的选定元素来创建新矩阵。

>> E=[8, 1, 6; 3, 5, 7; 4, 9, 2];
>> F=[3 7 3; 3 2 8; 9 2 1];
>> EF = [diag(E), diag(F)]
EF =
     8     3
     5     2
     2     1
10. 总结

通过对 MATLAB 线性代数操作和矩阵运算的详细介绍,我们可以看到 MATLAB 在处理线性代数问题方面具有强大的功能。从嵌入 MATLAB 函数块计算行列式和求解线性方程,到比较不同线性方程求解函数的准确性和效率,再到进行各种矩阵运算和使用标准矩阵生成器,MATLAB 提供了丰富的工具和函数。

在实际应用中,我们可以根据具体需求选择合适的方法和函数。例如,在求解线性方程时,如果对计算速度有较高要求,可以优先考虑使用 linsolve() mldivide 或反斜杠运算符 \ ;如果对计算精度要求较高,除了 RREF() BICG() LSQR() 函数外的其他方法通常能满足需求。在进行矩阵运算时,要注意矩阵的尺寸匹配和不同运算的区别,如矩阵乘法和数组乘法、矩阵幂和数组幂等。

同时,Simulink 环境为矩阵运算提供了直观的图形化界面,方便我们进行系统建模和仿真。通过合理利用 MATLAB 和 Simulink 的功能,我们可以更高效地解决各种线性代数问题,为科学研究和工程应用提供有力支持。

希望本文能够帮助读者更好地理解和掌握 MATLAB 中的线性代数操作和矩阵运算,在实际应用中灵活运用这些知识,提高工作效率和计算精度。

评论
成就一亿技术人!
拼手气红包6.0元
还能输入1000个字符  | 博主筛选后可见
 
红包 添加红包
表情包 插入表情
 条评论被折叠 查看
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值