matlab物理案例分析.ppt

《matlab物理案例分析.ppt》由会员分享，可在线阅读，更多相关《matlab物理案例分析.ppt（46页珍藏版）》请在课桌文档上搜索。

1、1,Matlab在物理中的应用案例研究,Matlab及其应用,2,一、Matlab 概述,MATLAB是美国Mathworks公司推出的一套高性能的数值分析和计算软件。最初的目的是为了解决矩阵的运算问题而开发的，经过不断完善，已成为将矩阵运算、数值分析、图形处理、编程技术结合在一起的功能强大的工具。它的优点在于快速开发计算方法，而不在于计算速度。,Matlab及其应用,3,Matlab的学习方法,必须做大量的练习，熟悉其中的函数联系和自己的课题相关，深刻的体会多看帮助文件，又一本好的参考书记住：Matlab可以做很多事情,Matlab及其应用,4,1.1 Matlab的出现,70年代中期，Cl

2、eve Moler和他的同事开发了LINPACK和EISPACK的Fortran子程序库 70年代末期，Cleve Moler 在新墨西哥大学给学生开线性代数，为学生编写了接口程序，这程序取名为MATLAB，即MATrix LABoratory1983年春天，工程师John Little与Moler、Steve Bangert一起开发了第二代专业版MATLAB1984年，MathWorks公司成立，MATLAB正是推向市场。,Matlab及其应用,5,1.2 Matlab的版本演化,Matlab 1.0Pc matlab-matlab 386Matlab3.5+simulinkMatlab 4

3、.0：simlink内嵌Matlab 5.0:全面的面向对象Matlab 5.15.3Matlab 6.0Matlab 6.5：购并了MATRIXxMatlab 7.0:2004,Matlab及其应用,6,1.3 Matlab的特点,1 高效方便的矩阵和数组运算2编程效率高（库函数非常丰富）3结构化、面向对象4方便的绘图工具5用户使用方便6功能强大的工具箱（核心部分中有数百个核心内部函数）总之，Matlab在信号处理、建模、系统识别、控制、优化等领域的问题处理上优势明显于其它高级语言。,Matlab及其应用,7,1.4 解决物理问题常用的软件,Mathmatica 数学计算功能。复杂符号计算，

4、主要在理论界MAGIC 电磁粒子仿真软件仿真软件：HFSS、CST等等,Matlab及其应用,8,二、MATLAB基础使用,2.1 MATLAB操作界面简介操作桌面的缺省外貌与先前版本相比，6.x版引入大量的交互工作界面，如通用操作界面、帮助界面、演示界面等指令窗、历史指令窗、工作空间浏览器和只有窗名的当前目录窗“开始”按钮(matlab6.5版本新加)通用操作界面,Matlab及其应用,9,通用交互界面介绍,内存数组编辑器(Array Editor)在工作空间浏览器中对变量操作时才启动交互界面分类目录窗（Lauch Pad)。可展开成树状结构，列着所有交互界面。M文件编辑/调试器（Edito

5、r/Debugger)当进行“打开文件”操作时才启动。帮助导航/浏览器(Help avigator/Browser)详尽展示由超文本写成的在线帮助。,Matlab及其应用,10,2.2 Matlab命令的两种执行方式,一.交互式的命令执行方式用户在命令窗口逐条输入命令，matlab逐条执行，这种方式操作简单直观，但速度慢，中间过程无法保留。通用操作界面二.M命令文件的程序设计方式将有关命令编成程序存储在一个文件(扩展名为.m)中，matlab自动依次执行，可调试复杂的程序(调试时只显示第一个错误)，是实际应用中主要的执行方式,Matlab及其应用,11,2.3 帮助的使用,MATLAB中有以下

