线性代数实践(教师班第9讲).ppt

《线性代数实践(教师班第9讲).ppt》由会员分享，可在线阅读，更多相关《线性代数实践(教师班第9讲).ppt（63页珍藏版）》请在三一办公上搜索。

1、第9章 线性变换及其特征,9.1 平面上线性变换的几何意义,例9.1 设x为二维平面上第一象限中的一个单位方块，其四个顶点的数据可写成把不同的A矩阵作用于此组数据，可以得到多种多样的结果yi=Ai*x。用程序ag911进行变换计算，并画出x及yi图形：x0,1,1,0;0,0,1,1;subplot(2,3,1),fill(x(1,:),0,x(2,:),0,r)A11,0;0,1,y1A1*xsubplot(2,3,2),fill(y1(1,:),0,y1(2,:),0,g),几种变换的行列式与特征值,看出的基本关系,可以看出，矩阵A1使原图对纵轴生成镜像，矩阵A2使原图在横轴方向膨胀，矩阵

2、A3使原图在纵轴方向压缩，矩阵A4使原图向右方剪切变形，矩阵A5使原图沿反时针方向旋转tpi/6。分别计算出这五个矩阵的行列式和特征值，对二维空间（平面），行列式的几何意义实际上是两个向量所构成的平行四边形的面积。一个变换所造成的图形的面积变化，取决于该变换的行列式。A1,A4和A5的行列式绝对值都是1，所以它们不会使变换后图形的面积发生改变。而A2和A3的行列式分别为1.5和0.2，,9.2 二维矩阵特征值的几何意义,二维矩阵的特征值表示该变换在原图形的特征向量的方向上的放大量。例如矩阵A1在第一特征向 量 方向的特征值为，即横轴 正方向的增益为1，其结果是把原图中横轴正方向的部分变换到新图

3、的负方向去了；A1在第二特 征向量 的方向的特征值为1(2)=1，即纵轴正方向的增益为1，因而保持了新图和原图在纵轴方向尺度不变。,用eigshow函数看特征值,对于比较复杂的情况，完全凭简单的几何关系去想像是困难的，应当用eigshow函数，联系x和Ax的向量图来思考。键入eigshow(A4)。绿色的x表示原坐标系中的单位向量，可以用鼠标左键点住x并拖动它围绕原点转动。图中同时出现以蓝色表示的Ax向量，它表示变换后的新向量。当两个向量处在同一条直线上时（包括同向和反向），表示两者相位相同，只存在一个（可正可负的）实数乘子，Axx,Eigshow(A4)产生的图形,eigshow(1,2;2

4、,2)的图形,将eigshow(1,2;2,2)粘贴到命令窗,A是对称实矩阵的情况,特别要注意A是对称实矩阵的情况，所谓对称矩阵是满足ATA的矩阵。，对22矩阵，只要求A(1,2)A(2,1)。例如令A=1,2;2,2 再键入eigshow(A)，这时的特点是：Axx出现在Ax椭圆轨迹的主轴上，所以两个特征值分别对应于单位圆映射的椭圆轨迹的长轴和短轴。此时A的特征值为-0.5616和 3.5616，可以和图形对照起来看。,例9.2 斜体字的生成,数据矩阵表示英文大写空心字母N的各个节点（1）用plot语句在子图1中画出其形状；（2）取 作为变换矩阵对x进行变换，并在子图2中画出其图形；画图的要

