MATLAB软件融入线性代数的三个应用案例

Three Application Cases of MATLAB Integrated with Linear Algebra

摘要: 本文探讨了将MATLAB软件应用于线性代数教学的重要性，并对其在教学过程中可能遇到的挑战进行了分析。通过三个具体的线性代数应用案例，展示了MATLAB软件在实际教学中的有效性，为如何在线性代数教学中更好地利用MATLAB提供了参考。

Abstract: This article discusses the importance of applying MATLAB software to the teaching of linear algebra and analyzes the challenges that may be encountered in the teaching process. Through three specific cases of linear algebra applications, it demonstrates the effectiveness of MATLAB software in teaching, providing a reference for how to better utilize MATLAB in the teaching of linear algebra.

1. 引言

线性代数是一门应用性很强、理论高度抽象的数学课程。在新工科的大背景下，线性代数越来越多地应用到工程、经济、信息等各个领域上，这种“需求牵引”使得它的重要性大大提高；另一方面，几十年来计算机软硬件的飞速发展，极大的提高了以线性代数为理论基础的科学计算能力，作为“技术推动”，用软件工具加强线性代数教学就显得十分重要[1]。

MATLAB是一款由MathWorks公司开发的多功能数值计算软件，具备强大的图形用户界面和丰富的内置函数库，能够进行高效的数值分析和可视化处理，广泛应用于工程计算、数据分析、算法开发、模型设计等领域，是科研和工程教育中不可或缺的工具。

近三十年来，在线性代数教学中引入MATLAB一直是该课程教学模式改革的一个重点课题[2]。“线性代数 + MATLAB”不仅可以有效减少繁琐的计算过程，其图形功能也极大降低了线性代数的抽象性，使学生更易于对知识点的理解与吸收。但是，国内在MATLAB教学中的应用大多围绕基本理论的学习，如矩阵的运算，线性方程组的求解，向量空间，特征值的求解等，缺少线性代数在科学生产中的实例问题。学生可能仅仅学会了如何使用MATLAB简化计算，而不能用其解决实际问题，缺少建模思想[3]。

2009年，西安电子科技大学负责的教育部“使用信息技术工具改造课程项目”——用MATLAB和建模实践改造工科线性代数子项目带动了全国18所高等院校进行教学改革，并取得了很好的教学效果和社会影响[3]。本文在前辈教师的工作基础上，引入三个线性代数的应用案例，帮助学生更好的了解利用“MATLAB + 线性代数”解决实际问题。

2. MATLAB在线性代数中的应用

下面分别举例说明MATLAB在线性代数中的应用案例。

2.1. 动漫技术的应用

借助计算机实现动漫技术是一项新兴的极具竞争力的产业。线性代数理论在该技术中起到非常重要的作用。比如把一个空间物体用多个顶点图形合成，把物体的动作分别用各个顶点的位置变化来表述，再把这些顶点在两个时间点之间的位移进行插补，使得动作可以用较小的步长连续实现，最后把立体形象投影到屏幕上。在这个过程中，很多环节都和线性代数知识相关。实际问题很复杂，以下给出一个简单平面运动的插补案例[4]。

例1 (1) 构造并绘制大写字母A的图形，然后对其逆时针旋转$\frac{3}{4}\pi$ ，向上平移30，再向右平移20，构造平移和旋转矩阵，并绘制移动后的图形。

(2) 构造微变化矩阵，通过动画描述(1)中大写字母A的运动过程。

解：(1) 用矩阵X来描述平面坐标系中的一个闭合图形(即刚体)，X的第i个列向量表示图形第i个顶点的坐标。为了实现刚体的平移及旋转运算，给矩阵X添加元素值都为1的一行，使矩阵X的形状为3 × n (图形共有n-1个顶点)。

若有矩阵：$M=\left[\begin{array}{ccc}1& 0& c_1\\ 0& 1& c_2\\ 0& 0& 1\end{array}\right]$ ，$R=\left[\begin{array}{ccc}\cos t& -\sin t& 0\\ \sin t& \cos t& 0\\ 0& 0& 1\end{array}\right]$ ，且：$Y_1=MX$ ，$Y_2=RX$

