matlab求解常微分方程

用matlab 求解常微分方程

在MATLAB 中，由函数dsolve ()解决常微分方程（组）的求解问题，其具体格式如下：

r = dsolve('eq1,eq2,...', 'cond1,cond2,...', 'v')

'eq1,eq2,...'为微分方程或微分方程组，'cond1,cond2,...',是初始条件或边界条件，'v'是独立变量，默认的独立变量是't'。

函数dsolve 用来解符号常微分方程、方程组，如果没有初始条件，则求出通解，如果有初始条件，则求出特解。

例1：求解常微分方程1dy dx x y =

+的MATLAB 程序为：dsolve('Dy=1/(x+y)','x')

，

注意，系统缺省的自变量为t ，因此这里要把自变量写明。 其中：Y=lambertw(X)表示函数关系Y*exp(Y)=X 。

例2：求解常微分方程的MATLAB 程序为：

2

'''0yy y −=Y2=dsolve('y*D2y-Dy^2=0','x')

Y2=dsolve('D2y*y-Dy^2=0','x')

我们看到有两个解，其中一个是常数0。

例3：求常微分方程组

2

5

3

t

t

dx

x y e

dt

dy

x y e

dt

⎧

++=

⎪⎪

⎨

⎪−−=

⎪⎩通解的MATLAB程序为：

[X,Y]=dsolve('Dx+5*x+y=exp(t),Dy-x-3*y=exp(2*t)','t')

例4：求常微分方程组

2

210cos,

24,

t

t

t

dx dy

x t x

dt dt

dx dy

y e y

dt dt

=

−

=

⎧

+−==

⎪⎪

⎨

⎪++==

⎪⎩

2

通解的MATLAB程序为：

[X,Y]=dsolve('Dx+2*x-Dy=10*cos(t),Dx+Dy+2*y=4*exp(-

2*t)','x(0)=2,y(0)=0','t')

以上这些都是常微分方程的精确解法，也称为常微分方程的符号解。但是，我们知道，有大量的常微分方程虽然从理论上讲，其解是存在的，但我们却无法求出其解析解，此时，我们需要寻求方程的数值解，在求常微分方程数值解方面，MATLAB具有丰富的函数，我们将其统称为solver，其一般格式为：

[T,Y]=solver(odefun,tspan,y0)

该函数表示在区间tspan=[t0,tf]上，用初始条件y0求解显式常微分方程。 '(,y f t y =)solver 为命令ode45，ode23，ode113，ode15s ，ode23s ，ode23t ，ode23tb 之一，这些命令各有特点。我们列表说明如下：

求解器 特点

说明

ode45

一步算法，4,5阶Runge-Kutta

方法累积截断误差

3

()x Δ大部分场合的首选

算法

ode23

一步算法，2,3阶Runge-Kutta

方法累积截断误差3

()

x Δ使用于精度较低的

情形 ode113多步法，Adams 算法，

高低精度均可达到

36

10~10−−计算时间比ode45

短 ode23t 采用梯形算法 适度刚性情形 ode15s

多步法，Gear’s 反向 数值积分，精度中等 若ode45失效时， 可尝试使用 ode23s

一步法，2阶Rosebrock

算法， 低精度。

当精度较低时， 计算时间比ode15s

短

odefun 为显式常微分方程中的'(,y f t y =)(,)

f t y tspan 为求解区间，要获得问题在其他指定点上的解，则令012,,,t t t 012[,,,,]f tspan t t t t = （要求i t 单调递增或递减），y0初始条件。

例5：求解常微分方程2

'222y y x =−++x .5，00x ≤≤，(0)1y =的MATLAB 程序如下：

y=dsolve('Dy=-2*y+2*x^2+2*x','y(0)=1','x') x=0:0.01:0.5;

yy=subs(y,x);

fun=inline('-2*y+2*x*x+2*x');[x,y]=ode15s(fun,[0:0.01:0.5],1);ys=x.*x+exp(-2*x);

plot(x,y,'r',x,ys,'b')

例6：求解常微分方程222

(1)0,(0)1,'(0)0d y dy y y y y dt

dt μ−−+===的解，并画出解的图形。 分析：这是一个二阶非线性方程(函数以及所有偏导数军委一次幂的是现性方程，高于一次的为非线性方程)，用现成的方法均不能求解，但我们可以通过下面的变换，将二阶方程化为一阶方程组，即可求解。

令：1x y =，

2dy

x dt =

，7μ=，则得到：

1

21221212

,(0)17(1),(0)0

dx x x dt dx x x x x dt ⎧==⎪⎪⎨

⎪=−−=⎪⎩

解：