劉白茹,劉俊利,呂 潘

(西安工程大學(xué) 理學(xué)院， 陜西 西安 710048)

0 引 言

在自然界中，捕食現(xiàn)象是最為普遍的一種，也是當(dāng)前研究的熱點問題。隨著捕食者-食餌模型研究的深入，不同類型的捕食者-食餌模型被提出[1-4]。在研究捕食者-食餌模型時，考慮物種的階段結(jié)構(gòu)十分必要[5-7]。種群在不同的階段會表現(xiàn)出不同的生理特征。一般來說，研究者們會將食餌或捕食者分成幼年和成年2個階段進行考慮。種群從幼年到成年的過程，會經(jīng)歷一段時間間隔，稱之為種群從幼年到成年存在成熟時滯。2016年，宋燕等研究了一類具有階段結(jié)構(gòu)及B-D功能反應(yīng)的捕食者-食餌模型，討論了在一定條件下，當(dāng)τ通過臨界值τ0時，正平衡點將改變其穩(wěn)定性，產(chǎn)生Hopf分支[8]。

除了種群的階段結(jié)構(gòu)和時滯對模型的影響，當(dāng)捕食者出現(xiàn)在食餌面前時，食餌會對捕食者產(chǎn)生恐懼，食餌的出生率會受到影響，從而影響食餌種群規(guī)模的大小。此時，捕食者對食餌的影響程度甚至超過直接捕殺的效果[9-11]。2011年，ZANETTE等通過歌雀實驗，證明了僅僅由于對捕食者產(chǎn)生的恐懼，導(dǎo)致成年歌雀的繁殖數(shù)量減少了40%[12]；2016年，WANG等首次提出了具有恐懼效應(yīng)的捕食者-食餌模型，理論分析及數(shù)值模擬結(jié)果表明，恐懼效應(yīng)會使系統(tǒng)更加穩(wěn)定[13]。2017年，WANG等還提出了一個具有恐懼代價和適應(yīng)躲避捕食者代價的捕食模型[14]。數(shù)學(xué)分析表明，恐懼效應(yīng)與時滯的相互作用，共同決定了種群的長期動力學(xué)行為：當(dāng)時滯處于中間值時，正平衡態(tài)可能失去穩(wěn)定性；而當(dāng)時滯較大時，正平衡可能恢復(fù)其穩(wěn)定性。2018年，王蓉等考慮了一類帶有恐懼效應(yīng)的捕食者-食餌模型的共存態(tài)問題[15]。利用譜分析和分支理論的方法，以捕食者的死亡率m為分支參數(shù)，討論了系統(tǒng)在半平凡解附近出現(xiàn)分支現(xiàn)象。

在上述文獻的基礎(chǔ)上，本文考慮食餌具有階段結(jié)構(gòu)，將食餌從幼年到成年的成熟時滯與恐懼效應(yīng)相結(jié)合，并考慮捕食者對食餌的功能反應(yīng)，以幼年食餌的出生率r和恐懼因子k作為分支參數(shù)，研究一類時滯微分方程模型。

1 模型的建立

(1)

式中:x1(t)、x2(t)、y(t)分別表示幼年食餌、成年食餌和捕食者在t時刻的密度；r為幼年食餌的出生率；τ為幼年食餌成長為成年食餌的成長期；d1、d2、d3分別為幼年食餌、成年食餌和捕食者的死亡率；a為成年食餌內(nèi)部競爭導(dǎo)致的死亡率；n為捕食者捕食食餌后轉(zhuǎn)化為自身生物量的轉(zhuǎn)化率。r、d1、d2、d3、a、m、n、b為正常數(shù)，k、τ為非負常數(shù)。

一般來說，成年食餌會將幼年食餌保護起來并喂養(yǎng)它們，因此幼年食餌受到捕食者的攻擊可以忽略不計，故在模型(1)中假設(shè)捕食者種群僅以成年食餌為食。模型(1)的初始條件為：

x1(θ)=φ1(θ)≥0，x2(θ)=φ2(θ)≥0，

y(θ)=φ3(θ)≥0，θ∈[-τ,0]

(2)

設(shè)φ=(φ1(θ),φ2(θ),φ3(θ))∈C+=C([-τ,0]；R3，R=[0,+∞))，C+表示從[-τ,0]到R3,R=[0,+∞)的連續(xù)函數(shù)所組成的Banach空間。對于連續(xù)函數(shù)u：[-τ:σφ)→R3，σφ>0,當(dāng)t≥0時，定義

ut(θ)=u(t+θ)，?θ∈[-τ,0]

定義集合

下面給出模型解的存在性、非負性和有界性。

引理1[16]考慮時滯微分方程

I′(t)=cI(t-τ)-dI(t)-hI2(t)