可以证明，矩阵Y₁是刚体X沿x轴正方向平移c₁，沿y轴正方向平移c₂后的结果；矩阵Y₂是刚体X以坐标原点为中心逆时针旋转t弧度的结果[5]。

表1给出大写字母A的11个顶点坐标值，其中最后一列与第一列相同，表示一个完整的闭合图形。

Table 1. The vertex coordinates of the letter A

表1. 字母A的顶点坐标

x	0	4	6	10	8	5	3.5	6.1	6.5	3.2	2	0
y	0	14	14	0	0	11	6	6	4.5	4.5	0	0

构造刚体矩阵$X=\left[\begin{array}{cccccccccccc}0& 4& 6& 10& 8& 5& 3.5& 6.1& 6.5& 3.2& 2& 0\\ 0& 14& 14& 0& 0& 11& 6& 6& 4.5& 4.5& 0& 0\\ 1& 1& 1& 1& 1& 1& 1& 1& 1& 1& 1& 1\end{array}\right]$ 旋转矩阵$R=\left[\begin{array}{ccc}\cos \frac{3\pi }{4}& -\sin \frac{3\pi }{4}& 0\\ \sin \frac{3\pi }{4}& \cos \frac{3\pi }{4}& 0\\ 0& 0& 1\end{array}\right]$ 及平移矩阵$M=\left[\begin{array}{ccc}1& 0& 20\\ 0& 1& 30\\ 0& 0& 1\end{array}\right]$ 。

在MATLAB的M文件编辑器中编写程序11-1.m

% 刚体的平面运动 close all X=[0,4,6,10,8,5,3.5,6.1,6.5,3.2,2,0;0,14,14,0,0,11,6,6,4.5,4.5,0,0;ones(1,12)]; %构造原刚体矩阵 M=[1,0,20;0,1,30;0,0,1]; %构造平移转矩阵 R=[cos(3*pi/4),-sin(3*pi/4),0;sin(3*pi/4),cos(3*pi/4),0;0,0,1]; %构造旋转矩阵 Y=M*R*X; %通过线性变换计算移动后的刚体矩阵 plot(X(1,:),X(2,:)); %绘制原来刚体 hold on axis equal fill(Y3(1,:),Y3(2,:),'black'); %绘制旋转及平移后刚体 grid on,hold off

在MATLAB命令窗口中输入：11-1。

绘制图形如图1所示。

Figure 1. The images of the letter A before and after moving

图1. 字母A移动前后图片

(2) 空心A图形(X矩阵)经过N次微运动变化到实心A图形(Y矩阵)，可以理解为用微变化矩阵K对X进行了N次左乘，即：$Y=\left(K^N\right)\ast X$ ，则有$K^N=M\ast R$ ，则微运动变化矩阵可以用MATLAB命令求得：$K={\left(M\ast R\right)}^{1/N}$ ，即K是对矩阵(M*R)开N次方的结果。

在MATLAB的M文件编辑器中编写程序l1-2.m

clear,close all X=[0,4,6,10,8,5,3.5,6.1,6.5,3.2,2,0;0,14,14,0,0,11,6,6,4.5,4.5,0,0;ones(1,12)]; %构造原始图形矩阵 M=[1,0,20;0,1,30;0,0,1]; %构造平移转矩阵 R=[cos(3*pi/4),-sin(3*pi/4),0;sin(3*pi/4),cos(3*pi/4),0;0,0,1]; %构造旋转矩阵 Y=M*R*X; %通过线性变换计算移动后的刚体矩阵 plot(X(1,:),X(2,:)) %画出原始图形 axis([-10,30,0,36]) %设置显示范围 hold on,grid on axis equal N=8; %设置分解次数 K=(M*R)^(1/N); %构造微变换矩阵 Qs=X; for i=1:N-1 Qs=K*Qs; %求出一次微转动后图形矩阵 F=fill(Qs(1,:),Qs(2,:),'r');pause(1) %画出微转动后图形 %set(F,{'LineStyle'},{'none'}) %F=fill(Qs(1,:),Qs(2,:),'w'); %消隐图形 %set(F,{'LineStyle'},{'none'}) end fill(Y(1,:),Y(2,:),'k') %画出最后图形 hold off

