《Matlab在线性代数中的应用1.ppt》由会员分享,可在线阅读,更多相关《Matlab在线性代数中的应用1.ppt(44页珍藏版)》请在淘文阁 - 分享文档赚钱的网站上搜索。
1、目标要求会给矩阵赋值会进行矩阵的基本运算,包括:加、减、数乘,乘法,转置,幂等运算会用命令inv计算矩阵的逆会用命令det计算行列式;会用命令rank计算矩阵的秩;会用命令rref把矩阵变为行最简型;会用命令rref计算矩阵的逆会用命令rref解方程组的解会用命令rref找出向量组的最大无关组会用命令null计算齐次线性方程组的基础解系会用左除运算计算非齐次方程组的特解会用命令orth把向量组正交规范化会用命令eig计算矩阵的特征值和特征向量会用命令eig把二次型标准化会用命令eig判断二次型的正定性1 矩阵赋值赋值语句一般形式变量=表达式(或数)如:输入a=1 2 3; 4 5 6; 7 8
2、 9 显示a = 1 2 3 4 5 6 7 8 9 输入x=-1.2 sqrt(3) (1+2+3)/5*4 显示x=-1.2000 1.7321 4.8000规则:矩阵元素放在方括号中,元素之间以空格或逗号分规则:矩阵元素放在方括号中,元素之间以空格或逗号分隔,不同行以分号分隔,语句结尾用回车或逗号将显示结隔,不同行以分号分隔,语句结尾用回车或逗号将显示结果果 基本赋值矩阵ones(m,n), zero(m,n), magic(n), eye(n), rand(m,n), round(A)如:输入 f1=ones(3, 2) 显示 f1= 1 1 1 1 1 1 输入 f2=zero(2,
3、 3) 显示 f2= 0 0 0 0 0 01 矩阵赋值 输入 f3=magic(3) 显示 f3= 8 1 6 3 5 7 4 9 2 输入 f4=eye(2) 显示 f4= 1 0 0 12 矩阵的基本运算 矩阵算术运算书写格式与普通算术相同,包括加、减、乘、除。可用括号规定运算的优先级。 Matlab将矩阵加、减、乘的程序编为内部函数,只要用+,-*做运算符号就包含阶数检查和执行运算的全过程 两相加矩阵有一个是标量时,Matlab承认算式有效,自动把标量扩展为同阶等元素矩阵 如:键入 X=-1 0 1; Y=X-1 得 Y= -2 -1 0 矩阵除法 矩阵求逆 inv(A),如果det(
4、A)等于或很接近零,Matlab会提示出错 “左除”与“右除”,左乘或右乘矩阵的逆,A或/A 2 矩阵的基本运算 幂运算幂运算 A*A*A=A5 转置转置 理论学习中,A的转置表示为AT,在Matlab中用“”表示3 行列式与方程组求解 相关命令U=rref(A), 对矩阵A进行初等行变换,矩阵U为A的最简梯矩阵det(A), 计算矩阵A的行列式rank(A),计算矩阵A的秩B(: , i)=b, 把向量b赋给矩阵B的第i行A(i, j), 引用矩阵A中第i行j列的元素A, eye(5), 创建510矩阵,前5列为A,后5列为单位矩阵syms x, 定义x为符号变量 3 行列式与方程组求解 逆
5、矩阵各种求法:clearA=-7,-2,-6,4,6;1,3,-6,3,11;3,-11,9,5,-2;-3,0,-2,9,-3;7,30,-18,11,4;% 1.命令法:An1=inv(A)% 2.幂运算法:An2=A-1% 3.右除法:An3=eye(5)/A % eye(5)为5阶单位矩阵% 4.左除法:An4=Aeye(5)% 5.初等行变换法:B=rref(A,eye(5); % 对矩阵A , I 进行初等行变换 % B为矩阵A的最简行阶梯矩阵if(rank(B(:,1:5)=5) % 判断最简行阶梯矩阵B的前5列是否为单位阵 An5=B(:,6:10) % 取出矩阵的后5列,并显
6、示elsedisp(A不可逆);end思考:如何用求逆阵或初等变换法解方程组?思考:如何用求逆阵或初等变换法解方程组?3 行列式与方程组求解% 求解符号行列式方程clear % 清除各种变量syms x % 定义x为符号变量A=3,2,1,1;3,2,2-x2,1;5,1,3,2;7-x2,1,3,2 D=det(A) % 计算含符号变量矩阵A的行列式Df=factor(D) % 对行列式D进行因式分解 % 从因式分解的结果,可以看出方程的解X=solve(f) % 求方程“D0”的解2232113221051327132xx 解方程:4 向量组的线性相关性及方程组的通解 相关命令R, s=r
7、ref(A), 把矩阵A的最简梯矩阵赋值给R;s是一个行向量,它的元素由R的首非零元所在列号构成null(A, r), 齐次线性方程组Ax=0的基础解系x0=Ab, 非齐次线性方程组Ax=b的一个特解x0length(s), 计算s向量的维数end, 矩阵的最大下标,最后一行或最后一列find(s), 向量s中非零元素的下标sub(A, k, n), 将A中所有符号变量k用数值n代替 4 向量组的线性相关性及方程组的通解123451234512345123452441623626237364619232521943xxxxxxxxxxxxxxxxxxxx 求非齐次线性方程组的通解4 向量组的线
8、性相关性及方程组的通解% 求齐次线性方程组的通解clear A=2,4,-1,4,16;-3,-6,2,-6,-23;3,6,-4,6,19;1,2,5,2,19; % 输入系数矩阵Ab=-2;7;-23;43; % 输入常数列向量bR,s=rref(A,b); % 把增广矩阵的最简行阶梯矩阵赋给R % 而R的所有基准元素在矩阵中的列号构成了行向量sm,n=size(A); % 矩阵A的行数、列数赋给了变量m、nx0=zeros(n,1); % 将特解x0初始化为N维零向量r=length(s); % 矩阵A的秩赋给变量rx0(s,:)=R(1:r,end); % 将矩阵R的最后一列按基准元素
9、的位置给特解x0赋值disp(非齐次线性方程组的特解为:)x0 % 显示特解x0disp(对应齐次线性方程组的基础解系为:)x=null(A,r) % 得到齐次线性方程组Ax0的基础解系x4 向量组的线性相关性及方程组的通解 当k取何值时方程组有非零解?在有非零解的情况下,求出其基础解系1234123412341234(12 )33303(2)33033(2)30333(11)0k xxxxxk xxxxxk xxxxxk x 已知齐次线性方程组:4 向量组的线性相关性及方程组的通解clearsyms k % 定义符号变量kA=1-2*k,3,3,3;3,2-k,3,3;3,3,2-k,3;3
10、,3,3,11-k; % 给系数矩阵赋值D=det(A); % 算出系数矩阵的行列式Df=factor(D)kk=solve(f); % 解方程“D0”,得到解kk,即k值for i=1:4 AA=subs(A,k,kk(i); % 分别把k值代入系数矩阵A中fprintf(当k=); disp(kk(i); % 显示k的取值 fprintf(基础解系为:n);disp(null(AA) % 计算齐次线性方程组“Ax=0”的基础解系end平板稳态温度的计算平板稳态温度的计算 (1020)/4(2040)/4(1030)/4(4030)/4abcbadcaddbcxxxxxxxxxxxx10.2
11、50.2507.50.25100.25150.25010.251000.250.25117.5abcdxxxx整理为化学方程的配平化学方程的配平 确定x1,x2,x3,x4,使两边原子数相等称为配平,方程为 写成矩阵方程138223242()()()()x C Hx Ox COxH O1234301080020221xxxx 12343010080020002210 xxAxxx -电阻电路的计算电阻电路的计算 设定三个回路电流ia,ib,ic,回路压降的方程为:scbauiiiRRRRRRRRRRRRR001007655554333321信号流图模型 信号流图是用来表示和分析复杂系统内的信号
12、变换关系的工具。右图方程如下。 写成矩阵方程 或x=Qx Pu移项整理,可以得到求信号向量x的公式。ux1G1x2G212221 1xuG xxG x,1212120100 xGxuxGx 信号流图的矩阵解法( I Q ) x= Pu, x = inv( I Q )*Pu 定义系统的传递函数W为输出信号与输入信号之比x/u,则W可按下式求得:W=x/u = inv( I Q )*P因为 得到221110010101GGGGI -Q2111211()11GGG GI -Q121121/1/1xuxuGG G-1x/uI -QP复杂点的信号流图 按右面的信号流图,照上述方法列出它的方程如下: x1
13、 = -G4x3 + u x2 = G1x1-G5x4 x3 = G2x2 x4 = G3x3信号流图的矩阵方程1412152323434000100000000000 xGxxGGxxuxGxxGx 列出的矩阵方程为: 矩阵中的参数是符号而不是数,MATLAB的许多函数(特别是求逆)都可以处理符号,带来了极大的方便。只要在程序第一行注明哪些是符号变量: syms G1 G2 用符号运算工具箱求解 矩阵代数方法的最大好处是可用于任意高的阶次的信号流图,实现传递函数推导的自动化 如下题的MATLAB程序ag863syms G1 G2 G3 G4 G5Q=0,0,G4,0;G1,0,0,G5;0,
14、G2,0,0;0,0,G3,0,P=1;0;0;0W=inv(eye(4)Q)*Ppretty(W(4)运行结果为(4) G1*G2*G3(4)1+G2*G5*G3+G1*G2*G4xWu5 特征向量与二次型 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中特征值一一对应 V, D=schur(A), 矩阵D为对称阵A的特征值所构
15、成的对角阵,V的列为A的单位特征向量,与D中特征值一一对应5 特征向量与二次型 已知矩阵 求其特征值。A=2,-2,-20,-19;-2,16,-9,11;-8,4,-6,1;0,-8,-4,-7;% 1.符号变量法syms k % 定义符号变量kB=A-k*eye(length(A); % 构造矩阵B=(A-kI)D=det(B); % 计算行列式:|A-kI|lamda1=solve(D) % 求|A-kI|=0的符号形式的解% 2.特征多项式法P=poly(A); % 计算矩阵A的特征多项式, 向量P的元素为该多项式的系数lamda2=roots(P) % 求该多项式的零点,即特征值%
16、3.命令法lamda3=eig(A) % 直接求出矩阵A的特征值22201921691184610847A 5 特征向量与二次型 求矩阵的特征值和特征向量,判断是否可对角化,如可以则找出可逆矩阵V,使V-1AV=DA=1,2,3;2,1,3;1,1,2;V,D=eig(A)123213112A 5 特征向量与二次型 用正交变换法将二次型 化为标准型clearA=1,0,0;0,2,2;0,2,2; % 输入二次型的矩阵AV,D=eig(A); % 其中矩阵V即为所求正交矩阵 % 矩阵D为矩阵A的特征值构成的对角阵% 或:V,D=schur(A) % 结果和eig( ) 函数相同disp(正交矩
17、阵为:);Vdisp(对角矩阵为:);Ddisp(标准化的二次型为:);syms y1 y2 y3f=y1,y2,y3*D*y1;y2;y322212312323( ,)224f x x xxxxx x 平面上线性变换的几何意义 例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,
18、2), fill(y1(1,:),0,y1(2,:),0,g) 01100011x几种变换的行列式与特征值1234 1 0det()1,11 , 0 1 0 1 det()1.5,1.0 1.5 , 1 0 0 1det()0.2,0.2 1.0 , 1 0 1.0 1.0det()1,11 , 0. DDDD 111222333444ApApApAp5 0. 0.7071 0.7071 1,0.866 + 0.5i 0.8660.5i ,00.7071i 0 + 0.7071iD55p二维矩阵特征值的几何意义 二维矩阵的特征值表示该变换在原图形的特征向量的方向上的放大量。例如矩阵A1在第一特
19、征向 量 方向的特征值为,即横轴 正方向的增益为1,其结果是把原图中横轴正方向的部分变换到新图的负方向去了; A1在第二特 征向量 的方向的特征值为1(2)=1, 即纵轴正方向的增益为1,因而保持了新图和原图在纵轴方向尺度不变。 1(:,1)0 1p(1)1 10(:,2)1 1p用用eigshow函数看特征值函数看特征值 对于比较复杂的情况,完全凭简单的几何关系去想像是困难的,应当用eigshow函数,联系x和Ax的向量图来思考。 键入eigshow(A4) 。绿色的x表示原坐标系中的单位向量,可以用鼠标左键点住x并拖动它围绕原点转动。图中同时出现以蓝色表示的Ax向量,它表示变换后的新向量。
20、当两个向量处在同一条直线上时(包括同向和反向),表示两者相位相同,只存在一个(可正可负的)实数乘子, Ax x Eigshow(A4)产生的图形eigshow(1,2; 2,2)的图形将eigshow(1,2; 2,2)粘贴到命令窗A是对称实矩阵的情况 特别要注意A是对称实矩阵的情况,所谓对称矩阵是满足AT A的矩阵。,对22矩阵,只要求A(1,2) A(2,1)。例如令A =1,2;2,2 再键入eigshow(A),这时的特点是:Axx出现在Ax椭圆轨迹的主轴上,所以两个特征值分别对应于单位圆映射的椭圆轨迹的长轴和短轴。此时A的特征值为 -0.5616和 3.5616,可以和图形对照起来看
21、。人口迁徙模型 设在一个大城市中的总人口是固定的。人口的分布则因居民在市区和郊区之间迁徙而变化。每年有6%的市区居民搬到郊区去住,而有2%的郊区居民搬到市区。假如开始时有30%的居民住在市区,70%的居民住在郊区,问10年后市区和郊区的居民人口比例是多少?30年、50年后又如何? 这个问题可以用矩阵乘法来描述。把人口变量用市区和郊区两个分量表示,一年以后,市区人口为xc1 (10.06) xc00.02xs0,郊区人口xs1 0.06xc0 (10.02)xs0,问题的矩阵描述 用矩阵乘法来描述,可写成:从初始到k年,此关系保持不变,因此上述算式可扩展为,故可用程序ag981n进行计算:A0.
22、94,0.02;0.06,0.98, x00.3;0.7x1A*x0, x10A10*x0, x30A30*x0, x50A50*x0得到:110.94 0.020.3 0.29600.06 0.980.7 0.7040csxx 10 xAx2kkk-1k-20 x = Ax= A x= A x1103050 0.2960 0.2717 0.2541 0.2508, 0.7040 0.7283 0.7459 0.7492xxxx本题特征值和特征向量的意义 无限增加时间k,市区和郊区人口之比将趋向一组常数0.25/0.75。为了弄清为什么这个过程趋向于一个稳态值,我们改变一下坐标系统。在这个坐标
23、系统中可以更清楚地看到乘以矩阵A的效果,先求A的特征值和特征向量,得到 令它是特征向量的整数化,得到 0.9200 0 -0.7071 -0.3162, 0 1.0000 0.7071 -0.9487lamdae1211,13uu 0210.250.05(0.92)kkkxA xuu6 直线和平面的快速绘制程序 平面曲线的快速绘制程序平面曲线的快速绘制程序 ezplot( ,a,b) 引号中函数可以只有一个自变量,代表显函数 ezplot(f(x), a,b) 系统将在 a x b的范围内画出 f = f(x) 引号中的函数若有两个自变量,那就代表隐函数,其典型格式为 ezplot(f(x,y
24、), a,b) 系统将在 a x b的范围内画出 f(x,y)=0。 a,b的默认值为-2, 2曲线快速绘制举例 ezplot(x1+0.2*x23+1) 画出 在x1=-2, -2的曲线画多条曲线可按下列方法编程 s1=x1+0.2*x23+1% 方程1 s2=3*x1+2*x2+3 % 方程2 ezplot(s1),hold on% 画方程1,保持 ezplot(s2),grid on% 画方程2,加网格 x1,x2=solve(s1,s2)% 解联立方程1,23120.210 xx解的结果 图5.3.1 用ezplot在 同 一 图 中 画 两 根 曲 线 -1.1 -2.2171612
25、389003691410154884062240 .21716123890036914101548840622401x 02 1.8257418583505537115232326093360-1.8257418583505537115232326093360 x解的说明 在线性代数中,遇到的都是一次函数,所以不会出现曲线。我们故意用一个三次曲线来说明ezplot的用法,是为了使读者知道,这个命令不限于画直线。 MATLAB用solve命令解题是采用了符号运算工具葙,它的数字精度是32位十进制,而不是一般数值计算时的16位十进制。尽管在本题中有两有效数字就够了 平面的快速绘制命令ezmesh ezmesh (f(x,y),a,b,c,d)可以绘制很多函数的曲面。 其第一输入变元可以直接输入用MATLAB语句写出函数的形式,引号中的函数只能是显函数z=f(x,y)中的f(x,y)。它应该有两个自变量,注意要用单引号括起来。 第二输入变元为自变量的取值范围,默认情况下其x,y的取值范围都是-2,2。 用ezmesh快速绘制平面举例 可以用以下的程序ag501在一张图内画两个平面 clear, clf ezmesh(3*x1+2*x2+3) hold on ezmesh(x1-2*x2+1)
限制150内