其中c，d，h，τ>0，當(dāng)t∈[-τ,0]時,有I(t)>0，則有：

定理1?φ∈X，模型(1)有唯一的非負解u(t,φ)，u0=φ，t∈[0,+∞)。

證明給定φ∈X，定義G(φ):=(G1(φ),G2(φ),G3(φ))，其中

d1φ1(0)

?φ∈X，G(φ)在R×X中每個緊集上都是連續(xù)的，關(guān)于φ是李普希茲連續(xù)的。由文獻[17]中定理2.3可知，?φ∈X，模型(1)通過點(0,φ)的解在其最大存在區(qū)間[0,σφ)是唯一的。

對模型(1)中第1個式子求解得到

當(dāng)φ∈X，φi(0)=0時，有Gi(φ)≥0，i=2,3;根據(jù)文獻[18]中定理5.2.1可知，?t∈[0,σφ)，ui(t,φ)≥0，i=2,3。因此，?t∈[0,σφ)，x1(t)≥0，則模型(1)的解是非負的。

下面證明模型(1)滿足初始條件(2)的正解是最終有界的。

由模型(1)第2個方程可知

令

L(t)=nx1(t)+nx2(t)+y(t)

當(dāng)t>T1時，得

其中u*=min{d1,d2,d3},則當(dāng)t>T2時，由比較定理知

因此x1(t)、x2(t)、y(t)有界 ，則σφ=∞。

2 平衡點的存在性和穩(wěn)定性

令模型(1)右邊的表達式等于零，計算得：

1) 模型(1)總存在滅絕平衡點E0(0,0,0)。

2) 當(dāng)rexp(-d1τ)-d2>0時，模型(1)存在邊界平衡點E1(x11,x21,0)，其中：

令E(x1,x2,y)是模型(1)任意一平衡點，則點E對應(yīng)的特征方程為

(3)

定理2當(dāng)rexp(-d1τ)d2時，E0不穩(wěn)定。

(λ+d1)(λ-rexp(-(d1+λ)τ)+

d2)(λ+d3)=0

(4)

方程(4)有2個特征值λ1=-d1，λ2=-d3，其他特征根滿足方程

λ-rexp(-(d1+λ)τ)+d2=0

(5)

當(dāng)τ=0時，式(5)式變?yōu)棣?r+d2=0，λ3=r-d2。當(dāng)rd2時，λ3>0，特征方程有一正根，E0不穩(wěn)定。

當(dāng)τ>0時，令λ=iω，ω∈R，代入式(5)中得到

iω-rexp(-(d1+iω)τ)+d2=0

分離實部與虛部得到

將此式分別平方后相加，得

ω2=(rexp(-d1τ)-d2)(rexp(-d1τ)+d2)

(6)

當(dāng)rexp(-d1τ)

定理3當(dāng)rexp(-d1τ)>d2時，邊界平衡點E1存在。

ⅰ) 當(dāng)nm-d3≤0時，平衡點E1是局部漸近穩(wěn)定的;

證明將RE1代入式(3)中，得到邊界平衡點E1所對應(yīng)的特征方程

(λ+d1)(λ-rexp(-(d1+λ)τ)+d2+2ax21)·

(7)

1) 若nm-d3≤0，此時顯然λ2<0。

2) 若nm-d3>0，則當(dāng)

時，λ2<0; 當(dāng)

時，λ2>0。

特征方程(7)的其他根滿足方程

λ-rexp(-(d1+λ)τ)-d2+

2rexp(-d1τ)=0

(8)

當(dāng)τ=0，式(8)變?yōu)棣?r-d2=0，則λ3=d2-r<0，特征方程(7)的所有根均為負數(shù)。當(dāng)τ>0時，令λ=iω，ω∈R，代入式(8)中得到

iω-rexp(-(d1+iω)τ)-

d2+2rexp(-d1τ)=0

(9)

將式(9)實部與虛部分離并整理得到

ω2=(-rexp(-d1τ)+d2)·

(3rexp(-d1τ)-d2)<0

故方程(9)無實數(shù)解，即方程(8)的根均具有負實部。

綜上所述，在定理3的ⅰ)、ⅱ)條件下，?τ≥0，點E1所對應(yīng)的特征方程(7)只有負實部的根；在定理3的ⅲ)條件下，因為λ2>0，平衡點E1不穩(wěn)定，因此定理3成立。

則當(dāng)k4-k3>0時，正平衡點E*是局部漸近穩(wěn)定的；當(dāng)k4-k3<0時，當(dāng)τ從0增加時，E*可能失去穩(wěn)定性。

證明將RE*代入式(3)，得到正平衡點E*所對應(yīng)的特征方程

(10)

方程(10)有一個特征值λ1=-d1<0，其余特征值滿足方程：

