matlab下的矩阵运算.docx
《matlab下的矩阵运算.docx》由会员分享，可在线阅读，更多相关《matlab下的矩阵运算.docx（16页珍藏版）》请在冰豆网上搜索。
matlab下的矩阵运算
第二讲矩阵运算
§1基本矩阵运算元
我们在第一讲章已说明过MATLAB的运算是以数组（array）及矩阵（matrix）方式在做运算，而这二者在MATLAB的基本运算性质不同，数组强调元素对元素的运算，而矩阵则采用线性代数的运算方式。
我们就来说明矩阵运算的特点。
以下将数组及矩阵的运算符号及其意义列出
数组运算符号
矩阵运算符号
功能
+
+
加
-
-
减
.*
*
乘
./
/
左除
.\
\
右除
.^
^
次方
.'
'
转置
利用这些运算符号即可进行以下的矩阵运算。
>>A=[251;738;4521;16130];
>>A'%A的转置矩阵
A=
27416
53513
18210
>>A=[4-13];B=[-252];
>>dot_prod=sum（A.*B）%二个数组做内积
dot_prod=
-7
>>c=dot（A,B）%以dot函数也可做内积运算
c=
-7
>>A=[4;-1;3];
>>dot_prod=sum（A'.*B）;%如果A是行数组则先做转置，再做内积
>>F=[25-1];G=[01-3];
>>out_prod=F'*G;%二矩阵做外积
>>A=[2,5,1;0,3,-1];
>>B=[1,0,2;-1,4,-2;5,2,1];
>>C=A*B%矩阵相乘，注意二个矩阵的大小须相容
C=
222-5
-810-7
>>A=[21;43];
>>A^2%矩阵次方
ans=
85
2013
下面我们演示一个具体的例子。
假设我们把本地区的天气分为3种状态：
晴，阴，雨。
若今天天阴，则明天天晴的概率为1/2，阴的概率为1/4，下雨的概率为1/4。
如果今天天阴，或者今天下雨，则明天天气是其它情况的概率会是其它的值，将这些概率值列入下面的表中。
天气状态转移概率表
今天
明天
晴
阴
雨
晴
3/4
1/2
1/4
阴
1/8
1/4
1/2
雨
1/8
1/4
1/4
表中的每一列对应于今天天气状态，它的每一行对应于明天天气状态。
例如，第2行第3列（最后一列）的值为1/2,，这给出了今天下雨明天转阴的概率。
将上表内的概率数据用矩阵A表示
矩阵A中概率称为转移概率，矩阵A称为转移矩阵。
已知今天天气晴、阴、雨的概率，可以用转移矩阵A提供的数据来计算明天天气晴、阴、雨的概率。
记
、
、
分别为今天天气是晴、阴、雨的概率，
、
、
分别为明天天气是晴、阴、雨的概率，两个列向量
分别称为今日概率向量和明日概率向量，则有矩阵计算式
以当前状态预测未来状态的概率模型称为Markov链。
如果在清晨我们听到的天气预报为，今天阴或雨的概率都是1/2,那么，今日概率向量
。
利用上式计算明日概率向量的Matlab操作是：
A=[3/41/21/4;1/81/41/2;1/81/41/4];%输入矩阵A
P=[01/21/2]';%输入向量P
P1=（A*P）'
P1=
3/83/81/4
这里，P1就是按行向量的形式输出的结果。
明天天气为晴、阴的概率都是3/8,下雨的概率是1/4。
明日的概率向量可以用前面的Markov链求出，那么两天的概率向量可用公式
给出。
这里需要计算
，我们再用Matlab完成矩阵乘积的计算。
A2=A^2
A2=
21/329/161/2
3/161/49/32
5/323/167/32
用键盘敲入A2=A^2后再按下【Enter】键就得到了结果A2。
从本例可见，Matlab的操作并不复杂，它用于矩阵运算十分方便。
§2线性方程组的求解
由于矩阵的乘法运算一般不具有可交换性，所以AX与XA一般不相同，为求解X，Matlab提供了两种除法运算符：
斜杠“\”和“/”分别表示左除和右除。
左除“\”：
用X=A\B表示AX=B的解；右除“/”：
用X=B/A表示XA=B的解。
例解线性方程组
输入：
A=[21-51;1-30-6;02-12;14-76],b=[89-50]
A=
21-51
1-30-6
02-12
14-76
b=
89-50
方程组的解为
X=（A\b'）'
X=
3.0000-4.0000-1.00001.0000
§3多项式与矩阵多项式
一、多项式行向量的创建方法
在MATLAB里，多项式由一个行向量表示，它的系数是按降序排列。
例如，输入多项式x4－12x3＋0x2＋25x＋116
>>p=[1-12025116]
p=
1-12025116
注意，必须包括具有零系数的项。
除非特别地辨认，MATLAB无法知道哪一项为零。
给出这种形式，用函数roots找出一个多项式的根。
?
>>r=roots（p）
r=
11.7473
2.7028
-1.2251+1.4672i
-1.2251-1.4672i
因为在MATLAB中，无论是一个多项式，还是它的根，都是向量，MATLAB按惯例规定，多项式是行向量，根是列向量。
给出一个多项式的根，也可以构造相应的多项式。
在MATLAB中，命令poly执行这个任务。
>>pp=poly（r）
pp=
1.0e+002*
Columns1through4
0.0100-0.12000.00000.2500
Column5
1.1600+0.0000i
>>pp=real（pp）%throwawayspuriousimaginarypart
pp=
1.0000-12.00000.000025.0000116.0000
因为MATLAB无隙地处理复数，当用根重组多项式时，如果一些根有虚部，由于截断误差，则poly的结果有一些小的虚部，这是很普通的。
消除虚假的虚部，如上所示，只要使用函数real抽取实部。
二、多项式乘法
函数conv支持多项式乘法（执行两个数组的卷积）。
考虑两个多项式a（x）=x3＋2x2＋3x＋4和b（x）=x3＋4x2＋9x＋16的乘积：
>>a=[1234];b=[14916];
>>c=conv（a,b）
c=
162050758464
结果是c（x）=x6＋6x5＋20x4＋50x3＋75x2＋84x＋64。
两个以上的多项式的乘法需要重复使用conv。
三、多项式加法
对多项式加法，MATLAB不提供一个直接的函数。
如果两个多项式向量大小相同，标准的数组加法有效。
把多项式a（x）与上面给出的b（x）相加。
>>d=a+b
d=
261220
结果是d（x）=2x3＋6x2＋12x＋20。
当两个多项式阶次不同，低阶的多项式必须用首零填补，使其与高阶多项式有同样的阶次。
考虑上面多项式c和d相加：
>>e=c+[000d]
e=
162052819684
结果是e（x）=x6＋6x5＋20x4＋52x3＋81x2＋96x＋84。
要求首零而不是尾零，是因为相关的系数象x幂次一样，必须整齐。
如果需要，可用一个文件编辑器创建一个函数M文件来执行一般的多项式加法。
：
functionp=mmpadd（a,b）
%MMPADDPolynomialaddition.
%MMPADD（A,B）addsthepolynomialAandB
%Copyright（c）1996byPrenticeHall,Inc.
ifnargin<2
error（'Notenoughinputarguments'）
end
a=a（:
）.';%makesureinputsarepolynomialrowvectors
b=b（:
）.';
na=length（a）;%findlengthsofaandb
nb=length（b）;
p=[zeros（1,nb-na）a]+[zeros（1,na-nb）b];%addzerosasnecessary
现在，为了阐述mmpadd的使用，再考虑前一页的例子。
>>f=mmpadd（c,d）
f=
162052819684
它与上面的e相同。
当然，mmpadd也用于减法。
>>g=mmpadd（c,-d）
g=
162048697244
结果是g（x）=x6＋6x5＋20x4＋48x3＋69x2＋72x＋44。
四、多项式除法
在一些特殊情况，一个多项式需要除以另一个多项式。
在MATLAB中，这由函数deconv完成。
用上面的多项式b和c
>>[q,r]=deconv（c,b）
q=
1234
r=
0000000
这个结果是b被c除，给出商多项式q和余数r，在现在情况下r是零，因为b和q的乘积恰好是c。
五、多项式导数
由于一个多项式的导数表示简单，MATLAB为多项式求导提供了函数polyder。
>>g
g=
162048697244
>>h=polyder（g）
h=
6308014413872
六、多项式估值
根据多项式系数的行向量，可对多项式进行加，减，乘，除和求导，也应该能对它们进行估值。
在MATLAB中，这由函数polyval来完成。
>>x=linspace（-1,3）;%choose100datapointsbetween-1and3.
>>p=[14-7-10];%usespolynomialp（x）=x3＋4x2－7x－10
>>v=polyval（p,x）;
计算x值上的p（x），把结果存在v里。
然后用函数plot绘出结果。
>>plot（x,v）,title（'x^3+4x^2-7x-10'）,xlabel（'x'）
多项式估值
七、有理多项式
在许多应用中，例如富里哀（Fourier），拉普拉斯（Laplace）和Z变换，出现有理多项式或两个多项式之比。
在MATLAB中，有理多项式由它们的分子多项式和分母多项式表示。
对有理多项式进行运算的两个函数是residue和polyder。
函数residue执行部分分式展开。
>>num=10*[12];%numeratorpolynomial
>>den=poly（[-1;-3;-4]）;%denominatorpolynomial
>>[res,poles,k]=residue（num,den）
res=
-6.6667
5.0000
1.6667
poles=
-4.0000
-3.0000
-1.0000
k=
[]
结果是余数、极点和部分分式展开的常数项。
上面的结果说明了该问题：
这个函数也执行逆运算。
>>[n,d]=residue（res,poles,k）
n=
0.000010.000020.0000
d=
1.00008.000019.000012.0000
>>roots（d）
ans=
-4.0000
-3.0000
-1.0000
在截断误差内，这与我们开始时的分子和分母多项式一致。
residue也能处理重极点的情况，尽管这里没有考虑。
正如前面所述，函数polyder，对多项式求导。
除此之外，如果给出两个输入，则它对有理多项式求导。
>>[b,a]=polyder（num,den）
b=
-20-140-320-260
a=
116102328553456144
该结果证实：
八、小结
下列表格概括了在本讲所讨论的多项式操作特性。
多项式函数
conv（a,b）
乘法
[q,r]=deconv（a,b）
除法
poly（r）
用根构造多项式
polyder（a）
对多项式或有理多项式求导
polyfit（x,y,n）
多项式数据拟合
polyval（p,x）
计算x点中多项式值
[r,p,k]=residue（a,b）
部分分式展开式
[a,b]=residue（r,p,k）
部分分式组合
roots（a）
求多项式的根
九、多项式函数
函数polyvalm是以矩阵方式做多项式函数计算，有别于polyval是以数组方式计算函数值。
它的语法为polyvalm（a,X），其中X为一矩阵而a则是一多项式。
以下的例子可说明其用法。
>>X=[111;222;333];
>>a=[111];%注意a=X*X+X+I
>>f=polyvalm（a,X）
f=
877
141514
212122
§4矩阵函数
矩阵函数按其用法大致可分为矩阵构造函数、矩阵计算函数和矩阵操作函数，其指令名及功能见以下列表：
（矩阵计算函数见下讲）
矩阵构造函数
指令
功能
指令
功能
eye
产生单位矩阵
tril
产生下三角阵
ones
产生全1矩阵
rand
产生均匀分布随机阵
zeros
产生全0矩阵
randn
产生正态分布随机阵
diag
产生对角形数组
hilb
产生Hilbert矩阵
triu
产生上三角阵
magic
产生魔方矩阵
例部分函数演示
A=hilb（3）
A=
1.00000.50000.3333
0.50000.33330.2500
0.33330.25000.2000
diag（A）
ans=
1.0000
0.3333
0.2000
tril（A）
ans=
1.000000
0.50000.33330
●0.33330.25000.2000
矩阵操作函数
指令
功能
指令
功能
reshape
改变矩阵阶数
flipud
矩阵作上下翻转
repmat
按指定的行列数复制矩阵
fliplr
矩阵作左右翻转
rot90
逆时针旋转矩阵90°
例部分函数功能演示
h=1:
6,H=reshape（h,2,3）
h=
123456
H=
135
246
H1=reshape（H,3,2）
H1=
14
25
36
G=repmat（h,2,2）
G=
123456123456
123456123456
fliplr（H）
ans=
531
642
>>rot90（H）
ans=
56
34
12
升级会员 