matlab编程基础

MATLAB实验指导书

线性代数部分

一、基础知识

1.1 常见数学函数

函 数 名 abs（x） 数 学 计 算 功 能 实数的绝对值或复数的幅值 函 数 名 floor（x） 数 学 计 算 功 能 对x朝-∞方向取整 acos（x） 反余弦arcsinx acosh（x） 反双曲余弦arccoshx angle（x） 在四象限内求复数 x 的相角 asin（x） atan（x） 反正弦arcsinx 反正切arctanx asinh（x） 反双曲正弦arcsinhx atan2（x,y） 在四象限内求反正切 atanh（x） 反双曲正切arctanhx ceil（x） conj（x） cos（x） cosh（x） exp（x） fix（x） 对x朝+∞方向取整 求复数x的共轭复数 余弦cosx 双曲余弦coshx 指数函数 ex 对x朝原点方向取整 gcd（m，n） 求正整数m和n的最大公约数 imag（x） 求复数x的虚部 lcm（m，n） 求正整数m和n的最小公倍数 log（x） real（x） 自然对数（以e为底数） 求复数x的实部 log10（x） 常用对数（以10为底数） rem（m，n） 求正整数m和n的m/n之余数 round（x） 对x四舍五入到最接近的整数 sign（x） sin（x） sinh（x） sqrt（x） tan（x） tanh（x） 符号函数：求出x的符号 正弦sinx 反双曲正弦sinhx 求实数x的平方根：正切tanx 双曲正切tanhx x 如：输入 x=[-4.85 -2.3 -0.2 1.3 4.56 6.75]，则： ceil(x)= -4 -2 0 2 5 7 fix(x) = -4 -2 0 1 4 6

floor(x) = -5 -3 -1 1 4 6 round(x) = -5 -2 0 1 5 7

1.2 系统的在线帮助

1 help 命令：

1.当不知系统有何帮助内容时，可直接输入help以寻求帮助: >> help（回车）

2.当想了解某一主题的内容时，如输入：

>> help syntax （了解Matlab的语法规定）

3.当想了解某一具体的函数或命令的帮助信息时，如输入： >> help sqrt （了解函数sqrt的相关信息）

2 lookfor命令

现需要完成某一具体操作，不知有何命令或函数可以完成，如输入：

·1·

>> lookfor line （查找与直线、线性问题有关的函数）

1.3 常量与变量

系统的变量命名规则：变量名区分字母大小写；变量名必须以字母打头，其后

可以是任意字母，数字，或下划线的组合。此外，系统内部预先定义了几个有特殊意

义和用途的变量，见下表：

特殊的变量、常量 ans pi eps inf NaN i，j 取 值 用于结果的缺省变量名 圆周率π的近似值（3.1416） 数学中无穷小（epsilon）的近似值（2.2204e - 016） 无穷大，如 1/0 = inf （infinity） 非数，如 0/0 = NaN （Not a Number），inf / inf = NaN 虚数单位：i = j =?1 1 数值型向量（矩阵）的输入

1．任何矩阵（向量），可以直接按行方式输入每个元素：同一行中的元素用逗号（，）．．．

或者用空格符来分隔；行与行之间用分号（；）分隔。所有元素处于一方括号（[ ]）内；

例1：

>> Time = [11 12 1 2 3 4 5 6 7 8 9 10] >> X_Data = [2.32 3.43；4.37 5.98]

2．系统中提供了多个命令用于输入特殊的矩阵： 函数 compan diag hadamard hankel invhilb kron magic pascal 功 能 伴随阵 对角阵 Hadamard矩阵 Hankel矩阵 Hilbert矩阵的逆阵 Kronercker张量积 魔方矩阵 Pascal矩阵 函数 toeplitz vander zeros ones rand randn eye meshgrid Toeplitz矩阵 Vandermonde矩阵 元素全为0的矩阵 元素全为1的矩阵 元素服从均匀分布的随机矩阵 元素服从正态分布的随机矩阵 对角线上元素为1的矩阵 由两个向量生成的矩阵 功 能 上面函数的具体用法，可以用帮助命令help得到。如：meshgrid(x,y) 输入 x=[1 2 3 4]; y=[1 0 5]; [X,Y]=meshgrid(x, y)，则 X = Y =

1 2 3 4 1 1 1 1 1 2 3 4 0 0 0 0

·2·

1 2 3 4 5 5 5 5 目的是将原始数据x，y转化为矩阵数据X，Y。

2 符号向量（矩阵）的输入

1．用函数 sym定义符号矩阵：

函数sym实际是在定义一个符号表达式，这时的符号矩阵中的元素可以是任何的符号或者是表达式，而且长度没有限制。只需将方括号置于单引号中。 例2：

>> sym_matrix = sym（'[a b c；Jack Help_Me NO_WAY]'） sym_matrix =

[ a, b, c] [Jack, Help_Me, NO_WAY]

2．用函数syms定义符号矩阵