λ2+k1λ+k2exp(-λτ)λ+

k3exp(-λτ)+k4=0

(11)

式中：

當(dāng)τ=0時，式(11)變?yōu)?/p>

λ2+(k1+k2)λ+k3+k4=0

(12)

當(dāng)τ>0時，令λ=iω，ω∈R，代入式(11)并分離實部與虛部，得到

(13)

整理得

(14)

式中：

z=ω2

則

1) 當(dāng)k4-k3>0時，方程(11)沒有純虛根，?τ>0，正平衡點E*局部漸近穩(wěn)定。

2) 當(dāng)k4-k3<0時，方程(11)有一組純虛根±iω0，當(dāng)τ從0增加時，E*可能失去穩(wěn)定性。其中

是方程(14)唯一的正根。因此定理4成立。

定理5當(dāng)rexp(-d1τ)

證明設(shè)(x1(t),x2(t),y(t))是模型(1)滿足初始條件式(2)的正解。由模型(1)第2個方程得

綜上所述，E0是全局吸引的，結(jié)合E0局部漸近穩(wěn)定性，得平衡點E0是全局漸近穩(wěn)定的。

定理6當(dāng)rexp(-d1τ)>d2且定理3中條件ⅰ)或條件ⅱ)成立時，平衡點E1是全局漸近穩(wěn)定的。

證明由模型(1)第2個方程得

2) 當(dāng)條件ⅱ)成立時，有nm-d3>0，

當(dāng)t>t3+τ時，由模型(1)中第2個方程可知

由引理1和比較定理知

由ε1的任意性可知

于是對于上述ε1>0，?t4>t3+τ，當(dāng)t>t4時，x2(t)>x21+ε1。因此

由模型(1)中第1個方程，當(dāng)t>t4+τ時，有

rexp(-d1τ)(x21+ε1)-d1x1(t)

由比較定理和ε1的任意性可知

因此

綜上所述，E1是全局吸引的，結(jié)合E1的局部漸近穩(wěn)定性知，邊界平衡點E1是全局漸近穩(wěn)定的。

3 系統(tǒng)(1)的持久性

定義：

X0={φ=(φ1,φ2,φ3)∈X|φ3(0)>0且?θ0∈[-τ,0],使得φ2(θ0)>0}

X1={φ∈X|φ3(0)=0}

X2={φ∈X|φ3(0)>0,φ2(θ)≡0,?θ∈[-τ,0]}

則X=X0∪X1∪X2。

定理7假設(shè)nm-d3>0，

(nm-d3)(rexp(-d1τ)-d2)>ad3b

證明首先證明X0，X1，X2關(guān)于系統(tǒng)(1)是正不變的。因為

因此,若y(0)=φ3(0)>0，則y(t)>0，?t≥0；若y(0)=φ3(0)=0，則y(t)≡0，?t≥0。所以X1是正不變的。

由模型(1)第2個式子得

(15)

y(t)有界，因此若x2(t)是持久的，則x1(t)也是持久的。

令?X0=X1∪X2，模型(1)中具有初值(φ1,φ2,φ3)的解的ω極限集為ω(φ1,φ2,φ3)。記Ω=∪{ω(φ1,φ2,φ3)|(φ1,φ2,φ3)∈?X0},則在?X0上：

1a) 若rexp(-d1τ)

2) 若φ∈X2，則x2(t)≡0，?t≥0。因此當(dāng)t→+∞時，有x1(t)→0，y(t)→0，即

(x1(t),x2(t),y(t))→RE0=(0,0,0)

因此Ω={RE0,RE1}，若RE0、RE1是孤立的，則Ω也是一個孤立的非循環(huán)覆蓋。RE0與RE1的孤立不變性可以從式(16)的證明中得到，即

Ws(RE0)∩X0=?，Ws(RE1)∩X0=?

(16)

(17)

因為rexp(-d1τ)>d2，因此?ε>0，使得

(18)

由引理1、比較定理和式(18)得

這與式(17)矛盾，所以Ws(RE0)∩X0=?。

(19)

