y (t )
n > 0,乙方胜
n = 0,平局 n < 0,甲方胜
cy 2 − 2bx = n
n
=
cy
2 0
−
2bx 0
n>0 乙方胜
⎜⎜⎝⎛
y0 x0
⎟⎟⎠⎞ 2
>
2b cx 0
⎛⎜⎜⎝
y 0
x 0
⎞⎟⎟⎠ 2
>
2r p s x xx rs x y ry 0
设 x0=100, rx/ry=1/2, px=0.1, sx=1(km2), sry=1(m2)
数学建模讲义
主讲人:穆学文
西安电子科技大学数学系 Email:mxw1334@
2013/2/21
第四讲 微分方程模型
-------多种群增长模型
4.1 正规战与游击战 4.2 微分方程解析解 4.3 微分方程数值解 4.4 微分方程稳定性 4.5 捕鱼业的持续收获 4.6 军备竞赛 4.7 种群的相互竞争 4.8 种群的相互依存 4.9 种群的弱肉强食
( y0 / x0 )2 > 100
0
x(t) 乙方必须10倍于甲方的兵力
4.2 微分方程的解析解
求微分方程（组）的解析解命令: dsolve(‘方程1’, ‘方程2’,…‘方程n’, ‘初始条件’, ‘自变量’)
记号: 在表达微分方程时，用字母 D 表示求微分，D2、D3 等 表示求高阶微分.任何 D 后所跟的字母为因变量，自变量可以指 定或由系统规则选定为确省. 例如，微分方程 d 2 y = 0 应表达为：D2y=0.
用于设定误差限(缺省时设定相对误差10-3, 绝对误差10-6), 命令为：options=odeset（’reltol’,rt,’abstol’,at）, rt，at：分别为设定的相对误差和绝对误差.
注意:
1、在解n个未知函数的方程组时，x0和x均为n维向量， m-文件中的待解方程组应以x的分量形式写成.
⎪ ⎨
y
=
− bx
• 假设没有增援
⎪⎩ x ( 0 ) = x 0 , y ( 0 ) = y 0
正规战争模型
⎧ x = − ay
⎪ ⎨
y
=
−bx
⎪⎩ x(0) = x0 , y (0) = y0
为判断战争的结局，不求x(t), y(t) 而在相平面上讨论 x 与 y 的关系
dy = bx dx ay
（一）常微分方程数值解的定义
在生产和科研中所处理的微分方程往往很复杂且大多 得不出一般解。而在实际上对初值问题，一般是要求得 到解在若干个点上满足规定精确度的近似值，或者得到 一个满足精确度要求的便于计算的表达式。
因此，研究常微分方程的数值解法是十分必要的。
对常微分方程
：⎨⎧ ⎩
y' = y(x
f(x, 0) =
为自治系统或动力系统。
若方程或方程组f (x)=0有解xo,x = xo显然满足(3.28)。 称点xo为微分方程或微分方程组(3.28)的平衡点或奇点。
⎪⎧
y (0) i+1
=
yi
+ hf (xi , yi )
⎨ ⎪⎩
y (k +1) i+1
=
yi
+
h 2
[
f
(
xi
,
yi
)
+
f
(
xi+1
,
y (k ) i +1
)]
k = 0,1,2,"
对于已给的精确度 ε，当满足
y (k +1) i +1
−
y(k) i +1
< ε 时，取 yi+1
=
y ， （k +1） i +1
dx2 例 1 求 du = 1+ u 2 的通解.
dt 解 输入命令：dsolve('Du=1+u^2','t')
结果：u = tg(t-c)
例 2 求微分方程的特解.
⎪⎧d 2 ⎨ dx
y
2
+
4
dy dx
+
29
y
=
0
⎪⎩ y(0) = 0, y' (0) = 15
解 : 输入命令:
y=dsolve('D2y+4*Dy+29*y=0','y(0)=0,Dy(0)=15','x')
1、用差商代替导数 若步长h较小，则有
y'(x) ≈ y(x + h) − y(x) h
故有公式：
⎧ ⎨ ⎩
y y
i+1 = yi + 0 = y(x0 )
hf(xi,yi)i = 0,1,2,", n -1
此即欧拉法。
2、使用数值积分
对方程y’=f(x,y), 两边由xi到xi+1积分，并利用梯形公式，有：
∫ y(xi+1 ) − y(xi ) =
xi +1 xi
f
(t,
y(t ))dt
≈
xi+1 − 2
xi
[f
(xi , y(xi )) +
f
(xi+1 ,
y(xi+1 ))]
故有公式：
⎪⎧ ⎨
y
i
+1
=
yi
+
h 2
[
f
(
xi
,
yi
)
+
f (xi+1 , yi+1 )]
⎪⎩ y0 = y(x0 )
实际应用时，与欧拉公式结合使用：
• 甲方战斗减员率还随着甲方兵力的增加而增加 f(x, y)=−cxy, c~ 乙方每个士兵的杀伤率
c = ry py ry~射击率 py ~命中率
py=sry /sx sx ~ 甲方活动面积 sry ~ 乙方射击有效面积
g(x, y) = −dxy, d = rx px = rxsrx / sy
• 忽略非战斗减员 • 假设没有增援
y=simple(y)
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
2
2013/2/21
4.3 微分方程的数值解
y) ，其数值解是指由初始点x y0
0
开始
的若干离散的x值处，即对x0 < x1 < x2 < " < xn，求出准确值y(x1 ),
y(
x2
),",
y(
x
n
)
的相应近似值y 1
,
y
2
,",
y
。
n
（二）建立数值解法的一些途径
设 x i+1 − xi = h, i = 0,1,2,"n −1, 可用以下离散化方法求解微分方程： ⎧y'= f(x, y) ⎩⎨y(x0 ) = y0
2、使用Matlab软件求数值解时，高阶微分方程必须 等价地变换成一阶微分方程组.
3
2013/2/21
例4
⎪⎧ ⎨
d2x dt 2
−1000(1
−
x2
)
dx dt
−
x
=
0
⎪⎩ x(0) = 2; x'(0) = 0
解: 令 y1=x，y2=y1’
则微分方程变为一阶微分方程组：
⎧ ⎪ ⎨
y2
'
=
y1' = y2 1000(1− y12
0.6
0.4
0.2
0
-0.2
-0.4
-0.6
-0.8
-1
0
2
4
6
8
10
12
2、取t0=0，tf=12，输入命令： [T,Y]=ode45('rigid',[0 12],[0 1 1]);
plot(T,Y(:,1),'-',T,Y(:,2),'*',T,Y(:,3),'+')
3、结果如图
图中，y1的图形为实线，y2的图形为“*”线，y3的图形为“+”线.
• 每方战斗减员率取决于双方的兵力和战斗力 • 每方非战斗减员率与本方兵力成正比 • 甲乙双方的增援率为u(t), v(t)
⎧ x(t) = f ( x, y) − αx + u(t), α > 0
模型
⎨ ⎩
y (t)
=
g
(x,
y)
−
βy
+
v (t ),
β >0
f, g 取决于战争类型
正规战争模型 双方均以正规部队作战
4.1 正规战与游击战
第一次世界大战Lanchester提出预测战 役结局的模型.
• 战争分类：正规战争，游击战争，混合战争 • 只考虑双方兵力多少和战斗力强弱 • 兵力因战斗及非战斗减员而减少，因增援而增加 • 战斗力与射击次数及命中率有关
一般模型 x(t) ~甲方兵力，y(t) ~乙方兵力
模型 假设
4.4 稳定性问题
在研究许多实际问题时，人们最为关心的 也许并非系统与时间有关的变化状态，而是系 统最终的发展趋势。例如，在研究某频危种群 时，虽然我们也想了解它当前或今后的数量， 但我们更为关心的却是它最终是否会绝灭，用 什么办法可以拯救这一种群，使之免于绝种等 等问题。要解决这类问题，需要用到微分方程 或微分方程组的稳定性理论。在下两节，我们 将研究几个与稳定性有关的问题。