先定义矩阵中的每一个元素为一个符号变量，而后像普通矩阵一样输入符号矩阵。 例3：

>> syms a b c ； >> M1 = sym（'Classical'）； >> M2 = sym（' Jazz'）； >> M3 = sym（'Blues'）；

>> A = [a b c； M1， M2， M3；sym（[2 3 5]）]

A =

[ a, b, c] [Classical, Jazz, Blues]

[ 2, 3, 5]

1.4 数组（矩阵）的点运算

运算符：+（加）、-（减）、./（右除）、.\\（左除）、.^（乘方）， 例4：

>> g = [1 2 3 4]；h = [4 3 2 1]；

>> s1 = g + h, s2 = g.*h, s3 = g.^h, s4 = g.^2, s5 = 2.^h

1.5 矩阵的运算

运算符：+（加）、-（减）、*（乘）、/（右除）、\\（左除）、^（乘方）、’（转置）等； 常用函数：det（行列式）、inv（逆矩阵）、rank（秩）、eig（特征值、特征向量）、rref

（化矩阵为行最简形）

·3·

例5：

>> A=[2 0 –1;1 3 2]; B=[1 7 –1;4 2 3;2 0 1];

>> M = A*B % 矩阵A与B按矩阵运算相乘 >> det_B = det(B) % 矩阵A的行列式

>> rank_A = rank（A） % 矩阵A的秩 >> inv_B = inv（B） % 矩阵B的逆矩阵

>> [V,D] = eig(B) % 矩阵B的特征值矩阵V与特征向量构成的矩阵D >> X = A/B % A/B = A*B，即XB=A，求X >> Y = B\\A % B\\A = B-1*A，即BY=A，求Y

上机练习（一）：

1．练习数据和符号的输入方式，将前面的命令在命令窗口中执行通过；

2．输入A=[7 1 5；2 5 6；3 1 5]，B=[1 1 1; 2 2 2; 3 3 3]，在命令窗口中执行下列表达式，掌握其含义： A(2, 3) A(:,2) A(3,:) A(:,1:2:3) A(:,3).*B(:,2) A(:,3)*B(2,:) A*B A.*B A^2 A.^2 B/A B./A 3．输入C=1:2:20，则C（i）表示什么？其中i=1,2,3,…,10； 4．查找已创建变量的信息，删除无用的变量；

5．欲通过系统做一平面图，请查找相关的命令与函数，获取函数的帮助信息。

-1

二、编程

2.1 无条件循环

当需要无条件重复执行某些命令时，可以使用for循环： for 循环变量t=表达式1 : 达式2 : 表达式3 语句体 end

说明：表达式1为循环初值，表达式2为步长，表达式3为循环终值；当表达式2省略时则默认步长为1；for语句允许嵌套。

例6： 如：矩阵输入程序

生成3×4阶的Hiltber矩阵。 m=input(‘矩阵行数：m=’); for i=1 : 3 n= input(‘矩阵列数：n=’); for j=1 : 4 for i=1:m H（i，j）=1/（i+j-1）； for j=1:n

end disp([‘输入第’,num2str(i),’行，第’, num2str(j),’列元素’]) end A(i, j) = input (‘ ’) end end

·4·

2.2 条件循环

1) if-else-then语句

if-else-then语句的常使用三种形式为：

(1) if 逻辑表达式 (3) if 逻辑表达式1 语句体 语句体1 end elseif 逻辑表达式2 语句体2 (2) if 逻辑表达式1 elseif 逻辑表达式3 语句体1 ? else else

语句体2 语句体n

end end 2) while循环语句

while循环的一般使用形式为：

while 表达式 语句体 end 例7：

用二分法计算多项式方程x3?2x?5= 0在[0，3]内的一个根。 解：

a = 0；fa = -inf； b = 3；fb = inf； while b-a > eps*b x =（a+b）/2； fx = x^3-2*x-5； if sign(fx)== sign(fa) a =x；fa = fx； else

b = x；fb = fx； end end

x

运行结果为：x = 2.0945515148154233

2.3 分支结构

若需要对不同的情形执行不同的操作，可用switch 分支语句： switch 表达式（标量或字符串）

·5·

case 值1

语句体1 case 值2 语句体2 …

otherwise 语句体n end

说明：当表达式不是“case”所列值时，执行otherwise语句体。

2.4 建立M文件

将多个可执行的系统命令，用文本编辑器编辑后并存放在后缀为 .m 的文件中，若在．．．MATLAB命令窗口中输入该m-文件的文件名（不跟后缀.m!），即可依次执行该文件中的多个命令。这个后缀为.m的文件，也称为Matlab的脚本文件(Script File)。

注意：文件存放路径必须在Matlab能搜索的范围内。

2.5 建立函数文件

对于一些特殊用户函数，系统提供了一个用于创建用户函数的命令function，以备用户随时调用。 1．格式：

function [输出变量列表]=fun_name(输入变量列表)