因為(nm-d3)(rexp(-d1τ)-d2)>ad3b，該不等式等價于(nm-d3)x21>d3b，則存在充分小的η>0 (η

(20)

由方程組(1)中第3式得

4 數(shù)值模擬

為了了解不同的參數(shù)對系統(tǒng)動力學(xué)的影響，利用Matlab軟件進行數(shù)值模擬，驗證理論結(jié)果的正確性。選取一組參數(shù)值(Ⅰ)：

r=0.2,d1=0.02,d2=0.1,d3=0.21,m=0.3,

n=0.8,a=0.01,k=0.1,b=1,τ=0.1

圖1(a)、(b)、(c)分別表示系統(tǒng)(1)在滅絕平衡點E0(0,0,0)、邊界平衡點E1(0.199 0,9.960 0,0)以及在正平衡點E*(0.001 8,5,0.722 9)處的穩(wěn)定性。

(a) d2=0.2，其余參數(shù)值同(Ⅰ)

圖1(a)選取d2=0.2，其余參數(shù)值同(Ⅰ)，此時rexp(-d1τ)-d2=-0.000 399 6<0，滿足定理5中條件，系統(tǒng)在點E0處是全局漸近穩(wěn)定的；圖1(b)選取d3=0.24，其余參數(shù)值同(Ⅰ)，此時

rexp(-d1τ)-d2=0.099 6>0

nm-d3=0，滿足定理6中的條件，系統(tǒng)在點E1處全局漸近穩(wěn)定;圖1(c)選取參數(shù)值(Ⅰ)，此時

nm-d3=0.04>0

k4-k3=0.001 2

滿足定理4中局部穩(wěn)定性條件，即系統(tǒng)在點E*處是局部漸近穩(wěn)定的。

選取另一組參數(shù)值(Ⅱ)：

r=0.5,d1=0.02,d2=0.1,d3=0.18,m=0.4,

n=0.5,a=0.01,k=0.3,b=1,τ=6

通過數(shù)值模擬得到不同的圖像， 如圖2、3所示。 圖2、3分別表示了系統(tǒng)(1)在時滯τ及恐懼因子k變化時， 對食餌種群和捕食者種群數(shù)量的影響。

(a) τ=5，其余參數(shù)同(Ⅱ) (b) 參數(shù)值(Ⅱ)

(a) k=0.2，其余參數(shù)同(Ⅱ) (b) k=0.28，其余參數(shù)同(Ⅱ)

當(dāng)系統(tǒng)(1)中幼年食餌到成年食餌的成熟時滯τ發(fā)生改變時，系統(tǒng)的穩(wěn)定性也會隨之發(fā)生改變，如圖2(a)、(b)所示。當(dāng)τ=5時，系統(tǒng)(1)在正平衡點周圍出現(xiàn)周期性震蕩，存在周期解，系統(tǒng)不穩(wěn)定；當(dāng)τ=6時，隨著τ的增加，系統(tǒng)(1)逐漸趨于穩(wěn)定。

隨著成年食餌對于捕食者的恐懼程度k的變化，系統(tǒng)(1)的動力學(xué)行為也隨之改變，如圖3(a)、(b)所示。當(dāng)k=0.2時，系統(tǒng)(1)在正平衡點處出現(xiàn)周期性震蕩，存在周期解；而當(dāng)k=0.28時，系統(tǒng)(1)逐漸趨于穩(wěn)定。

為了研究當(dāng)系統(tǒng)(1)在正平衡點附近出生率r及恐懼因素k的影響，將出生率r和恐懼因子k分別作為分岔參數(shù)，得出關(guān)于r和k的分岔圖，如圖4、5所示。

圖4中取初始值為(1,5,2)，除變量r外，其他參數(shù)與參數(shù)組(Ⅱ)相同。從圖4可以看出:對于幼年食餌出生率r，當(dāng)r>0.51時，隨著幼年食餌出生率的增加，系統(tǒng)原本的穩(wěn)定性被破壞，食餌和捕食者的數(shù)量呈周期變化;當(dāng)r∈(0.48,0.51]時，幼年食餌和捕食者的數(shù)量隨著幼年食餌出生率r的增長而增加，而成年食餌種群趨于穩(wěn)定狀態(tài)。

(a) x1(t)

圖5中取初始值為(1,5,2)，除變量k外，其他參數(shù)同參數(shù)組(Ⅱ)。

(a) x1(t)

從圖5可知：對于恐懼因子k，當(dāng)k∈(0,0.28]時，系統(tǒng)(1)中幼年食餌、成年食餌和捕食者的數(shù)量呈周期變化，系統(tǒng)不穩(wěn)定；當(dāng)k>0.28時，隨著恐懼程度的增加，食餌種群趨于穩(wěn)定狀態(tài)，恐懼效應(yīng)不再對其產(chǎn)生影響。從生物學(xué)角度上，食餌對捕食者產(chǎn)生一定程度的恐懼后，由于習(xí)慣，恐懼效應(yīng)在長期內(nèi)不再會影響食餌種群，且捕食者種群隨著k的增加而減少。

5 結(jié) 語

時，模型存在唯一的正平衡點E*。當(dāng)定理4中的條件ⅰ)成立時，E*是局部漸近穩(wěn)定的；當(dāng)正平衡點E*存在且初值滿足一定條件時，模型是一致持久的。最后，通過數(shù)值模擬驗證了理論結(jié)果的正確性及時滯、出生率或恐懼因子改變時對系統(tǒng)動力學(xué)行為的影響。