5、点是要在给定的数据右方，补上第一点的坐标，使画出的图形封闭。,程序与图形结果,x00,0.5,0.5,6,6,5.5,5.5,0;0,0,6.42,0,8,8,1.58,8;xx0,x0(:,1);%把首顶点坐标补到末顶点后A1,0.25;0,1;yA*x;subplot(1,2,1),plot(x(1,:),x(2,:)subplot(1,2,2),plot(y(1,:),y(2,:)画出的两个图形如右：,平移运动不能用二维变换实现,刚体在平面上的运动要用两个平移和一个转动来描述，转动可以从上面的线性变换A5得到，但平移yxc却不是一个线性变换。因为：（1）设yaxac；ybxbc；则它们的

6、和为yyaybxaxb2cxc，可见，它对加法不封闭；（2）设yaxac；将它乘以常数k，ykyak(xac)kxakckxacxc，可见，它对乘法也不封闭；就是说，这不符合线性变换的规则，x和y 不属于同一个向量空间，无法用矩阵乘法来实现平移变换yxc。,平面运动模型的齐次坐标系,把平面问题映射到高维的空间来建立方程，有可能把x和y由扩展了的向量空间来覆盖。把原来通过原点的平面沿垂直方向提高一个单位，与原平面保持平行，于是原来的x就用三维向量来表示为：这样的坐标系称为齐次坐标系。,刚体平面运动用线性变换描述,此时可以把平移矩阵写成：因而平移运动y就可用x经线性变换实现了。这个方法在研究刚体平

7、面运动时非常有用。,同时有旋转和平移的情况,对象若同时有旋转和平移，则可以分别列出旋转矩阵和平移矩阵。不过此时的旋转矩阵也要改为33维，这可以把上述A5中增加第三行和第三列，置A(3,3)1，其余新增元素为零。这就是既包括平移，又包括转动的平面齐次坐标系内的变换矩阵。,例9.3 刚体平面运动描述,设三角形的三个顶点坐标为(1,1),(1,1),(0,2)，今要使它旋转30度，右移2，上移3，以试设计变换矩阵A，并画出变换前后的图形。解：程序的要点是：1。列出三角形的数据矩阵2。扩展为齐次坐标(第三行加1)3。平移和转动变换矩阵也要用三维的变换矩阵4。按变换次序左乘5。绘图,9.3 空间线性变换

8、的几何意义,三维空间线性变换最直接的几何意义和应用价值可以从飞行器的三维转动坐标中得到解释。飞行器在空中可以围绕三个轴旋转。假如它在向北飞行，机头正对北方，则它围绕铅垂轴的旋转角称为偏航角（Yaw），它描述了飞机左右的偏转，用u表示；围绕翼展轴的旋转角称为倾斜角（Pitch），它描述了飞机俯仰姿态，用v表示；围绕机身轴的旋转角称为滚动角（Roll），用w表示；u,v和w三个变量统称为欧拉角，它们完全地描述了飞机的姿态。,演示程序quatdemo,演示画面的说明,画面中。左方为飞行器在三维空间中的模型，其中红色的是飞行器。右上方为三个姿态角u,v,w的设定标尺和显示窗，右下方为在地面坐标系中的另

9、外的三个姿态角：方位角、俯仰角和倾侧角。左下方还有【静态】和【动态】两个复选钮，我们只介绍【静态】，读者可自行试用【动态】进行演示。用键入参数或移动标尺的方法分别给u,v,w赋值并回车后，就可以得出相应的飞行器姿态，同时出现一根蓝色的线表示合成旋转的转轴。,例9.4 程序的实现方法。,把飞行器的三维图像用N个顶点描述，写成一个3N的数据矩阵G。用plot3命令时按顶点连线能绘制出飞行器的外观。例如以下的程序ag904a即可画出一个最简单的飞行器立体图。Gw=4,3,0;4,3,0;0,7,0;4,3,0;%主翼的顶点坐标 Gt=0,3,0;0,3,3;0,2,0;0,3,0;%尾翼的顶点坐标

10、G=Gw,Gt%整个飞行器外形的数据集 plot3(Gw(1,:),Gw(2,:),Gw(3,:),r),hold on plot3(Gt(1,:),Gt(2,:),Gt(3,:),g),axis equal,围绕各个轴的旋转变换矩阵,飞行器围绕各个轴的旋转的结果，表现为各个顶点坐标发生变化，也就是G的变化。只要把三种姿态的变换矩阵Y,P和R乘以图形数据矩阵G即可。其中,综合旋转的变换矩阵,单独变化某个姿态角所生成的图形由G1Y*G，G2P*G，G3R*G算出，如果同时变化三个姿态角，则最后的图像数据成为GfY*P*R*GQ*G。这里假定转动的次序为：先滚动R，再倾斜P，最后偏航Y，由于矩阵乘

11、法不服从交换律，转动次序不同时结果也不同。用MATLAB实现的程序ag904b如下：syms u w vY=cos(u),sin(u),0;sin(u cos(u),0;0,0,1)R=1,0,0;0,cos(w),sin(w);0,sin(w),cos(w)P=cos(v),0,sin(v);0,1,0;sin(v),0,cos(v)Q=Y*P*R,空间的齐次坐标系,三维空间考虑了平移运动后，如同二维情况那样，也必须扩展一维，成为4N数据集G4，成为空间的齐次坐标系：在四维空间的44变换矩阵为：其中c1,c2,c3为在三个轴x1,x2,x3方向上的平移距离。这种方法在机器人运动学研究中很有用

12、处。,9.4 基变换与坐标变换,在线性空间中常常需要进行坐标变换。用下图可以形象地说明这点。按照左图的笛卡儿坐标，x向量应该表为(1,6)，这是x按标准基e1,e2度量的结果，在斜坐标纸上的x点坐标就成为沿b1方向为2个单位而沿b2方向3个单位，即(-2,3)了。这反映了不同的基对坐标值的影响。,基坐标变换的公式,设线性空间Rn中的两组基向量u 和v都是n维列向量，它们在基准坐标系中的n个分量都是已知的，因此u和v都可表示为nn矩阵。如果Rn中的一个向量w在以u为基的坐标系内的坐标为wu（n1数组），在以v为基的坐标系内的坐标为wv（n1数组），它们在基准坐标系内的坐标应分别为u*wu和v*w

13、v，这两者应该相等。u*wu v*wv（9.18）所谓基坐标的变换就是已知wu，求出wv。将上式左右均左乘以inv(v)，得到（9.19）可见，坐标变换矩阵P可由u和v求得：P(uv)v u（9.20）,基变换的算例9.5,已知R4空间的两组基向量u,v如下：试求把u变换为v的坐标变换矩阵P(uv)。解的方法为:输入u和v矩阵后 键入uv，得到给出某点w的u坐标wu,即可求其v坐标wv=P*wu,9.5 特征值的MATLAB求法,MATLAB提供了计算方阵的特征值和特征向量各步骤的函数。这三个步骤是：（1）用f=poly(A)可以计算方阵A的特征多项式系数向量f；（2）用lamda=roots

14、(f)可以求特征多项式f的全部根lamda（表示为列向量）；（3）用函数p=null(lamda*I-A)直接给出基础解p,将n个特征列向量p排在一起，就是的特征向量矩阵。【例9.6】求矩阵 的特征值和特征向量。,解题的程序ag906,按上面所说的步骤，解题的程序ag906为A=3,2,4;2,0,2;4,2,3;%输入系数矩阵f=poly(A),%求特征多项式(步骤1)r=roots(f),r=real(r),%求特征根(步骤2),并去掉的虚部B1=r(1)*eye(3)-A;%后面可能需要插入一条语句%B1=rref(B1,1e-12),%为保证矩阵B1奇异性而设的语句p1=null(B1

15、,r)%求第一个特征值，的基础解系B2=r(2)*eye(3)-A;p2=null(B2,r)%求第二个个特征值的基础解系B3=r(3)*eye(3)-A;p3=null(B3,r)%求第三个个特征值的基础解系,程序运行的结果,程序运行的结果为：f=1-6-15-8（特征多项式系数向量）r=8.0000（三个特征根即特征值,后两个是重根）-1.0000+0.0000i（微小虚数可用r=real(r)去除）-1.0000-0.0000i对第一个特征值8的特征向量为空矩阵 p1=Empty matrix:3-by-0，即无特征向量对第二、三个特征值（重根）-1的特征向量为 p2=-0.5000-1

16、.0000 1.0000 0 0 1.0000 对第一个特征值所以找不到特征向量是由计算误差引起的，计算误差使得本应为零的det(B1)成了一个很小的数（-8.7130e-013），B1*x=0就没有非零解了。,结果的讨论,解决的方法是将B1化为行阶梯形式，因为在rref函数中可以人为地设定容差，例如把容差设为110-12，令B1=rref(B1,1e-12),这样变换出的行阶梯形式与原矩阵B1是等价的，但它将会把一些小的数看作零，从而保证系数矩阵的奇异性。在上述程序中求B1和求p1的两条语句之间，插入这条语句后执行，就可求出对应的特征向量为 p1=1,0.5,1T。把上述三个特征向量并列，可

17、以得到特征向量矩阵（元素取整后）和特征值矩阵为：,用eig函数计算结果,实际上MATLAB已经把求特征根和特征向量的步骤集成化，其中也包括了处理计算误差的功能，所以一条命令就解决问题了。这个功能强大的子程序名为eig(特征值英文是 p,lamda=eig(A)输出变元中的lamda是特征值，p是特征向量。把例9.6的系数矩阵A代入，即可得到：解出的特征向量怎么与前面计算的不同？其实特征向量本来不是唯一的，这两个特征向量只是相当于基础解系的一组基。只要检查由两种方法算出的两组特征向量是否相关，如果相关，说明解都是正确的。,例9.7,设矩阵，求以下矩阵的特征值；a）A；b）本题只要求出特征方程和特

18、征根，它的计算程序ag907是比较简单的：A=1,-1,2;0,2,1;0,0,-1;B=2*A3+A-5*eye(3);C=inv(A)+eye(3)fa=poly(A),ra=roots(fa)fb=poly(B),rb=roots(fb)fc=poly(C),rc=roots(fc),程序ag907的运行结果,程序运行结果：这意味着三个特征多项式分别为：相应的各组特征根为：,9.6 对称矩阵与二次型主轴,对称矩阵的特点是所有元素关于主对角线对称，即AA。所以对称矩阵一定是方阵。前面曾要求读者特别注意A是对称矩阵时x与Ax的对应关系，其特点就是Ax呈椭圆形状，在椭圆的两个主轴方向，Ax与x

19、在一条直线上长度差倍，即Ax x。当Ax与x方向相同时，为正数；当Ax与x方向相反时，为负数；22变换有两个特征值，在相互正交的两个主轴方向，各有一个。作为22正交变换的一个应用，我们来看看它对二次型图形的影响。二次型本身已经不是线性范围，不属于线性代数的范畴。现在要研究的是基坐标的线性变换对二次型图形发生何种影响。,例9.8 二次型例,设A=5,-2;-2,5，则令A的二次型xT*A*x等于常数得到的是一个椭圆方程，其图形如下图(a)所示。如果做一个基坐标的旋转变换，让坐标轴转过45度，此椭圆的主轴就与新的坐标方向y1,y2相同，如图(b)所示，即令y1x1cosx2siny2x1sinx2

20、cos用矩阵乘法表为,线性变换后的二次型,其逆变换R为，因此 用此变换式代入二次型的表达式，有本题中,=45，代入P和R,可得于是得到,二次型主轴等价于矩阵对角化,所以从几何图形上寻找二次型主轴的问题，在线性代数中就等价于使矩阵经过正交变换或相似变换R（注意这又是一个几何名词，说明被变换的图形的形状和尺寸保持不变），使矩阵A对角化。图中的(c)和(d)表示了对另一种双曲线二次型（它的两个特征值一正一负）的坐标变换，求主轴的方法就是把矩阵A对角化。找其主轴的大小和方向，也就是找它的特征值lamda和特征向量e。,双曲线二次型的算例,根据列出程序A=1,-4;-4,-5lamda,e=eig(A)

21、或R=orth(A)得到把两个特征向量e并列起来,即正交矩阵。lamda就是对角化的矩阵D，故标准化的二次型方程为,高维空间的算例9.9,化二次型为标准型。解：可以看出系数矩阵A，程序ag909n为A1,1,3;1,2,1;3,1,5,Rorth(A),Dinv(R)*A*R得知二次型最后的标准型为其中,例9.10,设有对称实矩阵，试求。解：建模：矩阵指数只有在A是对角矩阵时才可按对角元素直接计算，即：注意它的非对角元素不符合指数运算规则。,解题的思路,因此首先要把A对角化。使A=pDp-1，将此式代入y的表示式中，可得：求出其特征根矩阵D和特征向量矩阵p。编MATLAB程序ag910如下：A

22、=-2,4,1;4,-2,-1;1,-1,-3,p,D=eig(A)运行此程序，得 p=-0.6572-0.2610 0.70710.6572 0.2610 0.70710.3690-0.9294 0.0000D=-6.5616 0 00-2.4384 00 0 2.0000,解题的结果,由此得到：算式中的最后一步可以用下列语句来实现。y=p*exp(D(1,1),0,0;0,exp(D(2,2),0;0,0,exp(D(3,3)*inv(p)得到MATLAB中有直接计算矩阵指数的函数expm(A)，可以用它来检验所得结果。,9.7 奇异值分解简介,奇异值分解(Singular Value D

23、ecomposition)对于矩阵的应用具有重要意义。比如矩阵的条件数（对计算精度）、矩阵广义逆（解超定方程）都与矩阵奇异值有关。美国多数教材都有，而我国传统的教材大多不提，大纲也不要求。本书折衷，对这个内容作一扼要介绍。定义 设矩阵，若存在非负实数和n维非零向量u，m维非零向量v，使得（9.7.1）则称为A的奇异值，分别称u和v为对应于奇异值的右奇异向量和左奇异向量。由式（9.7.1）可得（9.7.2）（9.7.3）,奇异值的定义,因此2为ATA(nn)的特征值，也是AAT(mm)的特征值，说明ATA与AAT的特征值均为非负实数；而分别是和对应于特征值2的特征向量。ATA与AAT的非零特征值

24、的个数（重特征值按重数计算）等于A的秩。设 是ATA的全部特征值，,称为矩阵A的全部奇异值。矩阵的奇异值分解定理：设A是mn矩阵，设1,2,n是A的奇异值，12rr+1=r+2=n=0，则A=USVT，其中U是m阶正交矩阵，V是n阶正交矩阵，而。此式也可以表示为：(9.7.4)其中，ui是矩阵U的第i列，vj是矩阵V的第j列。,奇异值分解定理证明(1),证 不失一般性，设mn，R(A)=r，v1,v2,vn是ATA对应于特征值12,22,n2的标准正交特征向量组，12r r+1=r+2=n=0。令ui=A*vi/i，其中i=1,2,n，故有因此，u1,u2,ur组成了一个标准正交向量组。如果r

25、m,可以把向量组u1,u2,ur扩展为向量组u1,u2,ur,um，从而构成Rm空间的一个基。我们将证明AV=US，其中S是主对角元素为奇异值1,2,r,0,0的mn对角矩阵，而U和V分别是V=v1,v2,vn和U=u1,u2,ur,um，它们分别是m阶和n阶的正交矩阵。,奇异值分解定理证明(2),上式中的矩阵S由m行n列组成，按假设mn，所以左上角的rr矩阵Sr决定了S的秩为r。也就是行列式det(Sr)=1 2r0，而det(Sr1)=12r1=0。可见r1要小到什么程度才能看成零的问题，会直接影响秩r的判别。命题得证。,奇异值分解的应用：数字秩,我们曾把矩阵的秩定义为其最简行阶梯形的非全

26、零行数。对于用浮点计算机算法进行计算的工程实际中的矩阵，不太可能见到真正的全零行。在进行奇异值分解时，也不大会出现r+1=r+2=n=0的情况，通常这些i都有微小的值。这时把A的秩看作多少，就取决于把多小的i值看成零。如果取A的秩为r，那就意味着Sr中只取了r个特征值1,2,r，而认为r+1及以后各个更小的奇异值都等于零，经过这样的数值处理后得到的秩称为数字秩(Numerical Rank)。多小的数该看成零呢？这取决于所用计算软件的精度。就MATLAB而言，其基本部分采用的是IEEE浮点双精度格式标准，它的精度(即eps)是2.2410-16，当数eps时，它与数1相加时会被略去1+eps=

27、1。普遍地说，当 时，MATLAB将无法分辨a+b与a的差别。,数字秩和条件数,把这个概念推广到矩阵，假设A是nn方阵：式中uiviT中所有的nn个元素都是小于1的。决定各项数值大小的主要因素是i。如果i+1比构成的主要分量的和小得多，以致加上它根本不起作用，那末可以认为后的各项都为零。A的秩只能取作r。用式中的最大特征值与阶数n的乘积近似表征的大小，则可以把设为一个门限。大于这个门限特征值的数目r就是A的数字秩，把小于等于这个门限的特征值都看成零。把最大与最小奇异值的比定义为矩阵的条件数(Condition Number)。矩阵的条件数愈大，意味着愈小，它就愈接近奇异；反之，条件数愈小，它就

28、离开奇异门限愈远。,例9.11,设，求它的各奇异值，及条件数。解：这样阶次的问题，只能用计算机来解了。程序为：A=2,7,9,-5,4;-9,-9,5,3,-2;-2,5,-1,-3,5;-4,9,0,9,-4U,S,V=svd(A),condA=S(1,1)/S(4,4)运行结果为：,运行结果,不难验证：U,V都是规范正交矩阵，都有A的四个奇异值为：16.5933 13.9809 11.2638 5.9432，A的条件数为：condA=16.5933/5.9432=2.7920,9.8.1 人口迁徙模型,设在一个大城市中的总人口是固定的。人口的分布则因居民在市区和郊区之间迁徙而变化。每年有6

29、%的市区居民搬到郊区去住，而有2%的郊区居民搬到市区。假如开始时有30%的居民住在市区，70%的居民住在郊区，问10年后市区和郊区的居民人口比例是多少？30年、50年后又如何？这个问题可以用矩阵乘法来描述。把人口变量用市区和郊区两个分量表示，一年以后，市区人口为xc1(10.06)xc00.02xs0，郊区人口xs1 0.06xc0(10.02)xs0，,问题的矩阵描述,用矩阵乘法来描述，可写成：从初始到k年，此关系保持不变，因此上述算式可扩展为，故可用程序ag981n进行计算：A0.94,0.02;0.06,0.98,x00.3;0.7x1A*x0,x10A10*x0,x30A30*x0,x

30、50A50*x0得到：,本题特征值和特征向量的意义,无限增加时间k，市区和郊区人口之比将趋向一组常数0.25/0.75。为了弄清为什么这个过程趋向于一个稳态值，我们改变一下坐标系统。在这个坐标系统中可以更清楚地看到乘以矩阵A的效果，先求A的特征值和特征向量，得到令它是特征向量的整数化，得到,9.8.2 产品成本的计算,某厂生产三种成品，每件产品的成本及每季度生产件数如表9.1及表9.2所示。试提供该厂每季度在每种产品上的成本表。解：应当用矩阵来描述此问题，列出成本矩阵为M，季度产量矩阵为P,本例矩阵相乘的变换意义,将M和P相乘，得到的矩阵设为Q，Q的第一行第一列元素为Q(1,1)0.14000

31、0.320000.1558001870不难看出，Q表示了夏季消耗的原材料总成本。从线性变换的角度来看，Q矩阵把以件数为单位的产品空间映射到了以元为单位的成本空间。,9.6.3 情报检索模型,假如数据库中包括了n个文件，而搜索所用的关键词有m个。可以把数据库表示为mn的矩阵A。比如有7本书，6个关键词x（初等，代数，矩阵，理论，线性，应用）：则A就是67的矩阵。书名中有此关键词的就将该对应元素置1。搜索结果可以表示为乘积yATx，它是n1列向量。于是y的各个分量就表示各书与搜索向量匹配的程度。y值最大的元素对应于匹配最好的书籍，是读者可能最需要的。可见变换有很广泛的意义。在本例中，它是从关键词子

32、空间变换为文献目录的子空间。,结束语顺口溜,“一、两、三，三、两、一，还要三句作对比！”“一、两、三，”一个中心，两个方向，三种类型一个中心就是消元法生成行阶梯形式：用gauss,rref,lu,rank,inv,det,cond 等函数两个方向是从行向(方程)和列向(变量)分析行阶梯形式；用到rref,枢轴变量和自由变量三种类型即欠定，适定和超定。x=Ab的特点:x=inv(A)*b,x=pinv(A)*b,x0=null(A)。,结束语顺口溜(续),“三、两、一，”三个视点，1.解方程x，已如前述2.向量张成空间，内积xT*y，向量的夹角和模，基向量和正交化等，用到norm(v),qr(A)等 3.像空间b的特点，线性变换的特征值和特征向量，坐标变换，矩阵的对角化。用到 eig(A),orth(A)等两个动力，需求牵引和技术推动。一个目标，后续课的应用需求为主。,结束语顺口溜(续),“还要三句作对比！”线性代数抽象吗？解决的办法是？线性代数冗繁吗？解决的办法是？线性代数枯燥吗？解决的办法是？,