6、几种方法获得帮助：帮助命令、帮助窗口、MATLAB帮助台、在线帮助页或用Web菜单直接链接到MathWorks公司的主页。,Matlab及其应用,12,怎样获得详尽的帮助信息,首先是MATLAB窗口中的help命令，它最常用、最快，往往也是最有效的。给出关键字就能找到相应的说明。Help命令简单易用，但是它要求准确给出关键字，如果记忆不准，就很难找到。这时就应该用第二种方案。Help Desk：这是MATLAB的HTML格式的帮助文档，有着很好的组织。可以在命令行键入helpdesk，进入“Matlab Help Desk”。Help Desk的搜索可以使用不完全关键字，这样即使遇到记忆模糊的

7、情况也可以很快查到。前面两种方法所得到的帮助信息常常是不够细致的，需要用第三种途径来补充，这就是电子版的使用手册。MATLAB里包含了大量关于MATLAB安装、使用、编程以及各种工具箱等的电子版手册，全部为PDF格式,具有搜索功能，因此也可以在指定的手册中搜索关键字，使用非常方便。,Matlab及其应用,13,2.4 Matlab变量,变量要求变量名、函数名区分字母的大小写。如SIN、Sin与sin是不同的对6.5以前的版本，变量名字符不超过31个，对6.5版本，变量名最多可包含63个字符。超出长度的字符被忽略变量必须以字母开头，之后可以是任意字母、数字或者下滑线，但不能是空格和标点Matla

8、b不支持汉字，汉字不能出现在变量名和文件名中,Matlab及其应用,14,特殊变量,Matlab及其应用,15,2.5 复数,不需要特殊的处理。复数可以表示为：a=10-9i a=10-9j%i=j=sqrt(-1)注意下面的例子中的表达式 a4=5+sin(.5)*i如果把sin(.5)*i中的*去掉，就会出错，sin(.5)i对MATLAB是无意义的，直接以字符i和j结尾，只适用于简单数值，不适用于表达式。复数运算不需要特殊处理，可以直接进行。另外real(z),imag(z),abs(z)和angle(z)等,Matlab及其应用,16,三、MATLAB的数组(矩阵)及其运算,3.1 数

9、组(矩阵)的创建 1.直接输入法整个输入数组必须以方括号“”为其首尾数组的行与行之间必须用”;”或者【Enter】隔离数组元素必须由逗号”,”或空格分隔2.一维数组的冒号生成法 x=a:inc:ba是数组的第一个元素，inc是间隔，b是最后一个元素a，inc，b之间的冒号须在英文状态下输入。中文中的冒号将导致matlab出错。inc缺省值为1，可以取正数或负数。,Matlab及其应用,17,语句生成矩阵,3.linspace命令a=linspace(1,10,5);%产生线性等间距的列向量 4.logspace命令 b=logspace(0,2,10);%产生对数等间距的列向量 5.矩阵连接c

10、=a b;6.数组生成函数 zeros ones eye rand randn,Matlab及其应用,18,3.2 矩阵运算,矩阵的转置用符号来表示和实现非共轭转置.矩阵的四则运算+-*/点乘：.*右除：./左除：.矩阵的乘方矩阵乘方.元素对元素的乘方,Matlab及其应用,19,四、流程控制,为了便于应用，MATLAB提供了一些流程控制的命令。这些命令对脚本编写带来了一些方便，但是需要注意的是，尽量不要使用这些流程控制命令，尤其是循环控制命令。1 for循环for循环允许一组命令以固定的次数重复，它的一般形式是for x=array commandendfor 和end之间的命令串按数组ar

11、ray的每一列执行一次，直到n次后终止。,Matlab及其应用,20,1 For循环应当注意的一些方面,for循环不能使用内部重新赋值循环变量而终止；for循环内部接受任何有效的MATLAB数组；for循环可以嵌套；只要有矩阵形式可以解决的问题，不要使用for循环。使用for循环的算法执行很慢，一个好的MATLAB算法部应当出现循环语句。循环可以使用break跳出，但只跳出所在的循环，不跳出整个嵌套结构。,Matlab及其应用,21,2 while循环,与for循环以固定的次数求一组指令相反，while循环以不定的次数求一组语句的值。While循环的一般形式为：while expression

