实验05基于微分方程对象建模及实现(一).pptVIP

实验要求 1、按时上下课，不迟到、不早退； 2、爱护设备，不作非法操作； 3、每次上课时带好教材和笔、纸； 3、上课期间遵守实验室纪律，不上网、不打游戏； 4、下课时摆好凳子，光好计算机； 5、每次实验课后以实验报告形式在下次实验课上课时上交作业。 注意 1、因计算机重启后，存储在C盘的文件将不会存在，故在编程、练习、做作业时最好在D盘建立一个目录，并在MATLAB中设置好路径。 2、每次下课时用可移动磁盘将自己所编的程序、所作的文件拷回去。 3、联系： E-mail:linzongbing@ QQ:8552178 Tel问题背景和实验目的 Matlab 解初值问题 dsolve 求解析解 dsolve 的使用 dsolve 举例 dsolve 举例 dsolve 举例 Euler 折线法 初值问题的Euler折线法 Euler 折线法举例 Euler 折线法源程序 解析解与用Euler 折线法求得的数值解作比较 K=dsolve(Dy=y+2*x/y^2,y(0)=1,x) X=0:0.01:2; Y=((5*exp(3*X))/3 - 2*X - 2/3).^(1/3); hold on plot(X,Y,b) hold off Euler折线法举例（续） Runge-Kutta 方法 Runge-Kutta 方法 四阶 R-K 方法源程序 解析解与用四阶龙格-库塔法求得的数值解作比较 K=dsolve(Dy=y+2*x/y^2,y(0)=1,x) X=0:0.1:2; Y=((5*exp(3*X))/3 - 2*X - 2/3).^(1/3); hold on plot(X,Y,b+) hold off Runge-Kutta 方法 Euler 法与 R-K法误差比较 Matlab提供的ODE求解器 数值求解举例 数值求解举例 Matlab 求解微分方程小结 思考与练习 1 车间空气的清洁 2 古尸年代鉴定问题 3 鸭子游动的轨迹 使用于精度较低的情形 单步法；2，3 阶 R-K 方法；累计截断误差为 (△x)3 非刚性 ode23 计算时间比 ode45 短 多步法；Adams算法；高低精度均可到 10-3～10-6 非刚性 ode113 适度刚性情形 采用梯形算法 适度刚性 ode23t 若 ode45 失效时，可尝试使用 多步法；Gear’s 反向数值微分；精度中等 刚性 ode15s 当精度较低时，计算时间比 ode15s 短 单步法；2 阶Rosebrock 算法；低精度 刚性 ode23s 当精度较低时，计算时间比ode15s短 梯形算法；低精度 刚性 ode23tb 说明 特点 ODE类型 求解器 大部分场合的首选方法 单步法；4，5 阶 R-K 方法；累计截断误差为 (△x)3 非刚性 ode45 fun=inline(-2*y+2*x^2+2*x,x,y); [x,y]=ode23(fun,[0,0.5],1); 注：也可以在 tspan 中指定对求解区间的分割，如： [x,y]=ode23(fun,[0:0.1:0.5],1); %此时 x=[0:0.1:0.5] 求初值问题 的数值解，求解范围为 [0,0.5] 例 7 ： 如果需求解的问题是高阶常微分方程，则需将其化为一阶常微分方程组，此时需用函数文件来定义该常微分方程组。 令 ，则原方程可化为 求解 Ver der Pol 初值问题 例 8： 先编写函数文件 verderpol.m function xprime=verderpol(t,x) global mu; xprime=[x(2); mu*(1-x(1)^2)*x(2) - x(1)]; 再编写脚本文件 vdpl.m，在命令窗口直接运行该文件。 clear; global mu;mu=7; y0=[1;0]; [t,x]=ode45(verderpol,[0,40],y0); plot(t,x(:,1),r-); 解: 令 y1=x，y2=y1’ 1、建立m-文件vdp1000.m如下： 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),-) 3、结果如图 解 1