用户自定义的函数体

2．函数文件名为：fun_name，注意：保存时文件名与函数名最好相同； ．．．．．．3．存储路径：最好在系统的搜索路径上。 4. 调用方法：输出参量=fun_name (输入变量) 例8：

计算s = n!，在文本编辑器中输入：

function s=pp(n);

s=1; for i=1:n s=s*i; end

s;

在MATLAB命令窗口中输入：s=pp(5)

结果为s = 120

上机练习（二）：

·6·

1．编写程序，计算1+3+5+7+?+(2n+1)的值（用input语句输入n 值）。

?x?2．编写分段函数 f(x)??2?x?0?0?x?11?x?2 的函数文件，存放于文件ff.m中，计算出其它f(?3)，f(2)，f(?)的值。

三、矩阵及其运算

3.1 矩阵的创建

1.加、减运算

运算符：“＋”和“－”分别为加、减运算符。

运算规则：对应元素相加、减，即按线性代数中矩阵的“十”、“一”运算进行。 例3-1 在Matlab编辑器中建立m文件：LX0701.m

A=[1, 1, 1; 1, 2, 3; 1, 3, 6] B=[8, 1, 6; 3, 5, 7; 4, 9, 2]

A＋B=A+B A-Ｂ=A-B

在Matlab命令窗口建入LX0701，则 结果显示：A+B=

9 2 7 4 7 10 5 12 8 A－B=

-7 0 -5

-2 -3 -4 -3 -6 4

2.乘法

运算符：*

运算规则：按线性代数中矩阵乘法运算进行，即放在前面的矩阵的各行元素，分别与放在后面的矩阵的各列元素对应相乘并相加。

（1）两个矩阵相乘

例3-2 在Mtalab编辑器中建立M文件：LX0702.m

X= [2 3 4 5 1 2 2 1]；

·7·

Y=[0 1 1 1 1 0 0 0 1 1 0 0]；

Z=X*Y 存盘

在命令行中建入LX0702，回车后显示： Z=

8 5 6

3 3 3

（2）矩阵的数乘：数乘矩阵 上例中：a=2*X

则显示：a =

4 6 8 10

2 4 4 2

（3）向量的点乘（内积）：维数相同的两个向量的点乘。 命令：dot 向量点乘函数 例：X=[-1 0 2]；

Y=[-2 -1 1]; Z=dot(X, Y) 则显示：Z = 4

还可用另一种算法：

sum(X.*Y) ans= 4 （4）向量叉乘

在数学上，两向量的叉乘是一个过两相交向量的交点且垂直于两向量所在平面的向量。在Matlab中，用函数cross实现。

命令 cross 向量叉乘函数

例3-3 计算垂直于向量(1, 2, 3)和(4, 5, 6)的向量。 在Mtalab编辑器中建立M文件：LX0703.m a=[1 2 3]; b=[4 5 6]; c=cross(a,b) 结果显示： c=

-3 6 -3

可得垂直于向量(1, 2, 3)和(4, 5, 6)的向量为±(-3, 6, -3) （5）混合积

·8·

混合积由以上两函数实现：

例3-4 计算向量a=(1, 2, 3)、b=(4, 5, 6)和c=(-3, 6, -3) 的混合积a?(b?c) 在Matlab编辑器中建立M文件：LX0704.m a=[1 2 3]; b=[4 5 6]; c=[-3 6 -3]; x=dot(a, cross(b, c)) 结果显示：x = 54

注意：先叉乘后点乘，顺序不可颠倒。 3.矩阵的除法

Matlab提供了两种除法运算：左除（\\）和右除（/）。一般情况下，x=a\\b是方程a*x =b的解，而x=b/a是方程x*a=b的解

例：a=[1 2 3; 4 2 6; 7 4 9]

b=[4; 1; 2]; x=a\\b 则显示：x=

-1.5000

2.0000

0.5000

如果a为非奇异矩阵，则a\\b和b/a可通过a的逆矩阵与b阵得到： a\\b = inv(a)*b b/a = b*inv(a) 4.矩阵乘方 运算符：^ 运算规则：

（1）当A为方阵，p为大于0的整数时，A^P表示A的P次方，即A自乘P次；p为小于0的整数时，A^P表示A的P次方。

p?d11（2）当A为方阵，p为非整数时，则A^P?V????-1

???V?1其中V为A的特征向?dpnn???d11量，??????为特征值矩阵 ?dnn?5.矩阵的转置 运算符：′

运算规则：与线性代数中矩阵的转置相同。

·9·

6.矩阵的逆矩阵

?1?例3-5 求A??2?3?2243??1?的逆矩阵 3??方法一：在Matlab编辑器中建立M文件：LX07051.m A=[1 2 3; 2 2 1; 3 4 3];

inv(A)或A^(-1) 则结果显示为 ans =

1.0000 3.0000 -2.0000 -1.5000 -3.0000 2.5000 1.0000 1.0000 -1.0000

