微分方程建模方法

返 回
对 于 已 给 的 精 确 度 0 ，当 满 足 y i 1 然 后 继 续 下 一 步 yi 2 的 计 算 ．
( k 1)
y i 1 时 ， y i 1 y i 1 ， 取
(k )
（ k 1）
此即改进的欧拉法．
3．使用泰勒公式 以此方法为基础，有龙格-库塔法、线性多步法等方 法． 4．数值公式的精度 当一个数值公式的截断误差可表示为O（hk+1）（其 中k为正整数，h为步长）时，称它是一个k阶公式． k越大，则数值公式的精度越高．
例1 求
解
du dt
d y dx
2
2
0 应表达为：D2y=0.
2
1 u
的通解.
输入命令：dsolve('Du=1+u^2','t')
结 果：u = tg(tc)
To MATLAB（ff1）
例 2 求微分方程的特解.
d2 y dy 4 29 y 0 2 dx dx y (0 ) 0 , y '(0 ) 1 5 解 输入命令: y=dsolve('D2y+4*Dy+29*y=0','y(0)=0,Dy(0)=15','x')
•欧拉法是一阶公式，改进的欧拉法是二阶公式．
•龙格-库塔法有二阶公式和四阶公式． •线性多步法有四阶亚当斯外插公式和内插公式． 返 回
（三）用MATLAB软件求常微分方程,ts,x0,options)
自变 量值 函数 值
ode45 ode23 ode11 3ode1 5sode 23s
注意:
1．在解含n个未知数的方程组时，x0和x均为n维向量， M文件中的待解方程组应以x的分量形式写出． 2．使用MATLAB软件求数值解时，高阶微分方程必须 等价地变换成一阶微分方程组．
例4
解: 令 y1=x，y2=y1’
d2x dx 2 x 0 2 1 0 0 0 (1 x ) dt dt x (0 ) 2; x '(0 ) 0
To MATLAB（ff2）
结 果 为 : y =3e-2xsin（5x）
例 3 求微分方程组的通解.
dx 2x 3 y 3z dt dy 4x 5 y 3z dt dz 4x 4y 2z dt
解 输入命令 ：
[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) To MATLAB（ff3） z=simple(z)
结 果 为：x = (c1-c2+c3+c2e -3t-c3e-3t)e2t y = -c1e-4t+c2e-4t+c2e-3t-c3e-3t+c1-c2+c3)e2t z = (-c1e-4t+c2e-4t+c1-c2+c3)e2t 返 回
微分方程的数值解
（一）常微分方程数值解的定义 在生产和科研中所处理的微分方程往往很复杂,且大 多得不出一般解（解析解）．而实际中的对初值问题， 一般是要求得到解在若干个点上满足规定精确度的近似 值，或者得到一个满足精确度要求的便于计算的表达 式． 因此，研究常微分方程的数值解法是十分必要的．
微 分 方 程
微分方程的解析解
求微分方程（组）解析解的命令:
dsolve(‘方程1’,‘方程2’,…,‘方程n’,‘初始条件’,‘自变量’) 记号: 在表达微分方程时，用字母 D 表示求微分，D2、D3 等 表示求高阶微分.任何 D 后所跟的字母为因变量，自变量可以指 定或由系统规则选定为缺省.
例如，微分方程
y' f ( x, y ) 对常微分方程 ： ，其 数 值 解 是 指 由 初 始 点 x 0 开 始 y ( x0 ) y0 的 若 干 离 散 的 x 处 的 值 ， 即 对 x 0 x1 x 2 x n， 出 准 确 值 y ( x1 ), 求 y ( x 2 ), , y ( x n ) 的 相 应 近 似 值 y1 , y 2 , , y n .
则微分方程变为一阶微分方程组：
y1 ' y 2 2 y 2 ' 1000 (1 y 1 ) y 2 y 1 y (0) 2, y (0) 0 1 2
2 1.5 1 0.5 0 -0.5 -1
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),'-')

x i 1 xi
f (t , y (t )) d t [ f ( x i , y ( x i ) ) f ( x i 1 , y ( x i 1 ) ) ]
x i 1 x i 2
故有公式：
h y i 1 y i [ f ( x i , y i ) f ( x i 1 , y i 1 )] 2 y0 y(x0 )
故有公式：
y i 1 y i h f ( x i , y i ) y0 y ( x0 ) i 0,1, 2, , n - 1
此即欧拉法．
2．使用数值积分 对方程y’=f(x,y), 两边由xi到xi+1积分，并利用梯形公式，有：
y ( x i 1 ) y ( x i )
-1.5 -2 -2.5 0 500 1000 1500
2000
2500
3000
3．结果如图
To MATLAB（ff4）
实 验 作 业
1． 一个小孩借助长度为a的硬棒拉(或推)某玩具．此小孩沿 某曲线行走，计算并画出玩具的轨迹． 2． 讨论资金积累、国民收入与人口增长的关系． （1）若国民平均收入x与人口平均资金积累y成正比，说明仅当总 资金积累的相对增长率k大于人口的相对增长率r时，国民平均收 入才是增长的． （2）作出k(x)和r(x)的示意图，分析人口激增会导致什么后果．
实际应用时，与欧拉公式结合使用：
0 y i( 1) y i hf ( x i , y i ) h ( k 1) (k ) y i 1 y i [ f ( x i , y i ) f ( x i 1 , y i 1 )] k 0 ,1, 2 , 2
返 回
（二）建立数值解法的一些途径
设 x i 1 x i h , i 0,1, 2, , n 1, 则 可 用 以 下 离 散 化 方 法 求 解 微分方程 y ' f ( x, y ) y ( x0 ) y0
1.用差商代替导数 若步长h较小，则有
y'(x) y(x h) y(x) h
由待解 方程写 成的M 文件名
ts=[t0,t f]，t0、 tf为自变
函数 的初 值
量的初值 和终值
ode23：组合的2/3阶龙格–库塔–费尔贝格算法 ode45：运用组合的4/5阶龙格–库塔–费尔贝格算法 用于设定误差限(缺省时设定相对误差10-3, 绝对误差10-6), 命令为：options=odeset（’reltol’,rt,’abstol’,at）, rt，at：分别为设定的相对误差和绝对误差．