第08课-Matlab机械振动建模分析方法

机械振动(MechanicalVibration)交通与车辆工程学院刚宪约2007年9月27日第二课导论：Matlab基础主要内容Matlab基本操作特征值问题积分变换线性常微分方程主要内容Matlab基本操作特征值问题积分变换线性常微分方程Matlab基本操作数组与矩阵操作求解线性方程组多项式运算绘图操作编程初步MATLAB的概况MATLAB（MatrixLaboratory）是由美国MathWorks公司开发的；MATLAB具有数值计算、符号计算和图形可视化三大功能；MATLAB的基本数据单位是矩阵；具有丰富的工具箱源程序的开放性程序解释执行机械振动(MechanicalVibration)交通与车辆工程学院刚宪约2007年9月27日MATLAB是一种高级程序语言MatrixLaboratory只需输入几行简单的命令，而不必用编程语言（C++、VB等）去编写一大串程序，就可以快速地得到结果并解决问题应用MATLAB解线性代数Ax=BA,B,x是矩阵。当维数都是1时，则成为一般的线性方程。A=3,B=6x=B/A=2A=[5,12,8;6,5,8;9,6,10]B=[7;11;7]x=A\B=[-4.7073;-1.2439;5.6829]矩阵的输入：行内用空格或逗号行间用分号A=[5128;658;9,6,10]基本矩阵ones(3),ones(2,5)zeros(3),zeros(2,5)eye(3),eye(3,4),eye(4,3)矩阵大小[r,c]=size(A),size(A)r、c分别为A的行数和列数命令的编辑修改方向键及控制键可以编辑修改已输入的命令r1=(1+sqt(5))/2sqrt命令超过一行，用…连接下一行s=1-1/2+1/3-1/4+1/5-1/6+1/7...-1/8+1/9-1/10+1/11-1/12;各符号间的空格可有可无，加上空格增加可读性。-1/12后的分号使结果不回显常量和变量•缺省变量名ans•piπ•i,j复数(-1)^(1/2)•NanNot-A-Number非数0/0•inf无穷大1/0复数及其运算•表示：c1=1-2i,c2=3+sin(.5)*j•运算：c=c2^2•极坐标：b1=abs(c1),b2=angle(c1)•直角坐标：b3=real(c1),b4=imag(c1)c2=3.0000+0.4794ic=8.7702+2.8766ib1=2.2361b2=-1.1071（弧度）b21=b2*180/pi(角度）b3=1b4=-2矩阵运算+-*A=[123;456;789]B=[123;456;789]C=A+BD=A-Bhelpfixfix(3.6)helprandrand(3,3)rand(3,3)A=fix(15*rand(4,3))B=fix(15*rand(4,4))C=A*BD=B*A矩阵运算'invA=fix((10+17i)*rand(3,3))E=A'A=[5,12,8;6,5,8;9,6,10]B=[7;11;7]x=A\B=[-4.7073;-1.2439;5.6829]x=inv(A)*B简单数组的访问形式•x=[00.1*pi0.2*pi0.3*pi0.4*pi0.5*pi]x=00.31420.62830.94251.25661.5708•y=sin(x)y=00.30900.58780.80900.95111.0000•x(3)ans=0.6283x=00.31420.62830.94251.25661.5708y=00.30900.58780.80900.95111.0000•x(1:5)ans=00.31420.62830.94251.2566•x(2:2:7)ans=0.31420.94251.5708•y([5261])ans=0.95110.30901.00000如何创建数组•X=[abcd];由指定元素建立行向量•x=初值：终值从初值到终值，增量为1•x=初值：增量：终值•x=linspace(初值,终值,n)指定元素个数，线性等距分布x=linspace(1,100,10)x=linspace(0,100,10)x=linspace(0,100,9)x=linspace(0,100,11)•x=logspace(初值,终值,n)x=logspace(0,2,11)指定元素个数，对数等距分布数组的格式•a=[1,2,3,4,5]:行向量•a=[12345]:行向量•b=[1;2;3;4;5]:列向量•c=[123;456;789]:矩阵•d=c’共轭转置纯量与数组的运算•d=[123;456;789];•2*d-1•ans=1357911131517数组之间的运算•b=[123]•c=[1;2;3]b*c与b.*c的区别b*c=14c*b=[123;246;369]b.*c???Errorusing==.*Matrixdimensionsmustagree.c.*b???Errorusing==.*Matrixdimensionsmustagree.b.*c'=[149]c.*b'=[1;4;9]数组操作•A(r,c)：指定某一元素•A(r,:)：指定某一行元素•A(:,c)：指定某一列元素A=[123;456;789];A(1,2)A(2,:)A(:,3)数组规模•[r,c]=size(A)r、c分别为A的行数和列数A=[123;456;789][r,c]=size(A)A(2,:)•n=length(b)n为向量b的长度n=length(A(2,:))多项式求根•多项式表示法•用行向量P表示•求已知多项式的根•r=roots(P)•求一组根所对应的多项式•P=ploy(r)•P=[1-12025116];表示多项式x4-12x3+25x+116•P=[1-12025116];•r=roots(P)•r=11.74732.7028-1.2251+1.4672i-1.2251-1.4672i•p1=poly(r)•p1=1.0000-12.0000025.0000116.0000多项式的基本运算•多项式相乘c=conv(a,b)•如:a(x)=x3+2x2+3x+4b(x)=x3+4x2+9x+16求c(x)=a(x)*b(x)=x6+6x5+20x4+50x3+75x2+84x+64•多项式相加、减，即矩阵相加，但规模需相同如a(x)=x3+2x2+x+1;b(x)=x+4求c(x)=a(x)+b(x)•a=[1234];•b=[14916];•c=conv(a,b)•a=[1211];•b=[0014];•c=a+b;多项式的基本运算•多项式项除•[q,r]=deconv(a,b)•两多项式项除，q为商、r为余项•多项式微分•ployder(p)44*072*169*248*320*46*51*6•[q,r]=deconv([121],[11])•q=11r=000•polyder([162048697244])•ans=6308014413872基本二维画图命令Plot的用法•plot(x,y)•plot(y)•plot(x1,y1,x2,y2……)x为横轴，y为纵轴x=linspace(0,2*pi,30);y=sin(x);z=cos(x);plot(y);pause;plot(x,y);pause;plot(x,y,x,z)01234567-1-0.8-0.6-0.4-0.200.20.40.60.81图形标注•Title(‘text’)：图形顶部加标题•xlabel(‘text’)：横轴标注•ylabel(‘text’)：纵轴标注•text(x,y,‘text’)：(x,y)处加标注•gtext(‘text’):鼠标点击处加标注图形标注举例x=linspace(0,2*pi,30);y=sin(x);plot(x,y)title('y=sin(x)')xlabel('x')ylabel('sin(x)')text(pi/2,0.8,'max')gtext('min')01234567-1-0.8-0.6-0.4-0.200.20.40.60.81y=sin(x)xsin(x)maxmin图例说明•Legend(‘s1’,’s2’……)建立一个图例说明盒在右上角•legendoffplot(x,y,'b',x,z,'r')legend('sin(x)','cos(x)')pauselegendoff01234567-1-0.8-0.6-0.4-0.200.20.40.60.81sin(x)cos(x)其它二维图形命令•bar条形图•例：x=-2.9:0.2:2.9;y=exp(-x.*x);bar(x,y)-3-2-1012300.10.20.30.40.50.60.70.80.91其它二维图形命令•stairs阶梯图•例：x=-2.9:0.2:2.9;y=exp(-x.*x);stairs(x,y)-3-2-1012300.10.20.30.40.50.60.70.80.91其它二维图形命令•stem火柴杆图•例：x=linspace(0,4*pi);y=exp(-0.3*x).*sin(x);stem(x,y);title('stemplot')02468101214-0.3-0.2-0.100.10.20.30.40.50.60.7stemplot基本三维曲线命令•plot3(x,y,z)x,y,z为同规模的向量或矩阵•plot3(x,y,x,’s’)字符串s指定线型、颜色、点标记•plot3(x1,y1,z1,’s1’,x2,y2,z2,’s2’,……)基本三维曲线命令举例•绘制三维圆柱螺线t=linspace(0,10*pi);x=sin(t)y=cos(t)z=tplot3(sin(t),cos(t),t)-1-0.500.51-1-0.500.51010203040曲面彩色网线图•mesh(x,y,z)x(n维)向量,y(m维)向量,z(mXn维)矩阵绘制网线图•例：x=-7.5:0.5:7.5;y=x;[x,y]=meshgrid(x,y);r=sqrt(x.^2+y.^2)+eps;z=sin(r)./r;mesh(x,y,z)2222/)sin(yxyxz曲面彩色曲面图x=-7.5:0.5:7.5;y=x;[x,y]=meshgrid(x,y);r=sqrt(x.^2+y.^2)+eps;z=sin(r)./r;mesh(x,y,z)figuresurf(x,y,z)2.3MATLAB语言流程控制循环结构•for结构•while结构【例2-12】用循环求解【例2-13】用循环求解求最小的m【例2-14】求转移结构【例2-15】用循环求解求最大的m2.3.3开关结构2.4.1MATLAB语言的函数的基本结构•nargin，nargout,varargin,varargout一个有趣的问题房间里面有序号分别为1，2，…，100的电灯，初始时都关闭；又有100个序号分别为1，2，…，100的学生，当他们依次从这100盏灯前走过时分别将序号是本人序号的整数倍的所有灯开关切换一次（原来为开的变为关，原来为关的变为开）。请问当所有的人都走过之后，哪些灯是亮的？主要内容Matlab基本操作特征值问题积分变换线性常微分方程特征值问题(1)线性代数上所学的特征值问题广义特征值问题xAxBxAxd=eig(A)returnsavectoroftheeigenvaluesofmatrixA.d=eig(A,B)returnsavectorcontainingthegeneralizedeigenvalues,ifAandBaresquarematrices.[V,D]=eig(A)producesmatricesofeigenvalues(D)andeigenvectors(V)ofmatrixA,s