?1?方法二：由增广矩阵B??2?3?2243131000100??0?进行初等行变换 1??在Matlab编辑器中建立M文件：LX07052.m

B=[1, 2, 3, 1, 0, 0; 2, 2, 1, 0, 1, 0; 3, 4, 3, 0, 0, 1]； C=rref(B) %化行最简形

X=C(:, 4:6)

在Matlab命令窗口建入LX07052，则显示结果如下：

C =

1.0000 0 0 1.0000 3.0000 -2.0000 0 1.0000 0 -1.5000 -3.0000 2.5000 0 0 1.0000 1.0000 1.0000 -1.0000 X =

1.0000 3.0000 -2.0000 -1.5000 -3.0000 2.5000 1.0000 1.0000 -1.0000 这就是A的逆矩阵。 7.方阵的行列式

命令： det 计算行列式的值 例3-6 计算上例中A的行列式的值

在Matlab编辑器中建立M文件：LX0706.m A=[1 2 3; 2 2 1; 3 4 3]; D=det(A)

则结果显示为 D = 2

·10·

3.2 符号矩阵的运算

1.符号矩阵的四则运算

Matlab 5.x 抛弃了在4.2版中为符号矩阵设计的复杂函数形式，把符号矩阵的四则运算简化为与数值矩阵完全相同的运算方式，其运算符为：加（＋），减（－）、乘（×）、除（/、\\）等或：符号矩阵的和（symadd），差（symsub），乘 (symmul)。 例3-7 A?sym(?[1/x,1/(x?1);1/(x?2),1/(x?3)]?);

B?sym(?[x,1;x?2,0]?)；

C=B-A D=a\\b 则显示：

C=

x-1/x 1-1/(x+1)

x+2-1/(x+2) -1/(x+3) D=

-6*x-2*x^3-7*x^2 1/2*x^3+x+3/2*x^2

6+2*x^3+10*x^2+14*x -2*x^2-3/2*x-1/2*x^3 2.其他基本运算

符号矩阵的其他一些基本运算包括转置（'）、行列式（det）、逆（inv）、秩（rank）、幂（^）和指数（exp和expm）等都与数值矩阵相同

3.符号矩阵的简化

符号工具箱中提供了符号矩阵因式分解、展开、合并、简化及通分等符号操作函数。 （1）因式分解

命令：factor 符号表达式因式分解函数 格式：factor(s)

说明：s为符号矩阵或符号表达式。常用于多项式的因式分解 例3-8 将x -1分解因式 在Matlab命令窗口建入 syms x

factor(x^9-1) 则显示：ans =

(x-1)*(x^2+x+1)*(x+x^3+1)

例3-9 问入取何值时，齐次方程组

?(1??)x1?2x2?4x3?0??2x1?(3??)x2?x3?0 ??x1?x2?(1??)x3?06

9

有非0解

·11·

解：在Matlab编辑器中建立M文件：LX0709.m

syms k A=[1-k -2 4;2 3-k 1;1 1 1-k];

D=det(A)

factor(D) 其结果显示如下：

D =

-6*k+5*k^2-k^3 ans =

-k*(k-2)*(-3+k)

从而得到：当k=0、k=2或k=3时，原方程组有非0解。 （2）符号矩阵的展开

命令 expand 符号表达式展开函数 格式：expand(s)

说明：s为符号矩阵或表达式。常用在多项式的因式分解中，也常用于三角函数，指数函数和对数函数的展开中

例3-10 将(x+1)、sin(x+y)展开

在Matlab编辑器中建立M文件：LX0710.m syms x y p=expand((x+1)^3) q=expand(sin(x+y)) 则结果显示为 p =

x^3+3*x^2+3*x+1 q =

sin(x)*cos(y)+cos(x)*sin(y) （3）同类式合并

命令：Collect 合并系数函数

格式：Collect(s,v) 将s中的变量v的同幂项系数合并。

Collect(s) s — 矩阵或表达式，此命令对由命令findsym函数返回的默认变量进行同类项合并。

（4）符号简化

命令：simple或simplify 寻找符号矩阵或符号表达式的最简型 格式：Simple(s) s — 矩阵或表达式

说明：Simple(s)将表达式s的长度化到最短。若还想让表达式更加精美，可使用函数Pretty。

格式：Pretty(s) 使表达式s更加精美 例3-11 计算行列式

·12·

3

1aa21bb21cc21dd2d4

a4b4c4的值。

在Matlab编辑器中建立M文件：LX0711.m

syms a b c d

A=[1 1 1 1;a b c d;a^2 b^2 c^2 d^2;a^4 b^4 c^4 d^4]; d1=det(A)

d2=simple(d1) %化简表达式d1

pretty(d2) %让表达式d2符合人们的书写习惯 则显示结果如下：

d1 =

