一类非线性时滞微分方程模型系统的稳定性和Hopf分支分析

1、第33卷第3期2009年6月 江西师范大学学报(自然科学版)JOURNALOFJIANGXINORMALUNIVERSITY(NATURALSCIENCE)Vo.l33No.3 Jun.2009文章编号:1000 5862(2009)03 0335 06一类非线性时滞微分方程模型系统的稳定性和Hopf分支分析刘华祥, 曾广洪, 吴庆初122(1.广东海洋大学理学院数学系,广东湛江 524088;2.江西师范大学数学与信息科学学院,江西南昌 330022)摘要:探讨了2个具有竞争相克作用的浮游植物种群的二维非线性时滞微分方程模型.分析了当离散时滞在临界值附近变动时,无时滞和有时滞系统的共生平衡点

2、由局部渐近稳定变为不稳定.应用Hopf分支理论得出了分支周期解的存在条件,揭示了时滞 分支参数对系统基本动力学性质的影响.对所得的结果进行了数值模拟和可视化验证.通过选择适当的时滞参数值,所建立的时滞微分方程模型能够捕捉到各种浮游植物物种通常所展示的周期性暴发的振荡特性,这对于有害浮游植物的预防和控制具有一定的指导意义.关键词:浮游植物;时滞;Hopf分支;数值模拟中图分类号:O175.14 文献标识码:A0 引言浮游植物是水生食物链的基本营养层,它们不仅提供着水生环境中基本的食物资源,而且通过光合作用释放出大量的O2,又从周围的环境中吸收着大量的CO2,因此浮游植物对大自然的滋养和调节作用不

3、言而喻.然而,最近几十年来,有害水生植物的暴发在全球范围内有漫延趋势,它们不仅严重危害着有益水生动植物的生长和生存,而且给环境带来巨大威胁.例如,有些有害浮游植物能产生大量毒素,造成水产食物的污染,毒死鱼类,有些则通过自身的大量孳生漫延而引起水生环境的缺氧症,造成有益水生生物的大规模死亡,严重破坏水生环境中生物的多样性.因此有害水生植物的暴发及其预防和控制问题已经引起了人们的广泛关注,不少研究人员从生物学和生态学的角度对此进行了大量的研究,得出了很多有价值的研究成果.例如,研究人员证实,水生植物暴发及其自然演替背后的关键因素是,某些水生植物物种释放出的大量有害物质,对其他水生植物产生相克或抑制

4、作用,从而影响有益水生植物的生存,甚至造成有益水生植物的消亡.人们一直想弄清楚有害水生植物周期性暴发以及其暴发呈递增趋势的背后原因和机制,从而为有效的治理和防治措施提供科学的决策依据.从种群动力学角度研究上述问题早就成为了数学生态学感兴趣的问题.根据生物学上植物间的相克作用的机理,即植物通过产生种间相克毒素或刺激剂改变着某个种群的密度,从而实现植物相克的作用,Mol ish提出了最初的植物相克概念1.后来Rice将这一概念表述为 植物相克就是一种植物物种通过向周围2环境中释放化合物质而对其他物种的生长产生影响 .现在这一概念已被数学生态学广泛地用于定义浮游植物群落间的相克作用.Maynard

5、Smithc以著名的Lotka Volterra竞争模型为基础,提出了第1个表示2个具有竞争相克作用的两物种的数学模型收稿日期:2008 12 22基金项目:广东海洋大学科研(0612163)和江西师范大学青年成长基金(20071978)资助项目.作者简介:刘华祥(1958 ),男,江西泰和人,教授,主要从事生物数学方面的研究.3:336江西师范大学学报(自然科学版)2009年dN1(t)=N1(t) 1- 1N1(t)-!1N2(t)-1N1(t)N2(t),dtdN2(t)=N2(t) 2- 2N2(t)-!2N1(t)-2N1(t)N2(t),dt其中N1(t),N2(t)分别表示2种竞

