Matlab笔记——数值计算—线代篇016

1、16. 数值计算线代篇一、行列式det(A)矩阵A的行列式； inv(A)矩阵A的逆；rank(A)矩阵A的秩；B(: , i)=b将向量b赋给矩阵B的第i行；A, eye(5)在矩阵A右端，拼接5阶单位矩阵；U,s=rref(A)对矩阵A作行变换，U返回A的最简行阶梯形矩阵，s为行向量存储U的各行首个非0元所在列号，length(s)即为A的秩；例1 用初等行变换法求矩阵的逆。代码：format short g % 省略小数位多余的0A=1 2 3; 2 2 1; 3 4 3;B=rref(A,eye(3) % 对矩阵A,I进行初等行变换,得到最简行阶梯矩阵Bif(rank(B(:,1:3)

2、=3) % 判断B的前3列是否为单位阵，若是取出后3列，即A逆 A1=B(:,4:6) elsedisp('A不可逆');end运行结果：B =1 0 0 1 3 -20 1 0 -1.5 -3 2.50 0 1 1 1 -1A1 = 1 3 -2 -1.5 -3 2.51 1 -1例2解方程代码：syms x;A=3 2 1 1;3 2 2-x2 1;5 1 3 2;7-x2 1 3 2;D=det(A) f=factor(D) % 对行列式D进行因式分解X=solve(D) % 求方程“D0”的解运行结果：D =-3*(x2 - 1)*(x2 - 2) f =-3*(x -

3、 1)*(x + 1)*(x2 - 2) X = -1 1 2(1/2)-2(1/2)二、向量组的线性相关性例3 向量组求它的秩和一个最大线性无关组，并用来表示其它向量。代码：A=2 -1 3 5;4 -3 1 3;3 -2 3 4; 4 -1 15 17;7 -6 -7 0 '% format rat; % 使用分数表示rref(A)运行结果：ans = 10021 010-35 0 014 -5 0 00 0 0可见，向量组的秩是3，是一个最大线性无关组；并且，注：也可以用R,s=rref(A); length(s)得到秩。三、线性方程组的通解null(A, r)返回齐次线性方程组

4、Ax=0的基础解系，选项r返回有理数解，否则按分数显示；x0=inv(A)*b若A-1存在，直接可以得到Ax=b的一个特解x0,否则只能按求解理论求解；subs(A, k, n), 将矩阵或式A中的k用n代替。例4求下列方程组的通解：代码：syms x;A=2,4,-1,4,16;-3,-6,2,-6,-23;3,6,-4,6,19;1,2,5,2,19;b=-2;7;-23;43;R,s=rref(A,b)m,n=size(A);x0=zeros(n,1); % 将特解x0初始化为零向量r=length(s); % 矩阵A的秩赋给变量rx0(s,:)=R(1:r,end) %矩阵R的最后一列

5、按基准元素的位置得到特解x0x=null(A,'r') % 得到对应齐次方程组Ax0的基础解系运行结果：R = 1202930010280000000 00 0 0 0s = 1 3x0 = 3 0 8 0 0x = -2 -2 -9 1 0 0 0 0 -2 0 1 0 0 0 1例5已知齐次线性方程组问当k取何值时方程组有非零解？在有非零解的情况下，求出其基础解系。代码：syms k;A=1-2*k,3,3,3;3,2-k,3,3;3,3,2-k,3;3,3,3,11-k;D=det(A)K=solve(D) % 解方程“D0”即要求的k值for i=1:4 Ak=doub

6、le(subs(A,k,K(i); % 用K(i)替换矩阵A中的k,得到的是符号矩阵，% 再用double函数转化为数值矩阵x=null(Ak,'r') end运行结果：D =2*k4 - 31*k3 + 30*k2 + 161*k + 98 K = -1-17/214 x = -1 -1 1 0 0 1 0 0x = -1 -1 1 0 0 1 0 0x = -0.5 -1 -1 1x = 0.2 0.4 0.4 1四、基变换设Rn中的两组基向量U和V（都是n×n矩阵），若向量w在以U为基的坐标系内的坐标为wu（n×1数组），在以V为基的坐标系内的坐标为w

7、v（n×1数组），则在基准坐标系内的坐标应分别为U*wu和V*wv，这两者应该相等，即U*wu=V*wv所谓基坐标的变换，就是已知wu，求出wv. 将上式两边均左乘V-1，得到wv=V-1*U*wu故坐标变换矩阵为P=V-1*U.例6 已知R4中的两组基向量为，求从U到V的坐标变换矩阵P.代码：U=1 1 -1 -1; 2 -1 2 -1; -1 1 1 0; 0 1 1 1;V=2 0 -2 1; 1 1 1 3; 0 2 1 1; 1 2 2 2;P=inv(V)*U % 从U到V的基变换矩阵wu=1 2 3 4' wv=P*wu% 已知某向量在U坐标系下坐标为wu，求它

8、在V坐标系下的坐标运行结果：P = 0 1 -1 1 -1 1 0 00 0 0 11 -1 1 -1wv = 314-2五、特征值与特征向量变换可以用矩阵表示。矩阵A的特征向量是指，经过矩阵A（左乘）变换后不发生方向改变的那些向量；特征值是指在经过这些变换后特征向量的伸缩的倍数。对于实对称矩阵来说，不同特征值对应的特征向量必定正交。矩阵（变换）A的所有特征向量组成了该矩阵（变换）的一组基，可以理解为坐标系的坐标轴，可以把这个坐标系扭曲、拉伸、旋转，称为基变换。例如，在主成分分析中，通过在拉伸最大的方向设置基，忽略一些小的量，可以极大的压缩数据而减小失真。矩阵（变换）的所有特征向量作为空间的基

9、之所以重要，是因为在这些方向上矩阵（变换）可以拉伸向量而不必扭曲和选择它，使得计算大为简单。2. Matlab实现orth(A)返回矩阵A的列向量组构成空间的标准正交基；P=poly(A)返回矩阵A的特征多项式，P是行向量，元素为多项式系数；roots(P)求多项式P的零点；r=eig(A) r为列向量，元素为矩阵A的特征值；V, D=eig(A)矩阵D为A的特征值（从小到大排列）所构成的对角矩阵，V的列向量为A的特征向量，与D中特征值一一对应（VTAV=D，VT=V-1），D也是矩阵A的相似对角化矩阵；V, D=schur(A)矩阵D为对称矩阵A的特征值所构成的对角阵，V的列为A的单位特征向

10、量，与D中特征值一一对应；例7求下列矩阵的特征值和特征向量：代码：format short gA=1 2 3; 2 1 3; 1 1 2;% 特征多项式法% P=poly(A); % lamda1=roots(P) lamda2=eig(A)V,D=eig(A)inv(V)*A*V运行结果：lamda2 = 5.0000-1.0000-0.0000V = 0.6396 0.7071 -0.5774 0.6396 -0.7071 -0.5774 0.4264 -0.0000 0.5774D = 5.0000 0 0 0 -1.0000 0 0 0 -0.0000ans = 5.0000 -0.0

11、000 0-0.0000 -1.0000 -0.0000 -0.0000 0.0000 0例8（人口迁徙模型）设某大城市的总人口是固定的。人口的分布则因居民在市区和郊区之间迁徙而变化。每年有6%的市区居民搬到郊区去住，而有2%的郊区居民搬到市区。假如开始时有30%的居民住在市区，70%的居民住在郊区，问：(1) 10年后市区和郊区的居民人口比例是多少？30年、50年后又如何？(2) 无限增加时间，该比例最终是否会趋于稳定？问题分析：该问题可以用矩阵乘法来描述。(1) 把人口变量用市区xc和郊区xs两个分量表示，一年以后，市区人口为xc1= (1-0.06) xc0+0.02xs0，郊区人口xs

12、1=0.06xc0+ (1-0.02)xs0，用矩阵乘法来表示：从初始到k年，此关系保持不变，故对于(2), 借助矩阵A的特征值和特征向量，因为A作用在其特征向量上只改变值大小（为其特征值的倍数）而不改变方向。为此先求出A的特征值1, 2及特征向量v1, v2; 将x0用特征向量v1, v2作为坐标系表示出来（即做以v1,v2为基坐标的基变换），x0=v1+v2从而，令k趋于，判断xk是否存在极限。代码：x0=0.3; 0.7;A=0.94 0.02; 0.06 0.98;x1=A*x0x10=A10*x0x50=A50*x0V,D=eig(A)V(:,1)=V(:,1)./V(2,1);V(

13、:,2)=V(:,2)./V(1,2) % 对V做初等变换化简k=inv(V)*x0 % 求x0用特征向量作为坐标系的表示% x0在原坐标系I中的表示为x0,% 则基坐标变换到在坐标系V的表示为 x=P*x0=inv(V)*I*x0syms n;xn=k(1)*D(1,1)n*V(:,1)+k(2)*D(2,2)n*V(:,2);% xn=-0.05*(0.92)n*v1+0.25*1n*v2limit(xn,n,inf)运行结果：x1 = 0.29600.7040x10 = 0.27170.7283x50 = 0.25080.7492V = -0.7071 -0.31620.7071 -0.

14、9487D = 0.92000 0 1.0000V = -1.0000 1.0000 1.0000 3.0000k = -0.0500 0.2500ans = 1/43/4最终得到：无限增加时间k，市区和郊区人口之比将趋向一组常数0.25/0.75.六、二次型例9用正交变换法将下列二次型化为标准型（保持几何形状不变）：代码：A=1 -2 0; -2 2 -2; 0,-2,3; % 输入二次型的矩阵AV,D=eig(A) % 或用V,D=schur(A), 结果相同; V为正交变换矩阵，即Y=VXsyms y1 y2 y3; f=y1,y2,y3*D*y1;y2;y3运行结果：V = -0.66

15、67 -0.6667 0.3333-0.6667 0.3333 -0.6667 -0.3333 0.6667 0.6667D = -1.0000 0 0 0 2.0000 0 0 0 5.0000f =- y12 + 2*y22 + 5*y32注：XTAX=YTVTAVY=YTDY七、用矩阵做图形变换图形变换是指对图形进行平移、旋转、缩放、投影（透视）等变换，其实质是改变图形的各个顶点的坐标。图形变换可以通过对表示图形坐标的矩阵进行运算来实现，称为矩阵变换法：注：也可以用变换矩阵左乘实现，上式两边转置即可，即：BTAT=CT.下面只介绍二维图形的矩阵变换，代码示例是类似的省略之。1. 基本变换

16、(1) 缩放变换若a=d则为等比例缩放；若a=1，则x方向不变；若a=0, 则图形压缩为y轴上的线段。例10 已知三角形的三个顶点坐标，绘制三角形并进行缩放变换。代码：ABC=4 4; 1 3; 3 1; % 三角形的三个顶点E=ABC(1,1),ABC(1,2); % 为了让图形封闭，在末尾补上起始点ABC=ABC;Esubplot(1,2,1)plot(ABC(:,1),ABC(:,2),'r'); % 将各个顶点连线，绘制三角形axis(0,5,0,5);% 放缩变换a=0.6;d=1.2;T=a 0; 0 d; % 变换矩阵ABC1=ABC*T; % 对原顶点做变换得到

17、新的顶点subplot(1,2,2)plot(ABC1(:,1),ABC1(:,2),'g'); % 绘制变换后的三角形axis(0,5,0,5);运行结果：(2) 对称变换关于原点的对称变换关于x轴的对称变换关于y轴的对称变换关于直线y=x的对称变换关于直线y=-x的对称变换(3) 错切变换（延x轴或y轴一个方向移动，另一个方向不变）延x轴方向错切变换后，平行于x轴的直线变换后仍平行于x轴；平行于y轴的直线变换后，y=0的点不动(不动点)，y0的点沿x方向平移了cy，形成与y轴夹角为的直线，且tgcy / yc.延y轴方向错切(4) 绕坐标原点的旋转变换变换矩阵为逆时针方向旋

18、转时角度取正值；顺时针方向旋转时角度取负值。2. 平移变换2×2变换矩阵是不能实现平移变换的，为此需要扩展一维。注：补上一维的常数，称为齐次坐标表示法；若常数=1，称为标准齐次坐标表示法。若坐标变换结果是非标准化齐次坐标表示，应将其化为标准齐次坐标表示，方法是所有项都除以齐次项：例11已知三角形的三个顶点坐标，绘制三角形并进行平移变换。代码：ABC=3 4; 1 3; 3 1; % 三角形的三个顶点E=ABC(1,1),ABC(1,2); % 为了让图形封闭，在末尾补上起始点ABC=ABC;Esubplot(1,2,1)plot(ABC(:,1),ABC(:,2),'r'); % 将各个顶点连线，绘制三角形axis(0,5,0,5);% 平移变换ABC1=ABC,ones(length(ABC),1); % 加上一列1扩维, 变成标准齐次坐标表示l=1;m=0.5;T=1 0 0; 0 1 0; l