b*c^2*d^4-b*d^2*c^4-b^2*c*d^4+b^2*d*c^4+b^4*c*d^2-b^4*d*c^2-a*c^2*d^4+a*d^2*c^4+a*b^2*d^4-a*b^2*c^4-a*b^4*d^2+a*b^4*c^2+a^2*c*d^4-a^2*d*c^4-a^2*b*d^4+a^2*b*c^4+a^2*b^4*d-a^2*b^4*c-a^4*c*d^2+a^4*d*c^2+a^4*b*d^2-a^4*b*c^2-a^4*b^2*d+a^4*b^2*c

d2 =

(-d+c)*(b-d)*(b-c)*(-d+a)*(a-c)*(a-b)*(a+c+d+b) (-d+c)(b-d)(b-c)(-d+a)(a-c)(a-b)(a+c+d+b)

?1?例3-12 设A??2?3?2243??2?1?，B???53???11???，C??23??3?3??0? 1??求矩阵X，使满足：AXB = C

在Matlab编辑器中建立M文件：LX0712.m A=[1 2 3;2 2 1;3 4 3]; B=[2,1;5 3]; C=[1 3;2 0;3 1]; X=A\\C/B

则结果显示如下： X =

-2.0000 1.0000 10.0000 -4.0000 -10.0000 4.0000

?cos(t)例3-13 计算??sin(t)?sin(t)?? cos(t)?5在Matlab编辑器中建立M文件：LX0713.m syms t

A =[cos(t) -sin(t); sin(t), cos(t)];

B=sympow(A, 5) %计算A的5次幂

·13·

C=simple(B) %化简 pretty(C)

则显示结果如下 B =

[ cos(t)*(cos(t)*(cos(t)*(cos(t)^2-sin(t)^2)-2*sin(t)^2*cos(t))-sin(t)*(sin(t)*(cos(t)^2-sin(t)^2)+2*cos(t)^2*sin(t)))-sin(t)*(sin(t)*(cos(t)*(cos(t)^2-sin(t)^2)-2*sin(t)^2*cos(t))+cos(t)*(sin(t)*(cos(t)^2-sin(t)^2)+2*cos(t)^2*sin(t))),

cos(t)*(cos(t)*(-2*cos(t)^2*sin(t)-sin(t)*(cos(t)^2-sin(t)^2))-sin(t)*(cos(t)*(cos(t)^2-sin(t)^2)-2*sin(t)^2*cos(t)))-sin(t)*(sin(t)*(-2*cos(t)^2*sin(t)-sin(t)*(cos(t)^2-sin(t)^2))+cos(t)*(cos(t)*(cos(t)^2-sin(t)^2)-2*sin(t)^2*cos(t)))]

[ sin(t)*(cos(t)*(cos(t)*(cos(t)^2-sin(t)^2)-2*sin(t)^2*cos(t))-sin(t)*(sin(t)*(cos(t)^2-sin(t)^2)+2*cos(t)^2*sin(t)))+cos(t)*(sin(t)*(cos(t)*(cos(t)^2-sin(t)^2)-2*sin(t)^2*cos(t))+cos(t)*(sin(t)*(cos(t)^2-sin(t)^2)+2*cos(t)^2*sin(t))),

sin(t)*(cos(t)*(-2*cos(t)^2*sin(t)-sin(t)*(cos(t)^2-sin(t)^2))-sin(t)*(cos(t)*(cos(t)^2-sin(t)^2)-2*sin(t)^2*cos(t)))+cos(t)*(sin(t)*(-2*cos(t)^2*sin(t)-sin(t)*(cos(t)^2-sin(t)^2))+cos(t)*(cos(t)*(cos(t)^2-sin(t)^2)-2*sin(t)^2*cos(t)))]

C =

[ cos(5*t), -sin(5*t)] [ sin(5*t), cos(5*t)]

[cos(5 t) -sin(5 t)] [ ] [sin(5 t) cos(5 t) ]

四、秩与线性相关性

4.1 矩阵和向量组的秩以及向量组的线性相关性。

矩阵A的秩是矩阵A中最高阶非零子式的阶数；向量组的秩通常由该向量组构成的矩阵来计算。

命令：rank

格式：rank(A) A为矩阵式向量组构成的矩阵 例4-1 求向量组?1?(1?2?4?(0623)，?5?(2?6233)，?2?(?24?13)，?3?(?1203)，

4)的秩，并判断其线性相关性。

在Matlab编辑器中建立M文件：LX0714.m A=[1 -2 2 3;-2 4 -1 3;-1 2 0 3;0 6 2 3;2 -6 3 4]; B=rank(A)

运行后结果如下： B =

·14·

3

由于秩为3 < 向量个数，因此向量组线性相关。

4.2 向量组的最大无关组

矩阵的初等行变换有三条：

1.交换两行 ri?rj （第i、第j两行交换） 2.第i行的K倍 kri

3.第i行的K倍加到第j行上去 rj?kri

通过这三条变换可以将矩阵化成行最简形，从而找出列向量组的一个最大无关组，Matlab将矩阵化成行最简形的命令是：