6、争浮游植物种群的密度; 1, 2表示2个物种各自的增殖速度; 1, 2表示2个物种各自的种内竞争系数;1,2分别是2个物种彼此释放的 种间化学物质 对对方生长的抑制速度.后来,Bandyopahyay考虑了由无毒和有毒竞争浮游植物物种的数学模型,以N1(t)和N2(t)分别表示无毒和有毒浮游植物物种的密度,并设有毒物种对无毒物种的相克作用项采取形式1(N1(t)N2(t),则这2个浮游植物种群的相互作用模型由下面的非线性常微分方程dN1(t)=N1(t) 1- 1N1(t)-!1N2(t)-1N1(t)N2(t),dtdN2(t)=N2(t) 2- 2N2(t)-!2N1(t)dt上述提到的浮

7、游植物种群系统,Chattopadhyay5 924(1)决定,其中1表示由有毒物种释放的毒素对无毒物种的抑制速度,而其他参数的生物学意义与前述相同.对等人进行了详尽的分析,得出了这些系统的较丰富的动力学性质,如共生平衡点的局部稳定性条件等.遗憾的是所有这些模型研究都未能捕捉到引起浮游植物种群周期性暴发的振荡特性.本文将在模型(1)的基础上,考虑在2个物种的种内竞争项中引入时滞参数 ,得到下述修改后的模型:dN1(t)=N1(t) 1- 1N1(t- )-!1N2(t)-1N1(t)N2(t),dtdN2(t)=N2(t) 2- 2N2(t- )-!2N1(t),dt(2)并设初始条件为:N1

8、(#)!0,N2(#)!0,- #0,其中考虑了受损资源的恢复需要有一定的时间,因此对2种浮游植物的种内竞争采用Hutchinson型时滞竞争项 iNi(t)Ni(t- ),以取代经典的种内竞争项 i(Ni(t),从而得到一类非线性时滞微分方程模型系统.常微分方程模型系统(1)被推广到时滞微分方程模型系统(2)具有重要的意义,因为在确定性环境中,时滞具有打破平衡的扰动作用,时滞的存在经常起到引起种群密度振荡性变化的支配作用.模型(1)未能捕捉到浮游植物密度的振荡性,因此自然期待模型(2)能抓住有毒和无毒浮游植物物种竞争模型中的振荡特性.21 正平衡点的稳定性及Hopf分支根据系统的生态学意义,

9、只对正平衡点感兴趣.在对模型系统(2)分析之前,先简要地对相应的无时滞模型系统(1)的正平衡点作必要的分析.令系统(1)的右端为零,容易计算出系统(1)的正平衡点E(X,Y),其中Y是下述二次方程1 2y+( 1 2-!1!2-1 2)y+ 1!2- 1 2=0的正解.而X=( 2- 2Y)/!2时,易见,当1!2- 1 2<0时,方程(3)有唯一正解,此时若还有Y< 2/ 2,则系统(1)有唯一正平衡点.记(H1)*2*(3)(4)X>0,Y>0, .1!2- 1 2<0-!1XY- 1X 1X-1X -!2Y !2X- 2系统(1)在正平衡点E(X,Y)处的线

10、性化矩阵为M=其特征方程为(21+122 2X1Y=22,第3期刘华祥,等:一类非线性时滞微分方程模型系统的稳定性和Hopf分支分析337因为r= ,2+ 1Y+1XY-!2X= 1X+1XY+ 2Y>0s=-!2X 1+ 1 ,2X+1 2XY-1!2XY=( 1 2- 1!2)X+ 21XY>0所以由Hurwitz判据,系统(1)的正平衡点是渐近稳定的,于是得到下述定理.定理1 设条件(H1)成立,则系统(1)存在唯一正平衡点E,且正平衡点是局部渐近稳定的.下面探讨时滞模型系统(2).在同样的参数限制条件下,系统(1)的正平衡点E也是系统(2)的正平衡点.将以时滞 作为控制参数

11、,来探讨它的变化对该时滞系统的正平衡点的稳定性的影响.已经知道当 =0时,E是局部渐近稳定的,为了对模型系统(2)的正平衡点进行局部稳定性分析,先借助于平移变换:N1(t)=x(t)+X,N2(t)=y(t)+Y,把系统(2)化为=ax(t)+by(t)- Xx(t- )+h.o.t,1dt=cx(t)- Yy(t- )+h.o.t,2dt2*22(6)其中a=-!,b=-!1X-<0,c=-!,h.o.t表示2次及1Y-21XY+ 1- 1X=-1XY<01X2Y<0以上的高阶项,则系统(6)的平凡平衡点对应着系统(2)的正平衡点.系统(6)的线性化系统为=ax(t)+by

12、(t)- 1Xx(t- ),dt=cx(t)- 2Yy(t- ).dt(7)的特征方程为g()=-a-bc+( 2Y+ 1X-a 2Y)e+ 1X 2Ye=0.(8)对系统(6)或(2),根据稳定性判据易知,如果当 =0时,(8)的特征根的实部均为负,并且当 !0时,对于所有的实数%,有g(%i)#0,则系统的相应平衡点在 !0时都是局部渐近稳定的.因此,为了讨论相反的情况,下面将假定对于某个 >0,存在实数%#0,使得g(%i)=0,即假定对某个正数 ,%i是上述特征方程的纯虚根.分离出g(%i)=0的实部与虚部可得%+bc= n(% )-a 1X%sin(% )+ 2Y%si2Yco

13、s(% )+ 1X 2Ycos(2% ),a%= 1X%cos(% )+ 2Y%cos(% )+a 2Ysin(% )- 1X 2Ysin(2% ).将(9)的两式分别平方后相加可得(%+bc)+(a%)=-2 1X 2Y%sin(% )-2 1X 2Y%sin(% )-2a 1X 2Ycos(% )+1X 2Y+ 1X%+2 1X 2Y%+a 2Y+ 2Y%.设%>0,则将上式右端适当放大后可得%+(- 1X- 2Y-2a 1X 2Y+2bc+a)%+(-2 0.记P=- 1X- 2Y-1X 2Y-2 1X 2Y)%- 1X 2Y-a 2Y-2|a| 1X 2Y+bc2a 1X 2Y+

14、2bc+a,Q=2 1X 2Y+2 1X 2Y,R=- 1X 2Y-a 2Y-2|a| 1X 2Y+bc,则上述不等式可改写为%+P%-Q%+R0.若(H2)成立,则不等式(10)无实数解.因此当(H3)P<0或R<Q/(4P)2422222222222222222222222222222222222242222222222222222222222222222222- -2 (7)(9)(10)P>0,R>Q/(4P)2成立时,不等式(10)有实数解,由此可推知存在实数%,使g(%i)=0.又由(9)式可得2a%sin(% )+ 1X 2Y-(%+bc)cos(% )=

15、a 2Y, n(% )+a%cos(% )=( ,1X 2Y+(%+bc)si1X+ 2Y)%(2(11)3382江西师范大学学报(自然科学版)2009年2Y%a+(%+bc)- 1X 2Y( 1X+ 2Y)%sin(% )=,a%+(%+bc)-( 1X 2Y)cos(% )=( 1X+ 2Y)%a- 2Ya(%+bc)+ 1X 2Y22222a%+(%+bc)-( 1X 2Y)2222222(12)于是%是由(12)式及(sin(% )+(cos(% )=1组成的方程组确定的正根.由此可定义Y%a+(%+bc)- 1X 2Y( 1X+ 2Y)&=arcsi2,22222a%+(%+

16、bc)-( 1X 2Y)则时滞参数的临界值可表为 = ,k=0,1,2,3,.设 是 k=(&+2k)/%(k=0,k=(&+2k)/%1,2,)中的任一正数,则已证明存在正数 及%,使得特征方程(8)有一对纯虚根%i%.下面再研究当时滞参数在临界值 k(>0)的小邻域内变化时特征方程(8)的根实部符号的变化情况.设=(+i%是(8)的根,将其代入(8)中并分离出实部和虚部,得 K1(,%, )=(-%-bc-a(-a 2Ye1X(e-(22-( *cos(% )+ 1X 2Ye-(-2(-2(cos(2% )+-(cos(% )+%e-(-(sin(% )+ 2Y(ec

17、os(% )+%e-(sin(% ),-(K2(,%, )=2(%-a%+a 2Ye1X(-(eK1(,%, )K1(,%, )%于是(K1(,%, ),K2(,%, )( ),%=%( ), &K1(,%, ) K2(,%, )-(sin(% )- 1X 2Ye-(sin(2% )+sin(% )+%ecos(% ),sin(% )+%ecos(% )+ 2Y(-(e已有K1(0,%, ,K2(0,%, ,计算可得k)=0k)=0=(0,%, k)K2(,%, )=-a+(a k+1) 2Y+ 1Xcos(% k)-(0,%, k)% n(% k( 1X+ 2Y)sik)-2 k 1

18、X 2Ycos(2% k),=-(0,%, k)K2(,%, )(=(0,%, k)-2%+(a n(% k+1) 2Y+ 1Xsik)+% k( 1X+ 2Y)cos(% k)-2 k 1X 2Ysin(2% k),K1(,%, )2=(0,%, k)+K1(,%, )>0,(0,%, k)由隐函数存在定理知,方程组K1(,%, )=0,K2(,%, )=0在(0,%, k)的邻域内定义了隐函数组(=(-)+ k, )()>0),使得( ,%( )=%,又计算可得k+)k)=0=( k)+a% 2Ysin(% k)-2 1X 2Y%sin(2% k)=1X+ 2Y)%cos(%(

19、0,%, k)2%a%- 1X 2Ysin(2% k),=-( 1X+ 2Y)%sin(% k)+a% 2Ycos(% k)-2 1X 2Y%cos(2% k)=(0,%, )k(13)2-%(%+bc)+ 1X 2Ycos(2% k),利用(9),(12)(14),再进行计算可得 - (K1,K2)=-(0,%, k)2(14)%p -2%+( 1X+ 2Y)sin(% k)+ kp-%q -a+( 1X+ 2Y)cos(% k)- kq%p -2%+( 1X+ 2Y)sin(% k)-%q -a+( 1X+ 2Y)cos(% k)222222222=-=%+a 2Y+( 1X 2Y)-bc

20、-a% 1X 2Ysin(2% k)+2% 1X 2Ycos(2% a 1X 2Ycos(% k)-2k)=%+ 2Y(a-2a 1Xcos(% k)+( 1X)-bc+Xk)2k)4222222第3期刘华祥,等:一类非线性时滞微分方程模型系统的稳定性和Hopf分支分析2339其中p=a%- 1X 2Ysin(2% k),q=(%+bc)+ 1X 2Ycos(2% k).所以当(H4)%+ a 1Xcos(% 2Y(a-2k)+( 1X)-bc+% 1X 2Y(2%cos(2% k)-asin(2% k)>0成立时,有=- (K1,K2)/ ( ,%)/ (K1,K2)/ (,%)d k

21、Hopf分支.综上所述,可以叙述关于系统(2)的稳定性及Hopf分支定理.定理2 对于系统(2),()若条件(H1)(H2)成立,则系统的正平衡点E对一切 !0都是局部渐近稳定的;()若条件(H1)、(H3)及(H4)成立,则对于时滞参数的每个临界值 k(>0),存在相应的正数) < 是渐近稳定的,当 k,使得当 k-)kk时,系统的正平衡点Ek< k+)k时,系统的正平衡点E是不稳定的;系统(2)在每个临界值 = k处发生Hopf分支.*(0,%, k)4222222>0,*这样便得到截断性条件,从而Hopf分支定理的条件都得到满足,因而在系统的正平衡点E的附近将发生

22、2 数值模拟与讨论上节定理2的意义在于:它提供了时滞参数临界值存在的条件以及利用时滞参数的临界值判别这类时滞系统存在周期解的条件,由于这些条件都容易判断并可利用计算机进行计算,所以是一个有实际价值的结果10.下面的数值计算就是应用定理2研究这类时滞模型系统的一个算例.dN1(t)=N1(t)2-0.6N1(t- )-0.05N2(t)-0.009N1(t)N2(t),dtdN2(t)=N2(t)1-0.7N2(t- )-0.015N1(t),dt考虑形如(2)的如下时滞模型系统其中 , 1=0.6,!1=0.05,.009, 2=1, .7,!.015.经过计算不难得到该系统的正1=21=02

23、=02=0平衡点及各有关量的值:E(3.155503145,1.360953504),a=-0.03865043755,b=-0.2473899581,c=-0.02041430256,P=-8.087947039,Q=10.26647629,R=-3.387444560.由于 .428571429>Y=1.360953504,P<0,R<Q/4P,所以1/ 2=2<40= 1/2, 2/ 2=1(H1),(H3)成立.进一步计算可得A=10.26647629%,B=-0.1328271745,C=8.098047626%+3.254642891-(%+0.005050

24、29345),%=2222*2.325416971, 0=0.358,*且有d(/d 0=*0.5507046174>0,这样定理2的条件都得到满足,所以当 在 0附近变化,且 < 0时,正平衡点E是稳定的(见图1),而当 > 0时,系统是不稳定的,因此E在附近分支出周期解(见图2),其中取(N1(0),N2(0)=(2,3).图1时滞为0.45时x(实线)及y(虚线)随时间的变化 图2时滞为0.89时x(实线)及y(虚线)随时间的变化340江西师范大学学报(自然科学版)2009年3 结论在本文中,分析了一类反映有毒和无毒浮游植物种群竞争相克作用的非线性时滞微分方程模型系统,

25、简要讨论了相应的无时滞系统的正平衡点的局部稳定性,着重对时滞模型系统的正平衡点的局部稳定性和Hopf分支进行分析,建立了该系统在正平衡点的局部Hopf分支存在的条件,证明了当时滞参数在临界值附近作微小变化时,正平衡点从局部稳定变化为不稳定,系统在正平衡点附近发生Hopf分支.研究表明,时滞参数的引入将对系统的稳定性产生扰动,并在一定的临界值附近导致2种群规模的振荡性变化,因此这类关于两竞争浮游植物种群的时滞微分方程模型能捕捉到自然环境中浮游植物通常所展示的周期性振荡特性,因此本文的研究结果有一定的实际意义.参考文献:1OkaichiT,AndersonDM,NemotoT(Eds),eta.l

26、Biology,environmentalscienceandtoxicologyM.NewYork:ElsevierPress,1989.2RiceE.AllelopathyM.NewYork:AcademicPress,1984.3Maynard Smith.ModelsinEcologyM.Cambridge:CambridgeUniversityPress,1974.4BandyopadhyayM.DynamicalanalysisofaallelopathicphytoplanktonmodelJ.JBiolSys,2006,14:205 218.5ChattopadhyayJ.Ef

27、fectsoftoxicsubstancesonatwo speciescompetitivesystemJ.EcolMode,l1996,84:287 289.6AbdllaouiAE,ChattopadhyayJ,ArinoO.Comparisons,bymodelsofsomebasicmechanismsactingonthedynamicsofthezooplanktontoxicphytoplanktonsystemsJ.M3AS,2002,12(10):1421 1451.7Anderson.Marinebiotoxinsandharmfulalgae:anationalplan

28、R.WoodsHole,Mass:WoodsHoleOceanographicInstitutionTechnicalReport,1993.8Chattopadhyay,SarkarRR,MondalS.Toxin producingphytoplanktonmayactasabiologicalcontrolforplanktonicblooms fieldstudyandmathematicalmodelingJ.JTheorBio,l2002,215:333 344.9SarkarRR,ChattopadhyayJ.Theroleofenvironmentalstochasticityinatoxicphytoplankton non toxicphytoplankton zooplanktonsystemJ.Environmetrics,2003,14:775 792.10曾广洪,刘华祥,吴庆初.应用计算机实现规范型的自动化简J.南昌大学学报:理科版,2008,32(2):117 118.11马知恩.种群生态学的数学建模与研究M.合肥:安徽教育出版社,2000.11StabilityandHopfBifurcationAnalysisforaNonlinearDifferentialEquationModelSystemwi