实验二微分方程差分方程模型Matlab求解 下载本文

内容发布更新时间 : 2024/3/29 7:28:57星期一 下面是文章的全部内容请认真阅读。

资料

实验二: 微分方程与差分方程模型Matlab求解

一、实验目的

[1] 掌握解析、数值解法,并学会用图形观察解的形态和进行解的定性分析; [2] 熟悉MATLAB软件关于微分方程求解的各种命令;

[3] 通过范例学习建立微分方程方面的数学模型以及求解全过程; [4] 熟悉离散 Logistic模型的求解与混沌的产生过程。

二、实验原理

1. 微分方程模型与MATLAB求解

解析解

用MATLAB命令dsolve(‘eqn1’,’eqn2’, ...) 求常微分方程(组)的解析解。其中‘eqni'表示第i个微分方程,Dny表示y的n阶导数,默认的自变量为t。 (1) 微分方程 例1 求解一阶微分方程 (1) 求通解 输入:

dsolve('Dy=1+y^2')

dy?1?y2 dx输出:

ans =

tan(t+C1)

(2)求特解 输入:

dsolve('Dy=1+y^2','y(0)=1','x')

指定初值为1,自变量为x 输出:

ans =

tan(x+1/4*pi)

.

资料

1x2y???xy??(x2?)y?04例2 求解二阶微分方程 y(?/2)?2

y?(?/2)??2/?原方程两边都除以x2,得y??? 输入:

dsolve('D2y+(1/x)*Dy+(1-1/4/x^2)*y=0','y(pi/2)=2,Dy(pi/2)=-2/pi','x')

ans =

- (exp(x*i)*(pi/2)^(1/2)*i)/x^(1/2) +

(exp(x*i)*exp(-x*2*i)*(pi/2)^(3/2)*2*i)/(pi*x^(1/2))

试试能不用用simplify函数化简 输入: simplify(ans)

ans =

2^(1/2)*pi^(1/2)/x^(1/2)*sin(x) (2)微分方程组

例3 求解 df/dx=3f+4g; dg/dx=-4f+3g。

(1)通解:

[f,g]=dsolve('Df=3*f+4*g','Dg=-4*f+3*g') f =

exp(3*t)*(C1*sin(4*t)+C2*cos(4*t)) g =

exp(3*t)*(C1*cos(4*t)-C2*sin(4*t))

11y??(1?2)y?0 x4x特解:

[f,g]=dsolve('Df=3*f+4*g','Dg=-4*f+3*g','f(0)=0,g(0)=1')

f =

exp(3*t)*sin(4*t)

.

资料

g =

exp(3*t)*cos(4*t)

数值解

在微分方程(组)难以获得解析解的情况下,可以用Matlab方便地求出数值解。格式为:

[t,y] = ode23('F',ts,y0,options)

注意:

? 微分方程的形式:y' = F(t, y),t为自变量,y为因变量(可以是多

个,如微分方程组);

? [t, y]为输出矩阵,分别表示自变量和因变量的取值;

? F代表一阶微分方程组的函数名(m文件,必须返回一个列向量,每个

元素对应每个方程的右端);

? ts的取法有几种,(1)ts=[t0, tf] 表示自变量的取值范围,(2)

ts=[t0,t1,t2,…,tf],则输出在指定时刻t0,t1,t2,…,tf处给出,(3)ts=t0:k:tf,则输出在区间[t0,tf]的等分点给出; ? y0为初值条件;

? options用于设定误差限(缺省是设定相对误差是10^(-3),绝对误差

是10^(-6)); ode23是微分方程组数值解的低阶方法,ode45为中阶方法,与ode23类似。

例4 求解一个经典的范得波(Van Der pol)微分方程:

u''?(u2?1)u'?u?0,u(0)?1,u'(0)?0

解 形式转化:令y1?u(t);y2?u'(t)。则以上方程转化一阶微分方程组:

??y2;y2??(1?y12)y2?y1。 y1编写M文件如下,必须是M文件表示微分方程组,并保存,一般地,M文件的名字与函数名相同,保存位置可以为默认的work子目录,也可以保存在自定义文件夹,这时注意要增加搜索路径(File\\Set Path\\Add Folder)

function dot1=vdpol(t,y); dot1=[y(2); (1-y(1)^2)*y(2)-y(1)]; 在命令窗口写如下命令:

.