命令：rref

格式：rref(A) A为矩阵 例4-2 求向量组a1=(1,-2,2,3)，a2=(-2,4,-1,3)，a3=(-1,2,0,3)，a4=(0,6,2,3)，a5=(2,-6,3,4)的一个最大无关组。

在Matlab编辑器中建立M文件：LX0715.m a1=[1 -2 2 3]'; a2=[-2 4 -1 3]'; a3=[-1 2 0 3]'; a4=[0 6 2 3]';

a5=[2 -6 3 4]';

A=[a1 a2 a3 a4 a5]

format rat %以有理格式输出 B=rref(A) %求A的行最简形 运行后的结果为

A =

1 -2 -1 0 2 -2 4 2 6 -6 2 -1 0 2 3 3 3 3 3 4 B =

1 0 1/3 0 16/9 0 1 2/3 0 -1/9 0 0 0 1 -1/3 0 0 0 0 0 从B中可以得到：向量a1 a2 a4为其中一个最大无关组

·15·

五、线性方程的组的求解

我们将线性方程的求解分为两类：一类是方程组求唯一解或求特解，另一类是方程组求无穷解即通解。可以通过系数矩阵的秩来判断：

若系数矩阵的秩r=n（n为方程组中未知变量的个数），则有唯一解

若系数矩阵的秩r

线性方程组的无穷解 = 对应齐次方程组的通解+非齐次方程组的一个特解；其特解的求法属于解的第一类问题，通解部分属第二类问题。

5.1 求线性方程组的唯一解或特解（第一类问题）

这类问题的求法分为两类：一类主要用于解低阶稠密矩阵 —— 直接法；另一类是解大型稀疏矩阵 —— 迭代法。

5.1.1 利用矩阵除法求线性方程组的特解（或一个解） 方程：AX=b 解法：X=A\\b

?5x1?6x2?x1?5x2?6x3?x2?5x3?6x4例5-1 求方程组??x3?5x4?6x5?x4?5x5??1?0?0的解 ?0?1解：在Matlab编辑器中建立M文件：LX0716.m A=[5 6 0 0 0

1 5 6 0 0 0 1 5 6 0 0 0 1 5 6 0 0 0 1 5]; B=[1 0 0 0 1]';

R_A=rank(A) %求秩 X=A\\B %求解 运行后结果如下 R_A = 5 X =

2.2662 -1.7218 1.0571 -0.5940 0.3188

·16·

这就是方程组的解。 例5-2 求方程组

?x1?x2?3x3?x4?1??3x1?x2?3x3?4x4?4 ??x1?5x2?9x3?8x4?0的一个特解

解：在Matlab编辑器中建立M文件：LX0717.m

A=[1 1 -3 -1;3 -1 -3 4;1 5 -9 -8];

B=[1 4 0]'; X=A\\B X = 0 0

-0.5333 0.6000

5.1.2 利用矩阵的LU、QR和cholesky分解求方程组的解

1.LU分解：

LU分解又称Gauss消去分解，可把任意方阵分解为下三角矩阵的基本变换形式（行交换）和上三角矩阵的乘积。即A=LU，L为下三角阵，U为上三角阵。

则：A*X=b 变成L*U*X=b

∴X=U\\(L\\b) 这样可以大提高运算速度。 命令 [L，U]=lu (A) 例5-3 求方程组

?4x1?2x2?x3?2??3x1?x2?2x3?10 ??8?11x1?3x2的一个特解

解

?4?A??3?11?2?13?1??2?0??b?[2,10,8]?

在Matlab编辑器中建立M文件：LX0718.m A=[4 2 -1;3 -1 2;11 3 0]; B=[2 10 8]';

D=det(A) [L,U]=lu(A)

X=U\\(L\\B)

在Matlab命令窗口建入LX0718，回车后显示结果如下： D =

·17·

0 L =

0.3636 -0.5000 1.0000 0.2727 1.0000 0 1.0000 0 0 U =

11.0000 3.0000 0 0 -1.8182 2.0000 0 0 0.0000

Warning: Matrix is close to singular or badly scaled.

Results may be inaccurate. RCOND = 2.018587e-017. > In D:\\Matlab\\pujun\\lx0720.m at line 4 X =

1.0e+016 * -0.4053 1.4862 1.3511

说明：结果中的警告是由于系数行列式为零产生的。可以通过A*X验证其正确性。 2.Cholesky分解

若A为对称正定矩阵，则Cholesky分解可将矩阵A分解成上三角矩阵和其转置的乘积，即：A?R??R 其中R为上三角阵。

方程 A*X=b 变成 R??R*X?b 所以 X?R\\(R?\\b) 命令 R=chol(A) 3.QR分解

对于任何长方矩阵A，都可以进行QR分解，其中Q为正交矩阵，R为上三角矩阵的初等变换形式，即：A=QR