12、 commondsend只要表达式expression里的所有元素为真，就执行命令串commands。通常表达式求值给一个标量值，单数组值也同样有效。,Matlab及其应用,22,3 if-else-end结构,很多情况下，命令的序列必须根据关系的检验有条件的执行，它由if-else-end结构提供。它的结构如下：if expression1 commands1else if expression2 commands2else if else if commandsend,Matlab及其应用,23,五、数据的可视化,二维图形 在二维曲线的绘制命令中，最重要、最基本的命令时plot，其他许多特

13、殊绘图命令都以它为基础形成。1 plot(x,s)x为实向量时，以该向量元素的下标为横坐标、元素值为纵坐标画出一条连续曲线。X为实矩阵时，则按列绘制每列元素值相对其下标的曲线，图中曲线数等于X阵列数。s用来指定线型。（正方形）,Matlab及其应用,24,2 plot(x,y,s)x,y为同维向量时，绘制以x、y元素为横、纵坐标的曲线。3 plot(x1,y1,s1,x2,y2,s2,)在此格式中，每个绘线“三元组”（x,y,s）的结构和作用，与plot(x,y,s)相同，不同的三元组之间没有约束。,Matlab及其应用,25,线型和标记,Plot命令中，在每一对数组数据后面，给plot一个附

