MatLab中的矩阵

MatLab中有好多函数可以产生不同的矩阵,下面就让我们产生两个3*3的矩阵,这一章中,我们的学习就靠她们了!!!

A = pascal(3)
A =
　　 1 1 1
　　 1 2 3
　　 1 3 6
B = magic(3)
B =
　　 8 1 6
　　 3 5 7
　　 4 9 2

还有一个3*2的随机矩阵:

C = fix(10*rand(3,2))
C =
　　 9 4
　　 2 8
　　 6 7

看看列矩阵，行矩阵，以及常数的表达：

u = [3; 1; 4]
v = [2 0 —1]
s = 7

产生的矩阵是：

u =
　　 3
　　 1
　　 4
v =
　　 2 　0 　—1
s =
　　 7

加减法

X = A + B
X =
　　 9 2 7
　　 4 7 10
　　 5 12 8
Y = X –A
Y =
　　 8 1 6
　　 3 5 7
　　 4 9 2

若二矩阵维数不统一，则会出错！

X = A + C

Error using ==> +
Matrix dimensions must agree.

向量的乘积与转置

x = v*u
x =
　　 2
X = u*v
X =
　　 6 0 —3
　　 2 0 —1
　　 8 0 —4
X = B'
X =
　　 8 3 4
　　 1 5 9
　　 6 7 2

如ｘ与ｙ均是列向量，则x*y无解,但下二表达式却可以:

x'*y
y'*x

称内积或点积.

下面的语句产生单位矩阵

eye(m,n)

若用eye(n)则产生n*n的方阵

解线性方程

情况一:

x = A\u
x =
　　10
　　 —12
　　 5

又如：

X = A\B
X =
　　 19 –3 —1
　　 —17 4　 13
　　 　6 　0 —6

情况二；y是不同时刻t时的观测值:

t = [0 .3 .8 1.1 1.6 2.3]';
y = [.82 .72 .63 .60 .55 .50]';

若函数形式是:y(t)=c1+c2*exp(t);
构造矩阵:

E = [ones(size(t)) exp(–t)]
E =
　　 1.0000 1.0000
　　 1.0000 0.7408
　　 1.0000 0.4493 　　
　　 1.0000 0.3329
　　 1.0000 0.2019 　　
　　 1.0000 0.1003

则可求得系数c1及c2

c = E\y
c =
　　0.4760 0.3413

表明：y(t)=0.4760+0.3413*exp(t)
画图如下:

T = (0:0.1:2.5)';
Y = [ones(size(T)) exp(–T)]*c;
plot(T,Y,'–',t,y,'o')

转置与行列式

若A是方阵,且是非奇异的,则:

d = det(A)
X = inv(A)
d =
　　1
X =
　　 3 　—3 　1
　　 —3 　5 　—2
　　 1 　—2 　1

若c不是方阵,则用pinv:

X = pinv(C)
X =
　　 0.1159 　—0.0729 　0.0171
　　 —0.0534 　0.1152 　0.0418

那么我们可以发现，下面３个命令具有同样的功效（Ａ是ｍ*ｎ的矩阵，ｍ＞ｎ）：

x = A\b
x = pinv(A)*b
x = inv(A’*A)*A’*b

LU.RQ.及Cholesky分解

MatLab求解线性方程建立在以下三个分解之上:
Cholesky分解
Guass(高斯)分解
正交分解

Cholesky分解

A=p*p'

让我们临时把A变一变:

A = pascal(6)
A =
　　 1 　1 　1 　1 　1　　1
　　 1 　2　 3 　4 　5　　6
　　 1 　3 　6 　10　15 　21
　　 1 　4 　10　20　35 　56
　　 1 　5 　15　35　70 　126
　　 1 　6 　21　56　126　252

A是二项式系数,每一项是其左方与上方系数之和,求其Cholesky分解系数有:

R = chol(A)
R =
　　1 1 1 1 1 1
　　0 1 2 3 4 5
　　0 0 1 3 6 10
　　0 0 0 1 4 10
　　0 0 0 0 1 5
　　0 0 0 0 0 1