方程 A*X=b 变形成 QRX=b 所以 X=R\\(Q\\b) 命令 [Q, R]=qr(A) 上例中 [Q, R]=qr(A)

X=R\\(Q\\B)

说明：这三种分解，在求解大型方程组时很有用。其优点是运算速度快、可以节省磁盘空间、节省内存。

5.2 求线性齐次方程组的通解

在Matlab中，函数null用来求解零空间，即满足A·X=0的解空间，实际上是求出解

·18·

空间的一组基（基础解系）。

格式：z = null % z的列向量为方程组的正交规范基，满足Z??Z?I

z?nul(lA,?r?) % z的列向量是方程AX=0的有理基 例5-4 求解方程组的通解：

?x1?2x2?2x3?x4?0??2x1?x2?2x3?2x4?0 ??x1?x2?4x3?3x4?0解：在Matlab编辑器中建立M文件：LX0719.m A=[1 2 2 1;2 1 -2 -2;1 -1 -4 -3];

format rat %指定有理式格式输出 B=null(A,'r') %求解空间的有理基 运行后显示结果如下：

B =

2 5/3 -2 -4/3 1 0 0 1 写出通解：

syms k1 k2

X=k1*B(:,1)+k2*B(:,2) %写出方程组的通解 pretty(X) %让通解表达式更加精美 运行后结果如下： X =

[ 2*k1+5/3*k2] [ -2*k1-4/3*k2] [ k1] [ k2]

% 下面是其简化形式 [2k1 + 5/3k2 ] [ ] [-2k1 - 4/3k2] [ ] [ k1 ] [ ] [ k2 ]

5.3 求非齐次线性方程组的通解

非齐次线性方程组需要先判断方程组是否有解，若有解，再去求通解。因此， 步骤为：第一步：判断AX=b是否有解，若有解则进行第二步

·19·

第二步：求AX=b的一个特解 第三步：求AX=0的通解

第四步：AX=b的通解= AX=0的通解+AX=b的一个特解。

例5-5 求解方程组

?x1?2x2?3x3?x4?1??3x1?x2?5x3?3x4?2 ??2x1?x2?2x3?2x4?3解：在Matlab编辑器中建立M文件：LX0720.m A=[1 -2 3 -1;3 -1 5 -3;2 1 2 -2]; b=[1 2 3]';

B=[A b]; n=4; R_A=rank(A) R_B=rank(B)

format rat

if R_A==R_B&R_A==n %判断有唯一解 X=A\\b

elseif R_A==R_B&R_A

C=null(A,'r') %求AX=0的基础解系 else X='equition no solve' %判断无解 end

运行后结果显示： R_A = 2 R_B =

3 X =

equition no solve 说明该方程组无解

例5-6 求解方程组的通解：

?x1?x2?3x3?x4?1??3x1?x2?3x3?4x4?4 ??x1?5x2?9x3?8x4?0解法一：在Matlab编辑器中建立M文件：LX07211.m A=[1 1 -3 -1;3 -1 -3 4;1 5 -9 -8];

b=[1 4 0]'; B=[A b]; n=4;

R_A=rank(A)

·20·

R_B=rank(B) format rat

if R_A==R_B&R_A==n X=A\\b

elseif R_A==R_B&R_A

C=null(A,'r')

else X='Equation has no solves' end

运行后结果显示为： R_A = 2 R_B = 2

Warning: Rank deficient, rank = 2 tol = 8.8373e-015. > In D:\\Matlab\\pujun\\lx0723.m at line 11 X = 0 0 -8/15 3/5 C =

3/2 -3/4 3/2 7/4 1 0 0 1 所以原方程组的通解为

?3/2???3/4??0???????3/27/40?+k2??+?? X=k1??1??0???8/15????????0??1??3/5?解法二：在Matlab编辑器中建立M文件：LX07212.m A=[1 1 -3 -1;3 -1 -3 4;1 5 -9 -8];

b=[1 4 0]'; B=[A b];

C=rref(B) %求增广矩阵的行最简形，可得最简同解方程组。

运行后结果显示为： C =

1 0 -3/2 3/4 5/4 0 1 -3/2 -7/4 -1/4

·21·

0 0 0 0 0 对应齐次方程组的基础解系为：

?3/2???3/4?????3/27/4?， ?2??? ?1???1??0??????0??1?非齐次方程组的特解为：

?5/4????1/4? ?*???0????0?所以，原方程组的通解为：

X=k1?1+k2?2+?*

六、特征值与二次型

工程技术中的一些问题，如振动问题和稳定性问题，常归结为求一个方阵的特征值和特征向量。

6.1 方阵的特征值与特征向量

设A为n阶方阵，如果数入和n维列向量x使得关系式

Ax??x

成立，则称?为方阵A的特征值，非零向量x称为A对应于特征值入的特征向量。 在Matlab中，用如下几种调用格式来求A的特征值和特征向量。

1.d=eig(A) % d为矩阵A的特征值排成的向量

