常微分方程的求解与定性分析实验报告
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
常微分方程的求解与定性分析实验报告
一、实验综述
1、实验目的及要求
●归纳和学习求解常微分方程(组)的基本原理和方法;
●掌握解析、数值解法,并学会用图形观察解的形态和进行解的定性分析;
●熟悉MATLAB软件关于微分方程求解的各种命令;
●通过范例学习建立微分方程方面的数学模型以及求解全过程;
●通过该实验的学习,使学生掌握微分方程(组)求解方法(解析法、欧拉法、
梯度法、改进欧拉法等),对常微分方程的数值解法有一个初步了解,同时学会使用MATLAB软件求解微分方程的基本命令,学会建立微分方程方面的数学模型。这对于学生深入理解微分、积分的数学概念,掌握数学的分析思维方法,熟悉处理大量的工程计算问题的方法是十分必要的。
2、实验仪器、设备或软件
电脑、matlab7.0
二、实验过程(实验步骤、记录、数据、分析)
实验内容:
根据实验内容和步骤,完成以下实验,要求写出实验报告(实验目的→问题→数学模型→算法与编程→计算结果→分析、检验和结论)
1.求微分方程的解析解,并画出它们的图形。
y '= y + 2 x, y (0) = 1, 0< x <1;
m=dsolve('Dy=y+2*x','y(0)=1','x')
ezplot(m,[0 1])
m =
3*exp(x) - 2*x – 2
0.1
0.2
0.3
0.4
0.50.6
0.7
0.8
0.9
1
11.522.533.5
4x
3 exp(x) - 2 x - 2
1.求微分方程⎪⎩
⎪
⎨⎧====-+]100[0)0(;0)0(01.03t u
u u u u 的数值解,要求编写求解程序。 function dy=vdp1000(t,y) dy=zeros(2,1); dy(1)=y(2);
dy(2)=-y(1)+0.1*y(1)^3;
[T,Y]=ode15s('vdp1000',[0 10],[0 0]); plot(T,Y(:,1),'-')
-1-0.8-0.6-0.4-0.200.20.40.60.81
3.Rossler 微分方程组:
⎪⎩
⎪⎨⎧-+=+=--=)(''
'c x z b z ay
x y z y x 当固定参数b =2,c =4时,试讨论随参数a 由小到大变化(如 a ∈(0,0.65))而方程解的变化情况,并且画出空间曲线图形,观察空间曲线是否形成混沌状? function r=rossler(t,x) global a; global b; global c;
r=[-x(2)-x(3);x(1)+a*x(2);b+x(3)*(x(1)-c)];
global a; global b; global c; b=2; c=4;
t0=[0,200];
for a=0:0.1:0.6
[t,x]=ode45('rossler',t0,[0,0,0]); subplot(1,2,1);
plot(t,x(:,1),'r',t,x(:,2),'g',t,x(:,3),'b');
title('x(红色),y(绿色),z(蓝色)随t 的变化情况');xlabel('t');
subplot(1,2,2);
plot3(x(:,1),x(:,2),x(:,3))
title('相图');xlabel('x');ylabel('y');zlabel('z'); pause end 结果显示: a=0:
x(红色),y(绿色),z(篮色)随t 变化情况t
相图
z
a=0.1:
-1
-0.8-0.6-0.4-0.200.2
0.4
0.6x(红色),y(绿色),z(篮色)随t 变化情况t
相图
z
a=0.2:
50
100150
200
-4
-3-2-1012
3
4x(红色),y(绿色),z(篮色)随t 变化情况t
-5
2
相图
z
a=0.3:
50
100150
200
-5-4-3-2-1012345x(红色),y(绿色),z(篮色)随t 变化情况t
-5
相图
z
a=0.4:
050
100150200
-6
-4
-2
2
4
6
x(红色),y(绿色),z(篮色)随t 变化情况t
-5
相图
z
a=0.5:
50
100150
200
-8-6-4-2024681012x(红色),y(绿色),z(篮色)随t 变化情况t
-5
相图
z
结果分析:从图像可以看出,当a=0时,微分方程的解(x,y,z)收敛与(0,0.5,0.5);当a=0.1时,(x,y,z)仍收敛与(0,0.5,0.5),只是收敛速度减慢;当a=0.2时,(x,y,z)已发散,周期性变化;随着a的增大,(x,y,z)接近其极限环的速度加快,空间曲线成混沌状。
4.炮弹发射角的确定
炮弹发射视为斜抛运动,已知初始速度为200m/s,若要击中水平距离360m、在垂直距离160m的目标,当忽略空气阻力时,发射角应为多大?此时炮弹的运行轨迹如何?
要求:
(1) 建立在忽略空气阻力情况下的描述炮弹发射轨迹的数学模型;
(2) 用Matlab 软件求解方程和微分方程;
(3) 结合实际对解的合理性进行分析。
进一步思考:
如果要考虑水平方向的阻力,且设阻力与(水平方向)速度成正比,系数为0.1(1/s ),结果又如何?此时炮弹的运行轨迹如何?
解:
(1)忽略空气阻力时,设发射角为a,炮弹的飞行时间为t
水平方向:Vx=V0*cos a=200*cos(a)
竖直方向:Vy=V0*sin a=200*sin(a)
得到:x=Vx*t=200*cos(a)*t=360
y=Vy*t-1/2*g*t^2=200*sin(a)*t-1/2*9.8*t^2=160
得到炮弹的路程为Y=360*tan(a)-4.9*(360/200/cos(a)).^2-160
编程为:
function Y=fun2(a)
Y=360*tan(a)-4.9*(360/200/cos(a)).^2-160;
function Y=fun3(a0,a1,n,tol)
a(1)=a0;
a(2)=a1;
b=1;
i=2;
while(abs(b)>eps*a(i))
a(i+1)=a(i)-fun2(a(i))*(a(i)-a(i-1))/(fun2(a(i))-fun2(a(i-1)));
b=a(i+1)-a(i);
i=i+1;
if(i>n) error('n is full');
end
end