minz??36060?(k)(k)??3??xi,1??xi,0?b--------(1)i?1i?1??60(k)??xi,0?11000-----------(2)?i?1k?1)k)?xi(,0?xi(,1......(i?0,1,..59)------(3)?(k)??xi,1?0............(i?60)--------(4) ?(k)(k)x??x?i,11i,0......(i?0)--------(5)?(k)(k)x??x,2,...,59)-----(6)2i,0......(i?1?i,1?x(k)??y(k)??y(k)---------(7)1122?0,060?k)xi(,0??(k)(k)?y1?y2?i?13---------(8)??2其中
式(1)表示目标函数,其值趋于0;
式(2)表示该动物在第k年初的数量,其值取11000;
式(3)和(4)表示该动物第i龄到了年底全部转化为(i+1)龄; 式(5)和(6)表示该动物各年龄段的变化; 式(7)表示该动物新生的幼仔数量;
式(8)表示该动物在13~60岁的雌性个体数量之和;
对b以10为间距,分别从50取到300头,利用MATLAB软件进行求解(源程序见附录8),得到以下表二:
被运走的数量(头) 被注射的数量(头) 50 60 70 80 90 100 110 120 130 140 370 361 351 342 332 322 313 304 294 285
被运走的数量(头) 被注射的数量(头) 150 160 170 180 190 200 210 220 230 240 275 266 256 247 237 228 218 209 199 190
被运走的数量(头) 被注射的数量(头) 250 260 270 280 290 300 180 171 161 152 142 133
问题四:
用人工手段可以使得该动物群落的数量减少,使其数量控制在一定范围内,但是,我们并不知道用了这种人工手段后,对于一次突发事件后要增加其数量时,该种群的复壮能力如何。因此,我们建立数学模型,分析研究其复壮能力。我们假设由于某种原因(例如由于一场大灾难导致该动物的死亡),使得该动物减少至M头,此时不得不停止注射避孕的方法。在初步模型中,我们先假设灾难发生后该种群虽然数量减少,但其年龄结构没有改变。建立该动物群落重新壮大的能力的模型如下: 模型七:
?d?(k)?t?60(k)60(k)??xi,1??xi,0-------(1)?i?1i?1?dt?60(1)??xi,0?M-------------(2)??i?1k?1)k)?xi(,0?xi(,1.....(i?0,1,..59)-------(3)?k)?xi(,1?0...........(i?60)---------(4) ?(k)(k)?xi,1??1xi,0......(i?0)---------(5)?(k)(k)x??x?2i,0.....(i?1,2,...,59)------(6)?i,1对M取不同的值,利用c语言进行编程求解,(源程序见附录9)可以得到以下表: M值 第t年 t=1 增量d1 665 594 526 456 388 322 250 179 113 10000 9000 8000 7000 6000 5000 4000 3000 2000 当前数量M1 t=2 增量d2 当前数量M2 t=3 10665 9594 707 638 8526 562 7456 489 7945 521 8466 557 6388 414 6082 444 7246 472 7718 507 8225 538 5322 342 5664 365 6029 390 6419 418 6837 445 7272 472 7744 507 8251 543 4250 267 4517 286 4803 305 5108 325 5433 348 5781 371 6152 401 6553 426 6979 456 7435 22.8061 3179 191 3370 205 3575 219 3794 235 4029 251 4280 268 4548 286 4548 307 4835 327 5142 16.4364 2113 122 2235 130 2365 140 2505 145 2650 156 2806 166 2972 178 3150 191 3341 204 3545 10.7636 11372 10232 9088 757 恢复 699 602 增量d3 当前数量M3 10932 9690 729 恢复 640 t=4 增量d4 当前数量M4 10330 9023 689 596 t=5 增量d5 当前数量M5 11019 9619 741 恢复 638 t=6 增量d6 当前数量M6 10257 8792 682 581 t=7 增量d7 当前数量M7 10929 9373 729 恢复 623 9996 665 t=8 增量d8 当前数量M8 t=9 增量d9 当前数量M9 10661 8794 707 581 t=10 增量d10 当前数量M10 11368 9375 46.6000 42.6857 38.8810 35.6424 28.6121 d?2?t?46.00增长加速度:D? 2dt00
该表对M进行不同的取值,分别表示某一年(t=0)灾难发生后该动物群落剩下的数量。我们用第一题的C程序即可算出该群落在下一年(t=1)的增量,把这个增量加上上一年群落的总数即为第一年(t=1)群落的总数M1,接着我们又用程序算出这一年的增量,加上M1即为第二年(t=2)的个体总数M2。以这样的算法推算下去,一直到群落数量恢复到11000左右或者前十年(t=10)为止。
将计算出来的各年增量按时间t(1,2??10)列在座标上,用MATLAB进行曲线拟合(具体程序见附录10),得到的图形如下所示:
从图中可以看出,每年的增量都是遵从一定的线性关系的,也就是说,该群落的增
d?2?t?长加速度是一定的,通过MATLAB的曲线拟合可以求出不同M值所对应的增2dt长加速度(程序见附录 ),数据如下: M值 10000 9000 46.6000 8000 42.6857 7000 38.8810 6000 35.6424 5000 28.6121 4000 22.8061 3000 16.4364 2000 10.7636 D(只/46.00年^2) 00
通过这些数据,我们可以知道,在灾难后群落数量减少至M时的重新壮大能力。
再用MATLAB求出M与D的关系,进行曲线拟合求出斜率k后画出图形(程序见附录11),如下所示
可知M与D是服从一定的线性关系的,由MATLAB算出的的数据可知
M=3.9010+0.0047D。
该模型的稳定性分析见本论文第四部分。
四. 模型的稳定性分析
1.对于问题一:我们利用C语言进行编程,求得了第k年时各年龄结构的分布情况,并利用MATLAB软件对得到的数据进行画图分析,从输出的图可看出各年该动物各个年龄的分布是围绕一定的轨迹进行分布的,从而说明了该模型是稳定的,其次,用统计软件SAS对该模型得出的数据进行相关系数的分析,得到相关系数R都接近于1,而P=0.0001<0.01则表示其相关性显著。故建立的模型一是很稳定的,其各年的该动物的总数基本上都在11000左右。
2.对于问题二、问题三:我们从两个不同方面对问题进行分析,分别从?=新生的幼仔数量-自然死亡的数量~0,以及?=该动物群落在年底时的总数量与年初的数量加上被运走的数量b的差值~0两个不同方面进行模型的建立,从而把该动物群落的稳定性问题转化为求目标函数最小化问题。(这是本论文的一个亮点)。利用MATLAB对模型进行求解,得出的数据最后都能使得该动物群落维持稳定,可见对于问题二、问题三所建立的模型是很稳定的,结果也是很优的。