在MATLAB命令窗口中输入：l1-2。

MATLAB绘制出字母A经过8次微运动变化过程图，如图2所示。如果添加上消隐语句，可以得到更好的动画效果。

Figure 2. Micro-movement process of the letter A

图2. 字母A微移动过程

2.2. 太阳因子的计算

在目标与背景的红外辐射特性研究中，常常要分析太阳的辐射因素，而太阳因子的计算是一个非常重要的工作(太阳因子是面元的外法向与太阳光线夹角的余弦值)。以下给出一个利用线性代数知识和MATLAB软件计算建筑物表面太阳因子的应用案例。

例2 假设有一个半球形状的建筑物，计算其外表面各面元的太阳因子。

解：首先建立建筑物的几何模型，并进行网格划分，如图3所示。

Figure 3. Grid chart of hemispherical facade

图3. 半球型建筑表面网格划分图

建立固定在地面的坐标系，x轴为北方向，y轴为东方向，z为竖直向上方向。$\theta$ 值分别取0˚、10˚、20˚、……70˚、80˚九个不同的值；$\phi$ 值分别取0˚、10˚、20˚、……、340˚、350˚三十六个不同的值。把半球表面分为9 × 36 = 324个小面元，容易得到这324个面小元的外法向向量：

$\left(\cos \theta \cos \phi ,\cos \theta \sin \phi ,\sin \theta \right)$

若已知某时刻太阳的高度角$S\text{\_}\theta$ 和方位角$S\text{\_}\phi$ ，如图4所示。则可以算出太阳光线的方向向量为：

$\left(\cos S\text{\_}\theta \cos S\text{\_}\phi ,\cos S\text{\_}\theta \sin S\text{\_}\phi ,\sin S\text{\_}\theta \right)$

Figure 4. Schematic diagram of solar altitude and azimuth angles

图4. 太阳高度角和方位角示意图

最后，把每个面元的外法向向量和此时的太阳方向向量进行内积运算，从而求出该时刻每个面元的太阳因子。

在MATLAB的M文件编辑器中编写程序l2.m。

% 将9×36个面元的外法线方向向量放入三维数组v中。 for i=1:9 for j=1:36 v(i,j,1)=cos((i-1)*10*pi/180)*cos((j-1)*10*pi/180); v(i,j,2)=cos((i-1)*10*pi/180)*sin((j-1)*10*pi/180); v(i,j,3)=sin((i-1)*10*pi/180); end end s_theta=25; %给太阳的高度角赋值。 s_alpha=75; %给太阳的方位角赋值。 s_v=[cos(s_theta*pi/180)*cos(s_alpha*pi/180),cos(s_theta*pi/180)*sin(s_alpha*pi/180),sin(s_theta*pi/180)]; %计算此时太阳的方向向量s_v。 %根据两个向量的内积来计算9×36个面元的太阳因子，用2维数组d_g来存放。 for i=1:9 for j=1:36 d_g(i,j)=[v(i,j,1),v(i,j,2),v(i,j,3)]*s_v'; end end %让太阳因子为负值的点变为0(即太阳照射不到的位置) for i=1:9 for j=1:36 if d_g(i,j)<0 d_g(i,j)=0; end end end

在MATLAB命令窗口中输入：l2。在MATLAB命令窗口中输入：12。即可以计算出半球形建筑物表面所有面元该时刻的太阳因子，具体数据放入了数组d_g中。太阳因子数据的获得，为下一步建筑物表面红外辐射特性的研究奠定了基础。

2.3. 旅行线路的统计

