数组的建立
数组同一行的元素用逗号或空格分割,不同行之间通过分号分割,要注意都是在英文状态下。数组是有方向的,一维数组包括行向量和列向量。
A=[1 2 3 4 5]
B1=A(3) %获取第3个元素
B2=A(2:4) %获取第2个到第4个元素
B3=A(3:end) %获取从第3个开始到最后一位的元素
B4=A(3:-1:2) %获取第3个,第2个元素
B5=A(end:-1:1) %反序输出
B6=A([2 4]) %获取第2个、第4个元素
A=[1,2,3;4,5,6;7,8,9]
B=[1:3;4:6;7:1:9]
C=[A;B] %纵向拼接
D=[A B] %横向拼接
数组的运算
都是对应元素进行运算,比较简单。
数组信息的获取
函数 | 说明 |
---|---|
isempty(A) | 检测数组是否为空,是就返回1,否则返回0 |
isscalar(A) | 检测数组是否为单个元素的标量 |
isvector(A) | 检测数组是否为行向量或列向量 |
isrow() | 检测数组是否为列向量 |
iscolumn() | 检测数组是否为行向量 |
issparse() | 检测数组是否为稀疏矩阵(为0的元素数目远多于不为0的,且非0元素分布没有规律) |
size() | 获取数组的行数和列数 |
length() | 获取一维数组的长度,若是二维数组,返回行数和列数当中较大的一个 |
ndims() | 计算数组的维度(行数+列数) |
isnumeric() | 判断元素是否为数值型 |
isreal() | 判断数组元素是否为实数 |
isinteger() | 判断元素是否为整数 |
islogical() | 判断元素是否为逻辑型 |
whos | 获取数组的大小以及占用内存的多少 |
find() | 返回关系表达式为真的元素下标 |
sort() | 数组按升序排序 |
———————————————————————————————————
矩阵的创建
矩阵是按列存储的,先第一列,再第二列,依次类推,即第二个元素是第二列的第一个元素,第三个元素是第三列的第一个元素,依次往下类推。
对矩阵中元素赋值时,若行或列超出矩阵的大小,则自动扩充矩阵的大小,然后再进行赋值,扩充部分用0填充。
A=[1:4;5:8;9:12;13:16]
ind1=sub2ind(size(A),3,4) %将双下标转化为单下标,3行4列转化为15
A(ind1) %12
A(3,4) %12
[I,J]=ind2sub(size(A),10) %将单下标转化为双下标,I=2,J=3
ind2=sub2ind(size(A),I,J) %10
A =
1 2 3 4
5 6 7 8
9 10 11 12
13 14 15 16
矩阵的扩展
进行数组连接的函数:
函数 | 说明 |
---|---|
cat(DIM,A,B) | 在DIM维度上进行矩阵A和B的连接,返回值为连接后的矩阵 |
vertcat(A,B) | 在水平方向上连接数组A和B,相当于cat(1,A,B) |
horzcat(A,B) | 在垂直方向上连接数组A和B,相当于cat(2,A,B) |
矩阵的块操作
通过函数repmat()进行数据块的复制。
B=repmat(A,m,n):该函数产生大的矩阵B,把矩阵A看作单个元素,产生由m行和n列的矩阵A组成的大矩阵B。
B=repmat(A,m):该函数产生大的矩阵B,把矩阵A看作单个元素,产生由m行和m列的矩阵A组成的大矩阵B。
函数blkdiag()将多个矩阵作为对角块,产生新的矩阵。
如:Y=blkdiag(A,B):将矩阵A,B作为对角块,产生新的矩阵。
A=magic(3)
B=[1:2;3:4]
C=blkdiag(A,B)
D=blkdiag(B,A)
运行结果:
A =
8 1 6
3 5 7
4 9 2
B =
1 2
3 4
C =
8 1 6 0 0
3 5 7 0 0
4 9 2 0 0
0 0 0 1 2
0 0 0 3 4
D =
1 2 0 0 0
3 4 0 0 0
0 0 8 1 6
0 0 3 5 7
0 0 4 9 2
矩阵中元素的删除
利用空矩阵[]删除元素。
注意,X=[]和clear X不同,clear是将X从工作空间中删除,而空矩阵存在于工作空间中,只是维数为0。
A([1 3],:)=[] %删除第1行、第3行
A(:,end)=[] %删除最后一列
矩阵的转置
转置操作符 ’ ,A转置成A’,如果矩阵中含有复数,转置之后复数变为它的共轭复数。
矩阵的真正转置是A.’,即使是复数,也不转化为共轭复数。
也可以用函数transpose()来转化。
A=rand(2,4)
A1=A'
A2=transpose(A)
B=[2+3i,4+5i,3;2 4+i 5+3*i]
B1=B'
B2=B.'
B3=transpose(B)
B4=ctranspose(B)
运行结果:
A =
0.4218 0.7922 0.6557 0.8491
0.9157 0.9595 0.0357 0.9340
A1 =
0.4218 0.9157
0.7922 0.9595
0.6557 0.0357
0.8491 0.9340
A2 =
0.4218 0.9157
0.7922 0.9595
0.6557 0.0357
0.8491 0.9340
B =
2.0000 + 3.0000i 4.0000 + 5.0000i 3.0000 + 0.0000i
2.0000 + 0.0000i 4.0000 + 1.0000i 5.0000 + 3.0000i
B1 =
2.0000 - 3.0000i 2.0000 + 0.0000i
4.0000 - 5.0000i 4.0000 - 1.0000i
3.0000 + 0.0000i 5.0000 - 3.0000i
B2 =
2.0000 + 3.0000i 2.0000 + 0.0000i
4.0000 + 5.0000i 4.0000 + 1.0000i
3.0000 + 0.0000i 5.0000 + 3.0000i
B3 =
2.0000 + 3.0000i 2.0000 + 0.0000i
4.0000 + 5.0000i 4.0000 + 1.0000i
3.0000 + 0.0000i 5.0000 + 3.0000i
B4 =
2.0000 - 3.0000i 2.0000 + 0.0000i
4.0000 - 5.0000i 4.0000 - 1.0000i
3.0000 + 0.0000i 5.0000 - 3.0000i
矩阵的旋转
旋转可以采用转置方式,也可以用函数rot90()。
B=rot90(A):该函数将矩阵逆时针旋转90°;
B=rot90(A,k):该函数将矩阵逆时针旋转90°的k倍,k默认是1.
矩阵的翻转
矩阵左右翻转是,第一列和最后一列互换,第二列和倒数第二列互换,依次类推。函数是fliplr(A)。
矩阵上下翻转是,第一行和最后一行互换,第二行和倒数第二行互换,依次类推。函数是flipud(A)。
此外,还可以采用flipdim()进行矩阵的翻转。flipdim(A,k)在指定的方向进行矩阵的翻转,k=1,相当于flipud(A);k=2,相当于fliplr(A)。
A=rand(3,4)
A1=fliplr(A)
A2=flipud(A)
A3=flipdim(A,1)
A4=flipdim(A,2)
运行结果:
A =
0.9597 0.2238 0.5060 0.9593
0.3404 0.7513 0.6991 0.5472
0.5853 0.2551 0.8909 0.1386
A1 =
0.9593 0.5060 0.2238 0.9597
0.5472 0.6991 0.7513 0.3404
0.1386 0.8909 0.2551 0.5853
A2 =
0.5853 0.2551 0.8909 0.1386
0.3404 0.7513 0.6991 0.5472
0.9597 0.2238 0.5060 0.9593
A3 =
0.5853 0.2551 0.8909 0.1386
0.3404 0.7513 0.6991 0.5472
0.9597 0.2238 0.5060 0.9593
A4 =
0.9593 0.5060 0.2238 0.9597
0.5472 0.6991 0.7513 0.3404
0.1386 0.8909 0.2551 0.5853
矩阵尺寸的改变
在矩阵总元素保持不变的情况下,用函数reshape()改变矩阵的尺寸。
Y=reshape(X,m,n),将矩阵转化为m行n列的二维矩阵。矩阵的总元素数不变。
A=[1:4;5:8]
A1=reshape(A,1,8)
A2=reshape(A1,[4,2])
A3=reshape(A,size(A2))
运行结果:
A =
1 2 3 4
5 6 7 8
A1 =
1 5 2 6 3 7 4 8
A2 =
1 3
5 7
2 4
6 8
A3 =
1 3
5 7
2 4
6 8
矩阵的除法
左除,矩阵A和B的左除为X=A\B,表示方程组A*X=B的解;
右除,矩阵A和B的右除为X=B/ A,表示方程组X *A=B的解。
注意,矩阵的除法和矩阵元素的除法不同,元素的除法是:
左除./
右除.
即矩阵的对应元素相除。
A=[21,2,3;7 3 1;9 4 2]
B=[3,5,7;2,12,4;2,7,4]
C1=A\B
C2=inv(A)*B
D1=B/A
D2=B*inv(A)
E1=A^3
E2=A*A*A
运行结果:
A =
21 2 3
7 3 1
9 4 2
B =
3 5 7
2 12 4
2 7 4
C1 =
0.2286 1.6286 0.5143
0.4286 4.4286 0.7143
-0.8857 -12.6857 -1.7429
C2 =
0.2286 1.6286 0.5143
0.4286 4.4286 0.7143
-0.8857 -12.6857 -1.7429
D1 =
-0.3429 -10.3714 9.2000
-1.4857 -1.9429 5.2000
-0.7714 -4.0857 5.2000
D2 =
-0.3429 -10.3714 9.2000
-1.4857 -1.9429 5.2000
-0.7714 -4.0857 5.2000
E1 =
11181 1428 1648
4140 539 610
5516 724 813
E2 =
11181 1428 1648
4140 539 610
5516 724 813
矩阵元素的查找
i=find(A):查找矩阵中的非0元素,返回单下标。
[i,j]=find(A):查找矩阵中的非0元素,返回双下标和。
矩阵元素的排序
Y=sort(X):升序排列。当X为向量时,返回由小到大排序后的向量;当X为矩阵时,返回X中各列由小到大排序后的矩阵。
Y=sort(X,DIM):在给定的维数上进行升序排列后的结果。当DIM=1时,按列升序排序;当DIM=2时,按行进行升序排序。
Y=sort(X,DIM,‘MODE’):可以指定排序方式。默认值为MODE=‘ascend’,即升序排序;MODE='descend’时,按降序排序。
矩阵元素的求和
在MATLAB中,求和函数有:sum()和cumsum()。
Y=sum(X):对矩阵X的元素求和,返回矩阵中各列元素的和组成的向量。
Y=sum(X,DIM):返回在给定的维数DIM上的元素的和,当DIM=1时,计算矩阵各列元素的和;DIM=2时,计算矩阵各行元素的和。
函数cumsum()的调用格式和sum()差不多,不同在于cumsum()返回值为矩阵。
矩阵元素的求积
在MATLAB中,求积函数有:prod()和cumprod()。
Y=prod(X):对矩阵X的元素求积,返回矩阵中各列元素的积组成的向量。
Y=prod(X,DIM):返回在给定的维数DIM上的元素的积,当DIM=1时,计算矩阵各列元素的积;DIM=2时,计算矩阵各行元素的积。
函数cumprod()的调用格式和prod()差不多,不同在于cumprod()返回值为矩阵。
A=[1,2,3,0;7 3 1,4;3 4 2,1]
B=prod(A)
C=prod(A,2)
D=cumprod(A)
E=cumprod(A)
运行结果:
A =
1 2 3 0
7 3 1 4
3 4 2 1
B =
21 24 6 0
C =
0
84
24
D =
1 2 3 0
7 6 3 0
21 24 6 0
E =
1 2 6 0
7 21 21 84
3 12 24 24
矩阵元素的差分
在MATLAB中,用函数diff()来计算矩阵的差分。
Y=diff(X):计算矩阵各列的差分。
Y=diff(X,N):计算矩阵各列的N阶差分。
Y=diff(X,N):计算矩阵在方向DIM上的N阶差分。当DIM=1时,计算矩阵各列元素的差分;当DIM=2时,计算矩阵各行元素的差分。
A=[1,2,3,0;7 3 1,4;3 4 2,1]
B=diff(A) %每一列之间的差分
C=diff(A,2) %二阶差分,在一阶差分的基础上再差分
D=diff(A,1,1) %对矩阵的每一列求差分
E=diff(A,1,2) %对矩阵的每一行求差分
运行结果:
A =
1 2 3 0
7 3 1 4
3 4 2 1
B =
6 1 -2 4
-4 1 1 -3
C =
-10 0 3 -7
D =
6 1 -2 4
-4 1 1 -3
E =
1 1 -3
-4 -2 3
1 -2 -1
特殊矩阵的生成
一、全零矩阵
A=zeros(N):产生N行N列的全零矩阵;
A=zeros(M,N):产生M行N列的全零矩阵;
A=zeros(size(B)):产生和矩阵B维数相同的全零矩阵。
二、全1矩阵
函数ones()产生全1矩阵,调用格式和zeros()基本相同。
三、单位矩阵
eye(),调用格式和zeros()基本相同。
A=eye(3)
B=eye([3,2])
C=eye(size(B))
运行结果:
A =
1 0 0
0 1 0
0 0 1
B =
1 0
0 1
0 0
C =
1 0
0 1
0 0
四、0-1之间均匀分布的随机矩阵
rand(),调用格式和zeros()基本相同。
若要产生0-n之间均匀分布的随机矩阵,就用n*rand()。
五、标准正态分布随机矩阵
randn()产生均值为0,方差为1的标准正态分布随机矩阵,调用和前面一样。
六、魔方矩阵
魔方矩阵中,每行、每列及两条对角线上的元素之和都相等。对于n阶魔方矩阵,其元素由1,2,3,…,n2组成,共n2个整数。
在MATLAB中,用magic()来产生魔方矩阵。
七、范德蒙矩阵
范德蒙(Vandermode)矩阵最后一列全为1,倒数第二列为一个指定的向量,其他各列是其后一列与倒数第二列的点乘积(平方、立方……)。可以用一个指定向量生成一个范德蒙矩阵。
在MATLAB中,用vander(V)来产生以向量V为基础向量的范德蒙矩阵。
A=vander([1 2 3 4])
B=vander([1;2;4;6])
运行结果:
A =
1 1 1 1
8 4 2 1
27 9 3 1
64 16 4 1
B =
1 1 1 1
8 4 2 1
64 16 4 1
216 36 6 1
八、希尔伯特矩阵
函数hilb(N)产生N阶的希尔伯特(Hilbert)矩阵。第i行j列的元素等于1/(i+j-1)。希尔伯特矩阵是一种病态矩阵,矩阵中任意一个元素发生微小的变化,整个矩阵的值和逆矩阵都会发生巨大的变化。
函数invhilb()产生希尔伯特矩阵的逆矩阵。
A=hilb(3)
B=hilb(4)
C=invhilb(3)
D=invhilb(4)
A*C
运行结果:
A =
1.0000 0.5000 0.3333
0.5000 0.3333 0.2500
0.3333 0.2500 0.2000
B =
1.0000 0.5000 0.3333 0.2500
0.5000 0.3333 0.2500 0.2000
0.3333 0.2500 0.2000 0.1667
0.2500 0.2000 0.1667 0.1429
C =
9 -36 30
-36 192 -180
30 -180 180
D =
16 -120 240 -140
-120 1200 -2700 1680
240 -2700 6480 -4200
-140 1680 -4200 2800
ans =
1.0000 0.0000 -0.0000
0 1.0000 0
0.0000 -0.0000 1.0000
九、托普利兹矩阵
托普利兹(Toeplitz)矩阵除第一行和第一列外,其他每个元素都与左上角的元素相同。
toeplitz(x):用向量x生成一个对称托普利兹矩阵。
toeplitz(x,y):产生一个以x为第一列,y为第一行的托普利兹矩阵。x,y均为向量,注意,向量x和y的第一个元素必须相等。
A=toeplitz(3:6)
x=3:8
y=3:7
B=toeplitz(x,y)
运行结果:
A =
3 4 5 6
4 3 4 5
5 4 3 4
6 5 4 3
x =
3 4 5 6 7 8
y =
3 4 5 6 7
B =
3 4 5 6 7
4 3 4 5 6
5 4 3 4 5
6 5 4 3 4
7 6 5 4 3
8 7 6 5 4
十、伴随矩阵
compan(p):p为多项式的系数向量,高次幂系数排在前,低次幂系数在后。
p1=[2 3 4]
A=compan(p1)
p2=[2 3 4 5]
B=compan(p2)
运行结果:
p1 =
2 3 4
A =
-1.5000 -2.0000
1.0000 0
p2 =
2 3 4 5
B =
-1.5000 -2.0000 -2.5000
1.0000 0 0
0 1.0000 0
十一、帕斯卡矩阵
二次项展开后的系数随n的增大组成一个三角形表,称为杨辉三角形。由杨辉三角形表组成的矩阵称为帕斯卡(Pascal)矩阵。
pascal(n):生成一个n阶的帕斯卡矩阵。
A=pascal(4)
B=pascal(5)
运行结果:
A =
1 1 1 1
1 2 3 4
1 3 6 10
1 4 10 20
B =
1 1 1 1 1
1 2 3 4 5
1 3 6 10 15
1 4 10 20 35
1 5 15 35 70
方阵的行列式
把一个方针看作一个行列式,并对其按照行列式的规则求值,这个值就是矩阵所对应的行列式的值。
det(X):求方阵X的行列式的值。
A=magic(3)
B=[1:3;2 5 7;3 8 7]
y1=det(A)
y2=det(B)
运行结果:
A =
8 1 6
3 5 7
4 9 2
B =
1 2 3
2 5 7
3 8 7
y1 =
-360
y2 =
-4
特征值、特征向量和特征多项式
E=eig(A):求矩阵A的全部特征值,组成向量E。
[V,D]=eig(A):计算矩阵的特征值和特征向量,返回值是V和D两个方阵,方阵V每一列为一个特征向量,方阵D为对角矩阵,对角线上元素为特征值。
A=magic(3)
B=eig(A)
[V,D]=eig(A)
运行结果:
A =
8 1 6
3 5 7
4 9 2
B =
15.0000
4.8990
-4.8990
V =
-0.5774 -0.8131 -0.3416
-0.5774 0.4714 -0.4714
-0.5774 0.3416 0.8131
D =
15.0000 0 0
0 4.8990 0
0 0 -4.8990
p=[3 5 2 1] %3 5 2 1分别是多项式的系数
A=compan(p) %伴随矩阵
x1=eig(A) %伴随矩阵的特征值就是多项式的根
x2=roots(p) %roots直接求得多项式的根
运行结果:
p =
3 5 2 1
A =
-1.6667 -0.6667 -0.3333
1.0000 0 0
0 1.0000 0
x1 =
-1.3563 + 0.0000i
-0.1552 + 0.4708i
-0.1552 - 0.4708i
x2 =
-1.3563 + 0.0000i
-0.1552 + 0.4708i
-0.1552 - 0.4708i
对角阵
只有对角线上有非0元素的称为对角阵,对角线上的元素相等的对角矩阵是数量矩阵,对角线上的元素都为1的对角矩阵称为单位矩阵。
diag(A):提取矩阵A的主对角线元素,产生一个列向量。
diag(A,k):提取第k条对角线上的元素,组成一个列向量。
A=rand(3,4)
b1=diag(A) %主对角线也是第0条对角线
b2=diag(A,1)
b3=diag(A,2)
运行结果:
A =
0.6787 0.3922 0.7060 0.0462
0.7577 0.6555 0.0318 0.0971
0.7431 0.1712 0.2769 0.8235
b1 =
0.6787
0.6555
0.2769
b2 =
0.3922
0.0318
0.8235
b3 =
0.7060
0.0971
上三角阵和下三角阵
上三角阵:矩阵的对角线以下的元素全为0;
下三角阵:矩阵的对角线以上的元素全为0。
triu(A):返回矩阵A的上三角矩阵;
triu(A,k):返回矩阵A的第k条对角线以上的元素。
tril()求矩阵的下三角矩阵,用法和triu()一样。
A=rand(3,3)
b1=triu(A)
b2=tril(A)
运行结果:
A =
0.6948 0.0344 0.7655
0.3171 0.4387 0.7952
0.9502 0.3816 0.1869
b1 =
0.6948 0.0344 0.7655
0 0.4387 0.7952
0 0 0.1869
b2 =
0.6948 0 0
0.3171 0.4387 0
0.9502 0.3816 0.1869
矩阵的逆和伪逆
对于方阵A,如果存在一个与其同阶的方阵B,使得AB=BA=E,则称A和B互为逆矩阵。采用函数inv()求方阵的逆矩阵。
如果矩阵A不是一个方阵,或者是一个非满秩方阵时,没有逆矩阵,但可以找到一个与A的转置矩阵同型的矩阵B,使得A×B×A=A,B×A×B=B。此时,称矩阵B为矩阵A的伪逆,也称广义逆矩阵。求广义逆矩阵的函数是pinv()。
A=magic(3)
B=[1 3;2 6]
C=inv(A)
C*A
inv(B)
D=pinv(B)
B*D*B
运行结果:
A =
8 1 6
3 5 7
4 9 2
B =
1 3
2 6
C =
0.1472 -0.1444 0.0639
-0.0611 0.0222 0.1056
-0.0194 0.1889 -0.1028
ans =
1.0000 -0.0000 -0.0000
0 1.0000 0
0 0.0000 1.0000
警告: 矩阵为奇异工作精度。
> In test_one (line 5)
ans =
Inf Inf
Inf Inf
D =
0.0200 0.0400
0.0600 0.1200
ans =
1.0000 3.0000
2.0000 6.0000
矩阵的秩
矩阵的秩包括行秩和列秩,行秩和列秩相等,行秩为矩阵的行向量组成的极大无关组中行向量的个数,列秩为矩阵的行向量组成的极大无关组中列向量的个数。对于满秩矩阵,秩等于行数或列数,其各行向量或列向量都线性无关。
函数rank()求矩阵的秩。
矩阵的迹
矩阵的迹等于矩阵的对角线元素之和,也等于矩阵的特征值之和。
函数trace()求矩阵的迹。
矩阵的范数
求矩阵范数的函数是:norm(),该函数的调用格式是:
norm(X)或norm(X,2):计算矩阵的2-范数,返回矩阵的最大奇异值。
norm(X,1):计算矩阵的1-范数,返回矩阵的列向元素和的最大值。
norm(X,inf):计算矩阵的 ∞ -范数,返回矩阵的行向元素和的最大值。
norm(X,‘fro’):计算矩阵的Frobenius范数。
A=[1 2 3;3 5 7;2 5 8]
n1=norm(A,1)
n2=norm(A)
n3=norm(A,inf)
n4=norm(A,'fro')
n5=normest(A)
运行结果:
A =
1 2 3
3 5 7
2 5 8
n1 =
18
n2 =
13.7529
n3 =
15
n4 =
13.7840
n5 =
13.7529
矩阵的条件数
矩阵的条件数是用来判断矩阵病态的一个量,矩阵的条件数越大,表明该矩阵越病态,否则该矩阵越良态。Hilbert矩阵就是典型的病态矩阵。
通过函数cond()来求矩阵的条件数,调用格式为:
cond(X,1):计算矩阵X的1-范数下的条件数。
cond(X)或cond(X,2):计算矩阵X的2-范数下的条件数。
cond(X,inf):计算矩阵X的 ∞ -范数下的条件数。
A=[1 2 3;3 5 7;2 5 8]
n1=cond(A,1)
n2=cond(A)
n3=cond(A,inf)
运行结果:
A =
1 2 3
3 5 7
2 5 8
n1 =
9.7278e+18
n2 =
3.1503e+16
n3 =
5.6745e+18
矩阵的标准正交基
B=orth(A):矩阵B的列向量构成了矩阵A的一组标准正交基。
A=[1 2 3;3 5 7;2 5 8]
B=magic(3)
C=orth(A)
D=orth(B)
C'*C
eye(rank(A))
运行结果:
A =
1 2 3
3 5 7
2 5 8
B =
8 1 6
3 5 7
4 9 2
C =
-0.2721 -0.0075
-0.6606 -0.7256
-0.6997 0.6881
D =
-0.5774 0.7071 0.4082
-0.5774 0.0000 -0.8165
-0.5774 -0.7071 0.4082
ans =
1.0000 0.0000
0.0000 1.0000
ans =
1 0
0 1
矩阵的超越函数
sqrtm():计算矩阵的平方根。
logm():计算矩阵的自然对数。
expm():计算矩阵的指数。
funm():计算矩阵的超越函数值。
A=[1 2 3;3 5 7;2 5 8]
X1=expm(A)
Y1=logm(X1)
X2=logm(A)
Y2=expm(X2)
[V,D]=eig(A)
V*diag(exp(diag(D)))/V
运行结果:
A =
1 2 3
3 5 7
2 5 8
X1 =
1.0e+05 *
0.5456 1.1347 1.7237
1.3327 2.7715 4.2103
1.3953 2.9018 4.4083
Y1 =
1.0000 2.0000 3.0000
3.0000 5.0000 7.0000
2.0000 5.0000 8.0000
X2 =
-29.8095 6.3498 6.6107
59.0623 -11.2819 -9.8294
-28.6179 7.1326 6.9848
Y2 =
1.0000 2.0000 3.0000
3.0000 5.0000 7.0000
2.0000 5.0000 8.0000
V =
-0.2721 -0.4082 -0.0531
-0.6646 0.8165 -0.8263
-0.6959 -0.4082 0.5607
D =
13.5574 0 0
0 0.0000 0
0 0 0.4426
ans =
1.0e+05 *
0.5456 1.1347 1.7237
1.3327 2.7715 4.2103
1.3953 2.9018 4.4083
A=[1 2 3;3 5 7;2 5 8]
X1=funm(A,@sin)
X2=funm(A,@cos)
X3=funm(A,@exp)
Y=expm(A)
运行结果:
A =
1 2 3
3 5 7
2 5 8
X1 =
0.1002 0.1314 0.1625
0.7842 0.4320 0.0798
-0.2831 0.2248 0.7326
X2 =
0.9588 -0.0683 -0.0955
-0.2220 0.8081 -0.1618
0.0160 -0.1497 0.6845
X3 =
1.0e+05 *
0.5456 1.1347 1.7237
1.3327 2.7715 4.2103
1.3953 2.9018 4.4083
Y =
1.0e+05 *
0.5456 1.1347 1.7237
1.3327 2.7715 4.2103
1.3953 2.9018 4.4083
稀疏矩阵
在MATLAB中,对于矩阵的存储有两种方式:完全存储方式和稀疏存储方式。完全存储方式是将矩阵的全部元素按列进行存储。如果矩阵中的元素只有少数不是0,这种方式会浪费大量存储空间。
- 矩阵存储方式
稀疏矩阵只是存储方式不同,运算规则和普通矩阵是一样的。用户可以创建整型、双精度、复数类型和逻辑类型的稀疏矩阵。稀疏矩阵不能自动生成。定义在完全存储方式下的运算只能产生完全存储的矩阵,不论结果有多少个0。稀疏矩阵参与运算的结果还是稀疏矩阵。 - 产生稀疏矩阵
普通矩阵→稀疏矩阵
S=sparse(A):将矩阵A转换为稀疏矩阵S。当矩阵A是稀疏存储方式时,相当于调用S=A。
S=sparse(m,n):产生m行n列的元素全为0的稀疏矩阵。
稀疏矩阵→普通矩阵
B=full(A):把稀疏矩阵A转化为普通矩阵B
A=rand(15,10)>0.98
S=sparse(A)
whos
运行结果:
A =
15×10 logical 数组
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 1 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
S =
15×10 稀疏 logical 数组
矩阵中只有一个元素为1,位置13行8列
(13,8) 1
Name Size Bytes Class Attributes
A 15x10 150 logical
B 3x3 72 double
C 3x2 48 double
D 3x3 72 double
S 15x10 97 logical sparse
V 3x3 72 double
X1 3x3 72 double
X2 3x3 72 double
X3 3x3 72 double
Y 3x3 72 double
Y1 3x3 72 double
Y2 3x3 72 double
ans 3x3 72 double
n1 1x1 8 double
n2 1x1 8 double
n3 1x1 8 double
n4 1x1 8 double
n5 1x1 8 double
A=[0 0 0 2;0 0 3 0;0 0 0 0;4 0 0 0]
S1=sparse(A)
S2=sparse([4,2,1],[1,3,4],[4,3,2],4,4) %稀疏矩阵是4行4列的,有三个元素
B=full(S1)
运行结果:
A =
0 0 0 2
0 0 3 0
0 0 0 0
4 0 0 0
S1 =
(4,1) 4
(2,3) 3
(1,4) 2
S2 =
(4,1) 4
(2,3) 3
(1,4) 2
B =
0 0 0 2
0 0 3 0
0 0 0 0
4 0 0 0
函数nnz(S)计算稀疏矩阵S中非零值的个数。
函数spy()对矩阵中非零元素的分布进行图形化显示,获取稀疏矩阵中非零元素的密度(非零元素的个数/(稀疏矩阵的行数*列数))。
函数spalloc()为稀疏矩阵分配空间。
A=[0 0 0 2;0 0 3 0;0 0 0 0;4 0 0 0]
S=sparse(A)
n1=nnz(S)
n2=nonzeros(S)
n3=nzmax(S)
运行结果:
A =
0 0 0 2
0 0 3 0
0 0 0 0
4 0 0 0
S =
(4,1) 4
(2,3) 3
(1,4) 2
n1 =
3
n2 =
4
3
2
n3 =
3
v=[1 2 3 4 5]
X=diag(v,1)
S=sparse(X)
nnz(S)
d=nnz(S)/prod(size(S)) %计算稀疏矩阵中非零元素的密度
运行结果:
v =
1 2 3 4 5
X =
0 1 0 0 0 0
0 0 2 0 0 0
0 0 0 3 0 0
0 0 0 0 4 0
0 0 0 0 0 5
0 0 0 0 0 0
S =
(1,2) 1
(2,3) 2
(3,4) 3
(4,5) 4
(5,6) 5
ans =
5
d =
0.1389
- 特殊稀疏矩阵
单位矩阵只有对角线上元素是1,其余元素都为0,是一种具有稀疏特征的矩阵。产生稀疏存储方式的单位矩阵的函数是:speye()。
S=speye(n):产生一个n行n列的单位稀疏存储矩阵。
S=speye(m,n):产生一个m行n列的单位稀疏存储矩阵。
spones():将稀疏矩阵中的非零元素替换为1。
spconvert():将普通矩阵转化为稀疏矩阵。
A=sprandsym(8,0.1)
B=spones(A)
运行结果:
A =
(1,1) -1.1201
(2,1) 2.5260
(1,2) 2.5260
(3,2) 1.6555
(2,3) 1.6555
(7,6) 0.3075
(6,7) 0.3075
B =
(1,1) 1
(2,1) 1
(1,2) 1
(3,2) 1
(2,3) 1
(7,6) 1
(6,7) 1
A=[2 3 4 5;5 6 7 2;1 2 4 0;8 7 9 2]
S=spconvert(A)
运行结果:
A =
2 3 4 5
5 6 7 2
1 2 4 0
8 7 9 2
S =
(1,2) 4.0000 + 0.0000i
(2,3) 4.0000 + 5.0000i
(5,6) 7.0000 + 2.0000i
(8,7) 9.0000 + 2.0000i
矩阵的分解
矩阵分解是指根据一定的原理用某种算法将一个矩阵分解成若干个矩阵的乘积。下面是几种分解方式:
- Cholesky分解
对于正定矩阵,可以分解为上三角矩阵和下三角矩阵的乘积,这种分解称为Cholesky分解。能进行Cholesky分解的必须是正定矩阵,矩阵的所有对角元素必须是正的,同时矩阵的非对角元素不能太大。
通过函数chol()进行矩阵的Cholesky分解。在进行Cholesky分解之前,最好先用eig()查看矩阵的特征值,检验特征值是否为正。
A=pascal(4)
eig(A)
R=chol(A)
R'*R %验证
运行结果:
A =
1 1 1 1
1 2 3 4
1 3 6 10
1 4 10 20
ans =
0.0380
0.4538
2.2034
26.3047
R =
1 1 1 1
0 1 2 3
0 0 1 3
0 0 0 1
ans =
1 1 1 1
1 2 3 4
1 3 6 10
1 4 10 20
- LU分解
高斯消去法又称LU分解,将方阵A分解为下三角矩阵的置换矩阵L和上三角矩阵U的乘积。
通过函数lu()进行矩阵的LU分解,调用格式为:
[L1,U1]=lu(A):将方阵A分解为下三角矩阵的置换矩阵L1和上三角矩阵U1。
[L2,U2,P]=lu(A):将方阵A分解为下三角矩阵的置换矩阵L1和上三角矩阵U1以及置换矩阵P。
Y=lu(A):将下三角矩阵和上三角矩阵合并在矩阵Y中,矩阵Y的对角元素为上三角矩阵的对角元素。
A=[2 3 4;8 4 9;5 3 1]
[L1,U1]=lu(A)
[L2,U2,P]=lu(A)
Y1=lu(A)
L1*U1
Y2=L2+U2-eye(size(A))
运行结果:
A =
2 3 4
8 4 9
5 3 1
L1 =
0.2500 1.0000 0
1.0000 0 0
0.6250 0.2500 1.0000
U1 =
8.0000 4.0000 9.0000
0 2.0000 1.7500
0 0 -5.0625
L2 =
1.0000 0 0
0.2500 1.0000 0
0.6250 0.2500 1.0000
U2 =
8.0000 4.0000 9.0000
0 2.0000 1.7500
0 0 -5.0625
P =
0 1 0
1 0 0
0 0 1
Y1 =
8.0000 4.0000 9.0000
0.2500 2.0000 1.7500
0.6250 0.2500 -5.0625
ans =
2 3 4
8 4 9
5 3 1
Y2 =
8.0000 4.0000 9.0000
0.2500 2.0000 1.7500
0.6250 0.2500 -5.0625
- QR分解
矩阵的正交分解,就是QR分解。QR分解将一个m行n列的矩阵A分解为一个正交矩阵Q(m行m列)和一个上三角矩阵R(m行n列)的乘积。
通过函数qr()进行QR分解,调用格式为[Q,R]=qr(A)。
A=[2 3 4;8 4 9;5 3 1]
[Q1,R1]=qr(A)
Q1*R1
运行结果:
A =
2 3 4
8 4 9
5 3 1
Q1 =
-0.2074 0.9548 -0.2129
-0.8296 -0.2870 -0.4790
-0.5185 0.0773 0.8516
R1 =
-9.6437 -5.4958 -8.8141
0 1.9483 1.3136
0 0 -4.3112
ans =
2.0000 3.0000 4.0000
8.0000 4.0000 9.0000
5.0000 3.0000 1.0000
- SVD分解
即奇异值分解,通过函数svd()进行SVD分解(奇异值分解)。
s=svd(A):对矩阵A进行奇异值分解,返回由奇异值组成的列向量,奇异值按从大到小的顺序排列。
[U,S,V]=svd(A):U和V为酉矩阵,S为一个对角矩阵,对角线元素为矩阵的奇异值降序排列。
A=[2 3 4;8 4 9;5 3 1]
s=svd(A)
[U,S,V]=svd(A)
U*S*V'
norm(A)
运行结果:
A =
2 3 4
8 4 9
5 3 1
s =
14.5279
3.3393
1.6696
U =
-0.3497 0.2968 0.8886
-0.8701 0.2486 -0.4255
-0.3472 -0.9220 0.1713
S =
14.5279 0 0
0 3.3393 0
0 0 1.6696
V =
-0.6468 -0.6073 -0.4614
-0.3835 -0.2640 0.8850
-0.6592 0.7494 -0.0622
ans =
2.0000 3.0000 4.0000
8.0000 4.0000 9.0000
5.0000 3.0000 1.0000
ans =
14.5279
- Schur分解
对矩阵的Schur分解公式为A=USU’,矩阵A必须为方阵,U为酉矩阵,S为块对角矩阵。
通过函数schur()进行矩阵的Schur分解,调用格式是:
[U,S]=schur(A):返回酉矩阵U,块对角矩阵S;
S=schur(A):仅返回块对角矩阵S。
A=pascal(4)
[U,S]=schur(A)
运行结果:
A =
1 1 1 1
1 2 3 4
1 3 6 10
1 4 10 20
U =
0.3087 -0.7873 0.5304 0.0602
-0.7231 0.1632 0.6403 0.2012
0.5946 0.5321 0.3918 0.4581
-0.1684 -0.2654 -0.3939 0.8638
S =
0.0380 0 0 0
0 0.4538 0 0
0 0 2.2034 0
0 0 0 26.3047
- Hessenberg分解
对于任意一个n阶方阵可以进行Hessenberg分解,分解公式为:A=PHP’,其实P是酉矩阵,H的第一子对角线下的元素均为0,即H为Hessenberg矩阵。通过函数hess()进行Hessenberg分解,函数调用格式是:
H=hess(A):该函数对方阵A进行Hessenberg分解,返回Hessenberg矩阵。
[P,H]=hess(A):该函数对方阵A进行Hessenberg分解,返回值为P和H,满足A=PHP’。
A=[1 3 5 3;2 6 9 1;7 4 1 10;2 5 9 3]
H1=hess(A)
[P2,H2]=hess(A)
B=P2*H2*P2'
C=P2'*P2
运行结果:
A =
1 3 5 3
2 6 9 1
7 4 1 10
2 5 9 3
H1 =
1.0000 -6.2253 2.0602 -0.0327
-7.5498 9.7719 -9.0220 -2.4797
0 -11.8114 -1.3405 -4.7895
0 0 0.5086 1.5686
P2 =
1.0000 0 0 0
0 -0.2649 0.6443 -0.7174
0 -0.9272 -0.3746 0.0059
0 -0.2649 0.6667 0.6966
H2 =
1.0000 -6.2253 2.0602 -0.0327
-7.5498 9.7719 -9.0220 -2.4797
0 -11.8114 -1.3405 -4.7895
0 0 0.5086 1.5686
B =
1.0000 3.0000 5.0000 3.0000
2.0000 6.0000 9.0000 1.0000
7.0000 4.0000 1.0000 10.0000
2.0000 5.0000 9.0000 3.0000
C =
1.0000 0 0 0
0 1.0000 -0.0000 -0.0000
0 -0.0000 1.0000 -0.0000
0 -0.0000 -0.0000 1.0000