发布时间 : 星期一 文章Matlab讲义 - 图文更新完毕开始阅读fb1ad585a0116c175f0e48cd
4．彗星轨迹图
彗星轨迹图是最简单的动画，可以展现质点运动的过程，其使用格式为： comet(x,y) comet(x,y,p) comet3(x,y,z) comet3(x,y,z,p) 二维平面彗星轨迹，x,y的含义与plot命令中x,y含义相同，绘出的轨迹慧长为p*length(x)，0 ?n?1?设空气阻力与抛体运动速度的n次方成正比，F??bvv，b为正值常量．其运动微分方程为 ????b(x?2?y?2)(n?1)2x?，my????b(x?2?y?2)(n?1)2y? mx令 ?，y(3)?y，y(4)?y?，则化为 y(1)?x，y(2)?xby(2)[y(2)2?y(4)2](n?1)2 mb?(4)??9.8?y(4)[y(2)2?y(4)2](n?1)2 ?(3)?y(4)yym?(2)???(1)?y(2)，yy编写微分方程函数文件znxpfun.m如图39；再编写解微分方程的主程序znxp.m如图41；运行结果如图40． 图39 图40 图41 图42 注意例子中第4，5，6行中的命令．comet语句在例6图26中就已经用过了，从现在的图上看不出动态效果， 13 读者自己编程运行后就可以看到了． 把例9图35中第7行plot3改为comet3(x,y,z,p)就是三维彗星轨迹图了，如图42．程序中解微分方程的步长由0.1改为0.002，是为了使彗星轨迹图运行的速度减慢．由于不同计算机计算速度不同，可调整步长使图形运行速度适当． §3 一维非线性振动 一、一维振动 一维振动一般微分方程形式为 ???FR(x?)?F(x)?F0cos?t (3.1) mx二、一维线性振动 在恢复力和阻力满足线性关系条件下，一维线性振动的微分方程为 ????x??kx?F0cos?t mx引入2???m，?0?km，f?F0m，则 2???2?x???02x?fcos?t (3.2) x下面用数值计算的方法对一维线性振动做一点讨论．主要是通过x?t图和x?运动情况． 依据(3.2)式编写微分方程函数文件ywzdfun.m如图43；再编写解微分方程的主程序ywdz.m如图44．在主程序中改变参量，⑴令b=0，w0=1，f=0，w=2；运行，得简谐振动的x?t图和相图，如图45和fl.46所示．⑵令b=0.1，w0=1，f=0，w=2；运行，得弱阻尼振动的x?t图和相图，如图47和fl.48所示．⑶令b=0.1，w0=1，f=0.5，w=2；运行，得受迫振动的x?t图和相图，如图49和fl.50所示；再令w=0.99；运行，得共振状态的x?t图和相图，如图51和fl.52所示． ?图（即相图）x，去了解质点的 图43 图44 图45 图46 图47 图48 图49 图50 图51 图52 读者还可以继续改变参量，比如对阻尼振动情况，改变b的值就可以得到临界阻尼和过阻尼相应的图形． 14 三、用数值计算和相图研究大幅度单摆的运动 本节研究最大摆角不受限制的单摆的运动, 单摆由质点和刚性轻杆构成，它的运动微分方程为 ????2sin??0 (3.3) ?0这个非线性方程是可积的，它的解可用第一类椭圆积分表示；现在，用数值计算和相图方法研究它． 依据(3.3)式编写微分方程函数文件djdbfun.m如图53；再编写解微分方程并画x?t图的主程序djdb1.m如图54；运行，计算结果如图55． 图53 图54 图55 图56 图57 图58 编写解微分方程并画相图的主程序djdb2.m如图56；运行，计算结果如图57． 最后，编写解微分方程，作快速傅里叶变换并画功率谱的主程序djdb3.m如图58；运行，计算结果如图59． 由图55可见，振幅（最大摆角）越大，曲线与简谐振动正弦曲线的差别越显著．同时，周期随振幅增加而变大，当振幅接近?时，周期迅速增大（实际上，振幅??时，周期??）；图60为由数值计算画出的，振幅不太大时的“周期-振幅”曲线． 由图57可见，相轨迹分为两类，它们之间存在相应于总能量E?mgl的 分界线（以悬点为重力势能零点）；分界线之内的相轨迹是围绕中心的闭合曲线，相应于周期性的往复摆动，总能量满足0?E?mgl；分界线之外的相轨迹是 不闭合的，它代表另一种周期运动，即向一个方向的转动，相应总能量满足 E?mgl． 对应于振幅较大（即总能量较大）的往复摆动，谐频的影响亦不容忽视，它 的相轨迹虽然闭合，但已不是椭圆，与简谐振动的相图是不同的． 图59 15 图59为对数值解进行快速傅里叶变换（fft）得到的功率谱图（纵轴为对数坐标），清晰可见相应?的基频谱线，相应3?和5?的谐频谱线；隐约可见相应7?的谐频谱线．由于谐频振动的振幅随其级数的增加而迅速减少，所以更高倍的谐频谱线无法显示． 四、自激振动 范·德·波尔方程 ????(x02?x2)x???02x?0 (3.4) 图60 x依据(3.4)式编写微分方程函数文件zjzdfun.m如图61；再编写解微分方程并画相图的主程序zjzd.m如图62，参量取值为??0.3，x0?1，?0?1；运行，计算结果如图63?65． 图61 图62 图63 图64 图65 §4 能导致混沌的倒摆受迫振动 一、运动微分方程的建立 实验装置如图66，图67为示意原理图图67中，R为可绕O轴转动的均匀铝制圆盘，m为与盘固定的配重，并用m处指针标志盘的位置，总体上看是一个m在O轴上方的摆，故称为倒摆．螺旋形弹簧提供弹性回复力；D为电磁铁，提供线性阻尼；可调速电机E带动轮B旋转，轮B又经连杆L驱动弹簧A端，提供驱动力． 图66 图67 如图67，建立直角坐标系 Oxyz，以 ?描述系统位置．系统所受对 z轴的力矩有: 重力矩 ?，螺旋形弹簧线性恢复力矩?c?rmgsin??asin?，阻力矩?b?矩M，由于卷弹簧 A端被驱动而受周期性驱动力 cos?t．无驱动时系统有三个平衡位置，??0为不稳定平衡位置，其两侧还有O1和O2两个稳定平衡位置． 系统对z轴的转动惯量为I，由对z轴的角动量定理，得 ???asin??b???c??Mcos?t I? 16