Ｒ认识二项式系数．
这样对于线性方程便可化简：
A*x = b
R'*R*x = b
　　x = R\(R'\b)
复杂度由O(n^3)变为O(n^2);

LU分解

A = L U

其中,L时下三角阵,U是上三角阵,如:

[L,U] = lu(B)
L =
　　1.0000 0　　　　0
　　0.3750 0.5441 　1.0000
　　0.5000 1.0000 　0
U =
　　8.0000　1.0000　6.0000
　　0 　　　8.5000　—1.0000
　　0 　　　0 　　　5.2941

同样：

A*x = b可以解为
x = U\(L\b)

QR分解

正交阵有如下性质:

Q'Q = I

正交阵的好处在于,她保持了原阵的长度,角度,并且在计算的过程中不会扩大误差.

RQ分解如下:

A = Q R 　　　或　　　　　　A P = Q R

其中,Q是正交阵,R是上三角阵.

矩阵的幂与指数

若A是方阵,p是正数,则

X = A^2
X =
　　3 6 10
　　6 14 25
　　10 25 46

若Ａ是方阵，且是非奇异的，则X=A^(-P)将inv(A) P次方,如:

Y = B^(–3)
Y =
　　0.0053 　—0.0068 0.0018
　　—0.0034 0.0001 　0.0036
　　—0.0016 0.0070 　—0.0051

分数词幂将由Ａ的特征值决定．

若是对矩阵的每个元素进行幂，用.^,如

X = A.^2
A =
　　1 1 1
　　1 4 9
　　1 9 36

sqrtm(A)计算A^(1/2)，但要更精确，而
sqrt(A)则计算A.^(1/2),是一个元素一个元素的算．

dx/dt=Ax,可以表示为x(t)=exp(tA)*x(0);
下面来看看如何计算:--expm(A)

A =
　　0 —6　 —1
　　6 　2 　—16
　　—5 20　—10
x0 =
　　 1
　　 1
　　 1

计算如下：

X = [];
for t = 0:.01:1
　　X = [X expm(t*A)*x0];
end

作图有：

plot3(X(1,:),X(2,:),X(3,:),'–o')

特征值

Av=λv

若L是特阵值矩阵,则特征向量是V:　　

AV=VL;

如下:

A =
　　 0 —6　 —1
　　 6 　2　 —16
　　 —5 20　—10
lambda = eig(A)
lambda =
　　—3.0710
　　—2.4645+17.6008i
　　—2.4645-17.6008i

由exp(λt)可以看出exp(At)(见上小节)

若用二参数调用函数eig(),则返回特征向量及特征值矩阵:

[V,D] = eig(A)
V =
　　—0.8326 　—0.1203+ 0.2123i 　—0.1203– 0.2123i
　　—0.3553 　0.4691+ 0.4901i 　　0.4691– 0.4901i
　　—0.4248 　0.6249– 0.2997i　　0.6249+ 0.2997i
D =
　　—3.0710 　　0 　　　　　　　　0
　　0 　　　　　—2.4645+17.6008i　0
　　0 　　　　　　0 　　　　　　　—2.4645—17.6008i

对于下面的矩阵：

A =
　　6 　　12 　　19
　　—9 —20　 —33
　　4 　　9　　　15
V =
　　0.4741 　0.4082　 —0.4082
　　—0.8127 —0.8165 0.8165
　　0.3386 　0.4082 　—0.4082
D =
　　—1.0000 　0　　 　0
　　　0 　　　1.0000 　0
　　　0　　　 0　　　　1.0000

可以看出，有二特征值是一样的，其特征向量仅差一个符号，在Symbolic Math Toolbox中提供了Jordan标准型的函数,如下:

[X,J] = jordan(A)
X =
　　—1.7500 1.5000　 2.7500
　　3.0000 　—3.0000 —3.0000
　　—1.2500 1.5000 　1.2500
J =
　　—1 0 0
　　　0 1 1
　　　0 0 1