生物种群动态模拟

西北农林科技大学实验报告
学院名称：理学院 专业年级：2011级信计1班 姓 名： 学 号：2011014816 课 程：数学模型与数学建模 报告日期：2013年12月14日
1 实验题目: 生物种群动态模拟
2 实验问题陈述: 在这个问题中我们采用较复杂的增长模型
g （p ）=(r/c)(p-c)(1-p/k),
其中参数c 表示种群水平的最小值，低于这个水平出现的负增长，假设参数a=10∧-7，最小值的种群水平是蓝鲸3000条，长须鲸15000条。

(a)两种鲸鱼种群能否共存？建一个动力系统模型描述种群的动态。

(b)画这个模型的向量场。

确定每个平衡态的位置。

(c)确定状态空间每个平衡态是否稳定。

d)假设存在5 000条蓝鲸和70 000条长须鲸。

关于这两种鲸鱼的将来这个模型能预测些什么？
3 实验目的：(1)建立恰当的动力系统模型来描述待求解问题；
(2)画出模型向量场并确定平衡态位置； (3)确定状态空间中的平衡态是否稳定； (4)预测这两种鲸鱼的将来情况。

4 实验内容
对于单个种群的增长，有
B dB r B dt =，F dF r F dt
=. 当竞争种群存在时，有模型
()()B F F B dB
r C F B dt
dF r C B F dt
⎧=-⎪⎪⎨
⎪=-⎪⎩ （4.2.1） 这就是两个种群的竞争模型。

在鲸鱼种群竞争的问题中，已知：0.05B r =，0.08F r =；4
1510B K =⨯，44010F K =⨯；由假设5，B F C C α==，但α未知。

由此可得模型为
()()0.050.08dB
F B dt
dF B F dt
αα⎧=-⎪⎪⎨
⎪=-⎪⎩ （4.2.2） 这是一个有两个状态变量的微分方程组。

主要依赖于定性分析的方法了解这两个种群
的动态特征。

这是一个二维的自治系统，相空间应该是二维的平面(,)N B F =。

针对实际问题的要求，应该在这个平面的第一象限{(,):0,0}S B F B F =≥≥中讨论系统的动态。

为叙述方便，将模型（4.2.1）记为：()1',B f B F B =，()2',F f B F F =，其中，
()1,0.05B F f B F r C F F α=-=-，()2,0.08F B f B F r C B B α=-=-。

5 实验结果分析与讨论
5.1 分析Ⅰ——种群的动态特征
（一）平衡点：()*
00,0x =，*
1,F B B F r r x C C ⎛⎫
=
⎪⎝⎭
（二）等倾线：
（三）定性分析：
*0
x 不稳定。

关于*1x 有： 'F F B F F
B F B B F B r r r
C r B r B C B F F C C C C ⎛⎫⎛⎫⎛⎫⎛⎫
=+-++≈- ⎪
⎪⎪ ⎪⎝⎭⎝⎭⎝⎭⎝⎭
'B
F B
F B F B F B F F r r r C r F r F C B F B C C C C ⎛⎫⎛⎫⎛⎫⎛⎫
=+-++≈- ⎪
⎪⎪ ⎪⎝⎭⎝⎭⎝⎭⎝⎭
0,0,0B F p q r r ==-<∆>，在*
1
x 附近的轨线为同心圆。

（四）周期性：系统的解为在相平面12(,)x x 的第Ⅰ象限内的封闭轨线。

① 轨线方程
()()
B F F B B r c F dB dF F r c B -=
- F B B F r r c dB c dF B F ⎛⎫⎛⎫
-=- ⎪ ⎪⎝⎭⎝⎭
有解 (
)()F
B B F r c B r c F
B
e F e S
--=
② 系统的解是相平面(),B F 的第Ⅰ象限内的封闭曲线。

③ 周期：给定0t ，存在10t t >，使得()()()()0101,B t B t F t F t ==，且对于任何
01t t t <<都不成立上述关系，则称10T t t =-为种群波动的周期。

5.2 分析Ⅱ——人类控制活动的影响
① 周期平均：()x t 为周期为T 的函数，则称
()00
1
t T
t x t dt T
+⎰
为()x t 的周期平均。

由方程可知'/F B F F r c B =-，于是有
()()
()
000
1
0t T
F t T F B t F t r c B dt dF F
++-==⎰
⎰ ()00*1t T F F t B B r r B t dt B T c c +-===-⎰ ()
00*
1t T B t F
r F t dt F T c +==⎰ ② 人类控制行为的影响
按比率12,εε对,B F 进行捕杀，则有模型
()()1122''B F B F F B F B B r B c BF B r B c BF F r F c BF F r F c BF
εεεε=--=--=--=-- （5.2.1）
模型的周期平均（即平衡点）为
*22*1F F F
B B B
B B
F F
r r r B c c c r r F c c εεε-+-=
=>--=< （5.2.2）
6 实验程序（Matlab或者其它软件语言陈述）
6.1 绘制向量场
clear;clc
syms B F
alpha = 0.00001;%1*10^(-7);
DB = 0.05*B-alpha*B*F;
DF = 0.08*F-alpha*B*F;
[Bsteady,Fsteady] = solve(DB,DF);
disp('平衡点为：')
disp([Bsteady Fsteady])
M=10;
Bmin=0; Bmax=5000;
Fmin=0; Fmax=70000;
[X1,X2]=meshgrid(Bmin:(Bmax-Bmin)/M:Bmax,Fmin:(Fmax-Fmin)/M,Fmax); dX1 = 0.05*X1-alpha*X1.*X2;
dX2 = 0.08*X2-alpha*X1.*X2;
quiver(X1,X2,dX1,dX2,'r');
axis([Bmin Bmax Fmin Fmax]);
hold on
ezplot(DB,[Bmin Bmax Fmin Fmax]),hold on
ezplot(DF,[Bmin Bmax Fmin Fmax])
xlabel('蓝鲸');ylabel('长须鲸');
运行程序可得：
平衡点为：
[ 0, 0]
[ 8000, 5000]
4
蓝鲸
长须鲸
(2 F)/25 - (B F)/100000 = 0
使用MATLAB 绘制向量场 function zqjh() x=input('x='); y=input('y=');
[x,y]=meshgrid(0:1:15,30:1:45); u=x.*(0.05-0.01*y);%方程B'=Bf1(B,F) v=y.*(0.08-0.01*x);%方程F'=Ff2(B,F) quiver(x,y,u,v); grid; hold; x=2:0.01:6;
plot(x,50-(10/3).*x); hold on %水平集LB x=0:0.01:15;
plot(x,40-0.5.*x);%水平集LF end
在MA TLAB 界面中，输入 x=3.5 y=38.2
Current plot held
其中，为了方便，对模型表达式做了一些处理：4
4
10,10x B y F =⨯=⨯。

运行程序可得：
-2
0246810121416
283032343638404244
46
注：上图即为鲸鱼竞争种群动态的向量场。

其中，横坐标为B ，纵坐标为F ，较陡的直线是水平集B L 的部分，较平缓的直线是水平集F L 的部分。

图中的箭线指示该点相轨线前进的方向，箭线的长短表明该点相轨线前进速度的大小。

7．实验总结
分析不同竞争力下的蓝鲸与长须鲸的种群共存动态状况。

本次实习模型涉及到了平衡态、水平集等知识，并得到在竞争力为-710时的鲸鱼竞争种群动态的向量场图，使得两种群的竞争状态表现的更加直观，这类模型的建立有利于对珍惜生物品种的生存状况及发展前景的分析，可以更好的保护它们。