线性代数的矩阵运算常常可以解决生活中的很多问题，下面这个案例是一个关于旅行路线统计的问题。

例题3：李博士从西安要到a、b、c、d、e五个国家去旅行，规定在每一个国家只能停留在一个城市。图5给出李博士从西安到五个国家不同城市的路线图。从图中可以看出，a国有3个城市可以选择，从西安到a1城市有2种不同的途径，到a2有3种不同的途径，到a3有1种途径；b国有2个城市可以选择……问李博士从西安游遍这五个国家共有多少不同的旅行方式？

Figure 5. Map of routes from Xi'an to different cities in five countries

图5. 西安至五国不同城市路线图

解：构造两个相邻国家之间旅行路线矩阵。

从西安到a国的情况可以用下表或矩阵A描述：

	a1	a2	a3
西安	2	3	1

$A=\left(\begin{array}{ccc}2& 3& 1\end{array}\right)$

从a国到b国的情况可以用下表或用矩阵B描述：

	b1	b2
a1	3	1
a2	2	2
a3	0	3

$B=\left(\begin{array}{cc}3& 1\\ 2& 2\\ 0& 3\end{array}\right)$

从b国到c国的情况可以用下表或用矩阵C描述：

	c1	c2	c3	c4
b1	3	3	1	0
b2	0	0	2	2

$C=\left(\begin{array}{cccc}3& 3& 1& 0\\ 0& 0& 2& 2\end{array}\right)$

从c国到d国的情况可以用下表或矩阵D描述：

	d1	d2	d3
c1	2	1	0
c2	1	2	3
c3	0	0	3
c4	0	0	2

$D=\left(\begin{array}{ccc}2& 1& 0\\ 1& 2& 3\\ 0& 0& 3\\ 0& 0& 2\end{array}\right)$

从d国到e国的情况可以用下表或用矩阵E描述：

	e1	e2	e3	e4	e5
d1	2	4	5	0	0
d2	0	0	3	0	0
d3	0	0	2	2	2

$E=\left(\begin{array}{ccccc}2& 4& 5& 0& 0\\ 0& 0& 3& 0& 0\\ 0& 0& 2& 2& 2\end{array}\right)$

根据矩阵乘法规则，知$A\ast B=\left(\begin{array}{ccc}2& 3& 1\end{array}\right)\ast \left(\begin{array}{cc}3& 1\\ 2& 2\\ 0& 3\end{array}\right)=\left(\begin{array}{cc}12& 11\end{array}\right)$ ，12和11分别代表从西安到b国的b1和b2两个不同城市的路线总数，于是有A*B*C*D*E的乘积代表从西安到e国e1、e2、e3、e4和e5五个不同城市的路线总数。

M文件编辑器中编写程序11-3.m。

A=[2,3,1]; B=[3,1;2,2;0,3]; C=[3,3,1,0;0,0,2,2]; D=[2,1,0;1,2,3;0,0,3;0,0,2]; E=[2,4,5,0,0;0,0,3,0,0;0,0,2,2,2]; X=A*B*C*D*E sum(X)

在MATLAB命令窗口中输入：l1-3。

最后输出为：

X =

216 432 1372 508 508

ans =

3036

即从西安到E国家5个不同城市的路线总数由X给出，5个城市的路线总和为3036条路线。

3. 结论

上述实例涵盖了线性代数的核心概念，包括矩阵乘法、线性变换和向量空间等，它们是多个学科领域的基础，如计算机科学、红外物理学和运筹学等。这些实例不仅展示了线性代数在解决实际问题中的应用能力，而且更重要的是，它们启发了学生利用矩阵进行建模的思维方式，让学生意识到在自己的专业领域内，线性代数是一个强大的工具，配合MATLAB的计算能力，不仅能简化运算，还能够解决各种复杂问题。

基金项目

西安电子科技大学中央教改专项“线性代数课程智能化学习系统建设”(项目号：20901230005)。

NOTES

^*通讯作者。