2.[V, D]=eig(A) % D为A的特征值对角阵，V的列向量为对应特征值的特征向量

（且为单位向量）

3.[V, D]=eig(A，'nobalance?) % 当A中有小到和截断误差相当的元素时，用

nobalance选项，其作用是减少计算误差

??2?例6-1 求矩阵A??0??4?1211??0?的特征值和特征向量 3??解：在Matlab编辑器中建立M文件：LX0722.m A=[-2 1 1;0 2 0;-4 1 3]; [V,D]=eig(A)

运行后结果显示：

V =

-0.7071 -0.2425 0.3015

·22·

0 0 0.9045 -0.7071 -0.9701 0.3015 D =

-1 0 0 0 2 0

0 0 2

即：特征值-1对应特征向量(-0.7071 0 -0.7071)T

特征值2对应特征向量(-0.2425 0 -0.9701)和(-0.3015 0.9045 -0.3015)

??1?例6-2 求矩阵A???4?1?1300??0?的特征值和特征向量 2??T

T

解：在Matlab编辑器中建立M文件：LX0723.m A=[-1 1 0;-4 3 0;1 0 2];

[V,D]=eig(A)

运行后结果显示为 V =

0 0.4082 -0.4082 0 0.8165 -0.8165

1.0000 -0.4082 0.4082 D =

2 0 0 0 1 0 0 0 1

说明：当特征值为1 (二重根)时，对应特征向量都是k (0.4082 0.8165 -0.4082)k为任意常数。

T

6.2 正交矩阵及二次型

A为n阶方阵，且满足：A?A?E(即A'?A?1)，则A为正交矩阵

A为正交矩阵?A的各列(行)向量的长度为1，而且A的各列(行)向量两两正交。 6.2.1向量的长度（范数）

命令 norm

格式 norm(X) % 求X的范数 6.2.2求矩阵的正交矩阵

命令：orth

格式：orth(A) %将矩阵A正交规范化

·23·

?4?例如：求A??0?0?0310??1?的正交矩阵 3??在Matlab命令窗口键入 A=[4 0 0; 0 3 1; 0 1 3]; P=orth(A)

Q=P'*P

则显示结果为 P =

1.0000 0 0 0 0.7071 -0.7071 0 0.7071 0.7071 Q =

1.0000 0 0 0 1.0000 0 0 0 1.0000 6.2.3矩阵的schur分解

格式 [U, T]=schur(A) %U为正交矩阵，使得A?UTU?和U?U?E T = schur(A) %生成schur矩阵T。当A为实对称阵时，T为特征值对角阵

?4?例6-3 设A??0?0?0310??-1

1?，求一个正交矩阵P，使PAP=Λ为对角阵 3??解法1：在Matlab编辑器中建立M文件：LX07241.m A=[4 0 0;0 3 1;0 1 3]; [V,D]=eig(A)

运行后结果如下： V =

1.0000 0 0 0 0.7071 0.7071 0 0.7071 -0.7071 D =

4 0 0

0 4 0 0 0 2

这里，V就是所求的正交矩阵P，D就是对角矩阵Λ 解法2：在Matlab编辑器中建立M文件：LX07242.m A=[4 0 0;0 3 1;0 1 3]; [U,T]=schur(A)

·24·

运行后结果显示如下：

U =

1.0000 0 0 0 0.7071 0.7071 0 0.7071 -0.7071 T =

4 0 0

0 4 0 0 0 2

这里，U就是所求的正交矩阵P，T就是对角矩阵Λ 说明：对于实对称矩阵，用eig和schur分解效果一样。 例6-4 求一个正交变换X=PY，把二次型

f?2x1x2?2x1x3?2x1x4?2x2x3?2x2x4?2x3x4

化成标准形。

解：先写出二次型的实对称矩阵

?0?1A???1???110?111?101?1??1? 1??0?在Matlab编辑器中建立M文件：LX0725.m A=[0 1 1 -1;1 0 -1 1;1 -1 0 1;-1 1 1 0]; [P,D]=schur(A) syms y1 y2 y3 y4 y=[y1;y2;y3;y4]；

X=vpa(P,2)*y %vpa表示可变精度计算，这里取2位精度 f=[y1 y2 y3 y4]*D*y 运行后结果显示如下： P =

780/989 780/3691 1/2 -390/1351 780/3691 780/989 -1/2 390/1351 780/1351 -780/1351 -1/2 390/1351 0 0 1/2 1170/1351 D =

1 0 0 0 0 1 0 0 0 0 -3 0 0 0 0 1 X =

[ .79*y1+.21*y2+.50*y3-.29*y4]

·25·

[ .21*y1+.79*y2-.50*y3+.29*y4] [ .56*y1-.56*y2-.50*y3+.29*y4] [ .50*y3+.85*y4] f =

y1^2+y2^2-3*y3^2+y4^2

即 f = y12

22?y22?3y3?y4

·26·

本文来源：https://www.bwwdw.com/article/fu3x.html