




版权说明:本文档由用户提供并上传,收益归属内容提供方,若内容存在侵权,请进行举报或认领
文档简介
1、用 Matlab 求解微分方程,借助 Matlab 软件,可以方便地求出微分方程(组)的解析解和数值解。,微分方程(组)的解析解,求微分方程(组)解析解的命令为 dsolve(eqn1, eqn2, ., x) 其中“eqni”表示第 i 个方程,“x”表示微分方程(组)中的自变量,默认时自变量为 t。此外,在“eqni”表示的方程式中,用 D 表示求微分,D2、D3 等表示求高阶微分,任何 D 后所跟的字母表示因变量。,例 8.5.1 求解一阶微分方程 dy/dx = 1 + y2。,求通解 输入:dsolve(Dy=1+y2, x) 输出:ans = tan(x+C1) 求特解 输入:ds
2、olve(Dy=1+y2, y(0)=1, x) 输出:ans = tan(x+1/4*pi),例 8.5.2 求解下列微分方程的通解及 y(0) = 0 和 y (0) = 15 条件下的特解,求通解 输入:y=dsolve(D2y+4*Dy+29*y=0, x) 输出:y = C1*exp(-2*x)*sin(5*x)+C2*exp(-2*x)*cos(5*x) 求特解 输入:y=dsolve(D2y+4*Dy+29*y=0, y(0)=0, Dy(0)=15, x) 输出:y = 3*exp(-2*x)*sin(5*x),例 8.5.3 求解下列微分方程组,求通解 方式一 输入: x,y
3、,z=dsolve(Dx=2*x-3*y+3*z,Dy=4*x -5*y+3*z,Dz=4*x-4*y+2*z, t); 输出:x = C2*exp(-t)+C3*exp(2*t) y = C2*exp(-t)+C3*exp(2*t)+exp(-2*t)*C1 z = C3*exp(2*t)+exp(-2*t)*C1,方式二 输入:x,y,z=dsolve(Dx=2*x-3*y+3*z,Dy=4*x -5*y+3*z,Dz=4*x-4*y+2*z, t); x=simple(x) % 将x化简 y=simple(y) z=simple(z) 输出:x = C2/exp(t)+C3*exp(t)
4、2 y = C2*exp(-t)+C3*exp(2*t)+exp(-2*t)*C1 z = C3*exp(2*t)+exp(-2*t)*C1,求特解 输入:x,y,z=dsolve(Dx=2*x-3*y+3*z, Dy=4*x-5*y+3*z,Dz=4*x-4*y+2*z, x(0)=0, y(0)=1, z(0)=2, t); x=simple(x) % 将x化简 y=simple(y) z=simple(z) 输出:x = exp(2*t)-exp(-t) y = exp(2*t)-exp(-t)+exp(-2*t) z = exp(2*t)+exp(-2*t),微分方程(组)的数值解,事
5、实上,能够求得解析解的微分方程或微分方程组少之又少,多数情况下需要求出微分方程(组)的数值解。 Matlab中求微分方程数值解的函数有五个:ode45,ode23,ode113,ode15s,ode23s。调用格式为 t, x = solver (f, ts, x0, options),需要特别注意的是: solver 可以取以上五个函数之一,不同的函数代表不同的内部算法:ode23 运用组合的 2/3 阶龙格库塔费尔贝算法,ode45 运用组合的 4/5 阶龙格库塔费尔贝算法。通常使用函数 ode45; f 是由待解方程写成的m文件的文件名; ts=t0, tf,t0、tf为自变量的初值和终
6、值; x0为函数的初值;, options 用于设定误差限(可以缺省,缺省时设定为相对误差 103,绝对误差 106),程序为 options = odeset(reltol, rt, abstol, at) 其中rt和at分别为设定的相对误差和绝对误差; 在解 n 个未知函数的方程组时,x0、x 均为 n 维向量,m 文件中待解方程组应以 x 的分量形式写成; 使用 Matlab 软件求数值解时,高阶微分方程必须等价地变换成一阶微分方程组。,例 8.5.4 求解下列微分方程,解:令 y1 = x,y2 = y1,则微分方程变为一阶微分方程组:,(1) 建立 m 文件 vdp1000.m 如下
7、: function dy=vdp1000(t,y) dy=zeros(2,1); dy(1)=y(2); dy(2)=1000*(1-y(1)2)*y(2)-y(1); (2) 取 t0=0,tf=3000,输入命令: T,Y=ode15s(vdp1000,0 3000,2 0); plot(T,Y(:,1),-) 运行程序,得到如图的结果。,例 8.5.5 求解下列微分方程组,(1) 建立 m 文件 rigid.m 如下: function dy=rigid(t,y) dy=zeros(3,1); dy(1)=y(2)*y(3); dy(2)=-y(1)*y(3); dy(3)=-0.51
8、*y(1)*y(2); (2) 取 t0=0,tf=12,输入命令: T,Y=ode45(rigid,0 12,0 1 1); plot(T,Y(:,1),-,T,Y(:,2),*,T,Y(:,3),+),运行程序,得到如图的结果。图中,y1 的图形为实线,y2 的图形为“*”线,y3 的图形为“+”线。,例 8.5.6 导弹追踪问题 设位于坐标原点的甲舰向位于 x 轴上点 A(1, 0) 处的乙舰发射导弹,导弹头始终对准乙舰。如果乙舰以最大的速度 v0(是常数)沿平行于 y 轴的直线行驶,导弹的速度是 5v0,求导弹运行的曲线方程。又乙舰行驶多远时,导弹将它击中?,解:如图所示,假设导弹在
9、t 时刻的位置为P(x(t), y(t),乙舰位于 Q(1, v0t)。由于导弹头始终对准乙舰,故此时直线 PQ 就是导弹的轨迹曲线弧 OP 在点 P 处的切线,,于是有,即,又根据题意,弧 OP 的长度为 |AQ| 的 5 倍,于是,消去 t,得到导弹追踪模型如下:,下面求解这个初值问题。,解法一 解析解 利用微分方程初值问题的解析解法,得导弹的运行轨迹为:,参见下图。,根据题意,乙舰始终沿平行于 y 轴的直线 x = 1 行驶,且由上式知,当 x = 1 时 y = 5/24,故当乙舰航行到点 (1, 5/24) 处时被导弹击中。同时可求得被击中时间为:t = y/v0 = 5/24v0;若 v0 = 1,则在 t = 0.21 处被击中。,解法二 数值解 令 y1 = y,y2 = y1,将先前给出的导弹追踪模型化为一阶微分方程组,(1) 建立 m 文件 eq1.m function dy=eq1(x,y) dy=zeros(2,1); dy(1)=y(2); dy(2)=1/5*sqrt(1+y(1)2)/(1-x); (2) 取 x0=0,xf=0.9999,建立主
温馨提示
- 1. 本站所有资源如无特殊说明,都需要本地电脑安装OFFICE2007和PDF阅读器。图纸软件为CAD,CAXA,PROE,UG,SolidWorks等.压缩文件请下载最新的WinRAR软件解压。
- 2. 本站的文档不包含任何第三方提供的附件图纸等,如果需要附件,请联系上传者。文件的所有权益归上传用户所有。
- 3. 本站RAR压缩包中若带图纸,网页内容里面会有图纸预览,若没有图纸预览就没有图纸。
- 4. 未经权益所有人同意不得将文件中的内容挪作商业或盈利用途。
- 5. 人人文库网仅提供信息存储空间,仅对用户上传内容的表现方式做保护处理,对用户上传分享的文档内容本身不做任何修改或编辑,并不能对任何下载内容负责。
- 6. 下载文件中如有侵权或不适当内容,请与我们联系,我们立即纠正。
- 7. 本站不保证下载资源的准确性、安全性和完整性, 同时也不承担用户因使用这些下载资源对自己和他人造成任何形式的伤害或损失。
最新文档
- 2025至2030年中国鸳鸯梅行业发展研究报告
- 2025至2030年中国高毛彩带浴室地垫行业发展研究报告
- 2025至2030年中国链霉素针行业发展研究报告001
- 2025至2030年中国铁制佛帽行业发展研究报告
- 一年级道德与法治下册 第一单元 我的生活好习惯 第2课《爱清洁 讲卫生》教学设计1 苏教版
- 新闻文体知识2
- 小孩子防拐骗教育
- 2025至2030年中国精密周转盘行业投资前景及策略咨询报告001
- 2025至2030年中国空压高速夹头行业发展研究报告001
- 1WeallownEnglish-英语九年级下册
- 轮机工程船舶柴油机主动力推进装置培训课件
- 森林区划-组织森林经营类型(森林资源经营管理)
- 灰枣种植技术
- 医院污水处理技术指南(环发2023年197号2023年20实施)
- 求职登记表(标准模版)
- 高中历史 统编版 必修 《中外历史纲要》(下) 多元与交融:文明网络的初构《古代文明的产生与发展》第2课时 古代世界的帝国与文明的交流 课件
- 智慧景区建设实务 智慧景区视频监控系统
- 简单实用电子简历表格,个人简历模板word格式
- 备考2023新高考英语听力12(答案听力原文)
- GB/T 23724.1-2016起重机检查第1部分:总则
- GB 2707-2016食品安全国家标准鲜(冻)畜、禽产品
评论
0/150
提交评论