14、加参量，就可以指定所要的颜色和线型。plot(x,y,rd,x,z,k-,x,y,b-.,x,z,co)b蓝色,c青色，g绿色，k黑色，m紫红色，r红色，w白色，y黄色 加格线和标注，grid on命令在当前图形的单位标记处加格线。grid off命令取消格线。用text命令可以在图形的任意指定位置增加标记和其他文本信息，格式是text(x,y,string).,Matlab及其应用,26,plot(x,y,x,z)grid title(sin(x)and cos(x)xlabel(independent variable x)ylabel(dependent variable y and z

15、)text(2,0.8.sin(x),Matlab及其应用,27,plot应用举例画一个圆,程序如下：R=2;t=0:pi/20:2*pi;x=R*sin(t);y=R*cos(t)plot(x,y,b:),grid,Matlab及其应用,28,例1无阻力抛射体的运动,设目标相对于射点的高度为，给定初速，试计算物体在真空中飞行的时间和距离。建模 无阻力抛射体的飞行是中学物理就解决了的问题，本题的不同点是目标和射点不在同一高度上，用MATLAB可使整个计算和绘图过程自动化。其好处是快速地计算物体在不同初速和射角下的飞行时间和距离。关键是在求落点时间 时，需要解一个二次线性代速方程。,Matlab

16、及其应用,29,MATLAB程序,clear;y0=0;x0=0;%初始位置vMag=input(输入初始速度(m/s):);%输入初始速度的大小和方向vDir=input(输入初速方向(度)：);yf=input(输入目标高度(m):);%输入目标高度vx0=vMag*cos(vDir*(pi/180);%计算x,y方向的初始速度vy0=vMag*sin(vDir*(pi/180);wy=-9.81;wx=0;%重力加速度（m/s2)tf=roots(wy/2,vy0,y0-yf);%解方程，计算落点时间tf=max(tf);%去掉庸解t=0:0.1:tf;y=y0+vy0*t+wy*t.2

17、/2;%计算轨迹x=x0+vx0*t+wx*t.2/2;xf=max(x),plot(x,y);%计算射程grid,Matlab及其应用,30,初始速度50m/s，初速方向为40度，目标高度为8m时的输入结果,Matlab及其应用,31,例2有空气阻力抛射体的运动,步骤一：给定小球质量m,重力加速度g，拖曳系数C和小球面积A（半径为r），根据D=rou*C*A/2求出D步骤二：给定小球的初始速度v（1）（大小和方向角theta），小球的位置（x（1）和y（1）。设定时间间隔detat步骤三：选定最大的间隔数N（或者最大的时间）,Matlab及其应用,32,求解步骤,步骤四：当nN或者t 时循环

18、步骤五到步骤九步骤五：计算加速度分量 步骤六：打印或者画,Matlab及其应用,33,步骤七：计算新的速度分量 步骤八：计算新的位移坐标 步骤九：时间上增加一个间隔 步骤十：结束,Matlab及其应用,34,Matlab程序,g=9.8;m=0.145;C=0.5;r=0.0366;A=pi*r2;rou=1.2;D=(rou*C*A)/2;detat=0.01;theta=(35/180*pi);x(1)=0;y(1)=0;V(1)=50;Vx(1)=V(1)*cos(theta);Vy(1)=V(1)*sin(theta);N=700;,Matlab及其应用,35,for n=1:N V(

19、n)=sqrt(Vx(n)2+Vy(n)2)ax(n)=-(D/m)*V(n)*Vx(n);ay(n)=-g-(D/m)*V(n)*Vy(n);Vx(n+1)=Vx(n)+ax(n)*detat;Vy(n+1)=Vy(n)+ay(n)*detat;x(n+1)=x(n)+Vx(n)*detat+0.5*ax(n)*detat2;y(n+1)=y(n)+Vy(n)*detat+0.5*ay(n)*detat2;px(n)=Vx(1)*n*detat;py(n)=Vy(1)*n*detat-0.5*g*(n*detat)2;endplot(x,y,r,px,py,g);gridxlabel(x),

20、ylabel(y);title(Projectile motion with air resistance);,Matlab及其应用,36,结果,Matlab及其应用,37,例3有空气阻力的落体运动。,一个质量为75kg的人在某一高度静止下落，空气阻力F与v的关系满足F=-k，其中常数k=0.22kg/m，根据牛顿第二定律知 若令g=9.80m/,求下落过程中速度v(t)和位移x(t)的表达式，并画图。,Matlab及其应用,38,推导,Matlab及其应用,39,Matlab程序,v=dsolve(m*D2v+2*k*v*Dv=0,v(0)=0,Dv(0)=g,t)x=dsolve(m*g-

21、k*(Dx)2=m*D2x,x(0)=0,Dx(0)=0,D2x(0)=g,t)m=75;k=0.22;g=9.8t=0:0.5:30;v=tanh(t.*(g*m*k)(1/2)/m)*(g*m*k)(1/2)/k%x=(m/k)*log(cosh(t.*(k*g/m)1/2)plot(t,v),gtext(v(t)hold onplot(t,x),gtext(x(t),Matlab及其应用,40,结果,Matlab及其应用,41,例4小球碰撞问题,质量为m的小球以速度 正面碰撞质量为M的静止小球，假设碰撞是完全弹性的，即没有能量损失，求碰撞后两球的速度，及它们与两球质量比KM/m的关系。建

22、模 设碰撞后两球速度都与 同向，球m的速度为u，球M的速度为v，列出动量守恒和能量守恒方程，引入质量比K=M/m和无量纲速度ur=u/，vr=v/后，有,Matlab及其应用,42,推导过程,动量守恒动能守恒化为,（1）,（2）,（3）,（4）,Matlab及其应用,43,由（3）代入（4）展开并整理多项式(6)得可用roots命令求根,（5）,（6）,Matlab及其应用,44,MATLAB程序,clearK=logspace(-1,1,11);%设自变量数组K，从K=0.1 10，按等比取11个点for i=1:length(K)%对各个K循环计算ur1=roots(1+1/K(i),-2/K(i),(1/K(i)-1);%求两个解ur(i)=ur1(abs(ur1-1)0.001);%去掉在1邻近的庸解endvr=(1-ur)./K;%用(5)式求vr，用元素群运算K;ur;vr%显示输出数据semilogx(K,ur,K,vr),grid,gtext(ur),gtext(vr)%绘图,Matlab及其应用,45,输出结果,Matlab及其应用,46,结果的物理分析,从上图中可以看出当K1时，为负，即当静止球质量大于主动球质量时，主动球将反弹。当K=1时，=0，即主动球将全部动能传给静止球。当K1时，为正，说明主动球将继续沿原来方向运动。,