于智博, 呂堂紅, 周林華

(長春理工大學(xué) 數(shù)學(xué)與統(tǒng)計學(xué)院， 吉林 長春 130022)

研究Lotka-Volterra模型的目的是描述和分析捕食者和獵物之間的相互作用，反映真實的自然現(xiàn)象，揭示種群動力學(xué)的本質(zhì)。生態(tài)學(xué)中的概念，如物種遷移、功能響應(yīng)函數(shù)、種內(nèi)競爭、Allee效應(yīng)、食餌避難所、恐懼效應(yīng)和時間延遲等因素被添加到Lotka-Volterra方程中，以獲得對種群數(shù)量變化和動力學(xué)行為更準(zhǔn)確的描述。在現(xiàn)實世界中，食餌尋找避難所是一種可以提高存活率的策略。避難所的存在對捕食者種群和獵物種群的共存產(chǎn)生明顯的影響，被視為影響種群動態(tài)的關(guān)鍵性因素。研究表明，隨著躲入避難所的獵物數(shù)量增加，捕食者種群的數(shù)量會先增加后減小，直到滅絕[1-2]。早期的研究工作主要使用線性獵物避難所，即有一定比例的獵物被避難所保護(hù)。此時系統(tǒng)平衡點的穩(wěn)定性與存在性取決于庇護(hù)強度[3-5]。近年來，有學(xué)者考慮更為實際的情況，即躲入避難所的獵物數(shù)量受到獵物種群數(shù)量和捕食者種群數(shù)量的影響，因此提出非線性避難所，且以非線性避難所作為控制項研究系統(tǒng)的動力學(xué)性質(zhì)，并得出豐富的結(jié)論[6-7]。

食餌種群數(shù)量的減少，除了受到捕食者種群的直接捕殺外，還受到恐懼效應(yīng)的間接影響。2016年，Wang Xiao-ying等[8]研究了如下形式的具有恐懼效應(yīng)的捕食系統(tǒng):

(1)

將追捕獵物的時間、食餌的妊娠期、捕食者的進(jìn)食時間納入經(jīng)典的Lotka-Volterra模型更符合自然界的實際情況。在生物模型中引入時滯比常微分方程具有更復(fù)雜的動力學(xué)特性，因為它可以使平衡點的穩(wěn)定性發(fā)生改變，并產(chǎn)生極限環(huán)[10-11]。受到上述研究內(nèi)容的啟發(fā)，本文在模型(1)的基礎(chǔ)上，提出了一個包括非線性避難所與雙時滯的捕食-食餌模型，并研究其動力學(xué)性態(tài)。改進(jìn)后的系統(tǒng)模型如下：

(2)

1 正平衡點的存在性

系統(tǒng)(2)正平衡點E*記做

如果系統(tǒng)(2)存在平衡點，則應(yīng)滿足

(3)

由式(3)中的第二個方程可得

(4)

將式(4)代入式(3)中的第一個方程，有

η1(x*)4+η2(x*)3+η3(x*)2+η4x*+η5=0,

(5)

其中

η1=-c3p3δ3d2-d2c3p3δ2k，

η3=3c2p2dδ2r-3c2p2δ2dd1-3cpd2d2δ-2d1c2p2kδd-d2d2cpk-2c3p3δd-c2p3dk，

η4=3rcpd2δ-3cpd2d1δ-d2d3-d1d2cpk-cp2d2，η5=(r-d1)d3。

根據(jù)Descarte’s符號準(zhǔn)則，當(dāng)下列條件成立時，至少存在一個正平衡點E*：

(H1)η1<0,η2<0,η3<0,η4>0,η5>0或η1<0,η2<0,η3<0,η4<0,η5>0。

2 局部穩(wěn)定性與Hopf分支

系統(tǒng)(2)在E*處的Jacobi矩陣為

(6)

則有系統(tǒng)(2)的特征方程：

λ2+A1λ+A2e-λτ1+(A3λ+A4)e-λτ2+A5=0,

(7)

其中A1=-(a11+a22)，A2=-a12a21′，A3=-a22′，A4=a11a22′，A5=a11a22。

情形一τ1=τ2=0。

特征方程(7)變?yōu)?/p>

λ2+(A1+A3)λ+A2+A4+A5=0。

(8)

根據(jù)Routh-Hurwitz準(zhǔn)則，若

(H2)A1+A3>0，A2+A4+A5>0,

則系統(tǒng)(2)的正平衡點E*在沒有時滯的影響下是穩(wěn)定的。

情形二τ1=0,τ2>0。

特征方程(7)變?yōu)?/p>

λ2+A1λ+A2+A5+(A3λ+A4)e-λτ2=0。

(9)

令λ=iω1(ω1>0)是該方程的根，代入式(9),可得

(10)

簡化可得

(11)

(12)

于是有如下結(jié)論。

定理1 對于系統(tǒng)(2)，當(dāng)τ1=0,τ2>0與(H1)、(H2)成立，存在一個τ10，使得當(dāng)τ2∈(0,τ10)時，E*是穩(wěn)定的；當(dāng)τ2>τ10時，E*不穩(wěn)定；當(dāng)τ2=τ10時，系統(tǒng)(2)發(fā)生Hopf分支。

情形三τ1>0,τ2=0。

證明過程類似于情形二。有

(13)

(14)

對于系統(tǒng)(2)，當(dāng)τ1>0、τ2=0與(H1)成立，當(dāng)τ1未超過τ20時，平衡點E*是穩(wěn)定的；當(dāng)τ1>τ20時，平衡點E*不穩(wěn)定；當(dāng)τ1=τ20時，系統(tǒng)(2)發(fā)生Hopf分支。

情形四τ1=τ2=τ>0。

特征方程(7)簡化為

λ2+A1λ+A5+(A3λ+A2+A4)e-λτ=0。

(15)

令λ=iω3(ω3>0)是該方程的根，代入式(15)，可得

(16)

簡化可得

(17)

(18)

剩余證明過程與情形三類似，則對于系統(tǒng)(2)，當(dāng)τ1=τ2=τ與(H1)成立，當(dāng)τ未超過臨界值τ30時，平衡點E*是穩(wěn)定的；當(dāng)τ>τ30時，平衡點E*不穩(wěn)定；當(dāng)τ=τ30時，系統(tǒng)(2)發(fā)生Hopf分支。

情形五τ1>0,τ2>0。

考慮式(2)中τ1在穩(wěn)定的區(qū)間，τ2作為參數(shù)。設(shè)λ=iω4為式(7)的根，代入式(7)有

化簡得

(19)

若(H6)式(19)至少具有一個非零解成立，則將式(19)的非零解表示為ω40。顯然

k=0,1,2,…。

對式(7)關(guān)于τ2求導(dǎo)，令τ2=τ40，經(jīng)計算有

(20)

假設(shè)(H7)R≠0，則

定理2 對于系統(tǒng)(2)，固定τ1∈(0,τ20)，若(H1)、(H6)和(H7)成立，則當(dāng)τ2∈(0,τ40)時，平衡點E*是漸近穩(wěn)定的；若τ2>τ40，平衡點E*不穩(wěn)定；當(dāng)τ2=τ40，系統(tǒng)(2)發(fā)生Hopf分支。

3 Hopf分支方向及其穩(wěn)定性

本節(jié)在τ1=τ2=τ=τ30條件下，運用Hassard[12]理論，得到確定E*(x*,y*)附近Hopf分支性質(zhì)的顯示公式。

令u(t)=(u1(t),u2(t))T∈R2，其中u1(t)=x(τt)，u2(t)=y(τt)，τ=τ30+μ，μ∈R,則當(dāng)μ=0時是系統(tǒng)(2)的Hopf分支值，在C=C([-1,0],R2)上記做

(21)

其中Lμ:C→R2,F:R×C→R2為

(22)

F(μ,φ)=(τ30+μ)(F1(μ,φ),F2(μ,φ))T,

(23)

其中

φ=(φ1,φ2)∈C([-1,0],R2)，

c11=cp-2cpδy*，c12=-d，c13=-4cpδx*，

由Riesz表示定理，可得

η(θ,μ):[-1,0]→R2,

使得

(24)

這里

其中δ(θ)是Dirac-delta函數(shù)。

對于φ∈C1([-1,0],R2)，定義

于是，系統(tǒng)(21)可改寫為

(25)

這里u=(u1,u2)，ut(θ)=u(t+θ),θ∈[-1,0]。

對于Ψ∈C1([-1,0],(R2)*)，定義A=A(0)的伴隨算子A*：

對于φ∈C1([-1,0],R2)，Ψ∈C1([-1,0],(R2)*)，為了使算子A和伴隨算子A*的特征向量規(guī)范化，可以定義以下的雙線性內(nèi)積:

設(shè)q(θ)和q*(s)是對應(yīng)于特征根iω3τ30與-iω3τ30的特征向量。于是

A(0)q(θ)=iω3τ30q(θ),A*(0)q*(s)=-iω3τ30q*(s),

計算得

這里

令Xt為μ=0時方程(25)的解。定義

z(t)=〈q*,Xt〉,

(26)

在中心流型C0上，有

(27)

(28)

在由式(27)和(28)得

Xt(θ)=W(t,θ)+2Re{z(t)q(θ)},

綜合式(23)得

比較系數(shù)可得

為了確定g21,下面計算W20(θ)和W11(θ)。由式(25)和(26)得到

即

(29)

其中

(30)

結(jié)合式(26)和(27)得到

其中

(31)

通過比較系數(shù)得到

(2iω0-A)W20(θ)=H20(θ), -AW11(θ)=H11(θ)。

當(dāng)θ∈[-1,0)時

(32)

比較式(31)和(32)有

其中W20(θ)、W11(θ)的結(jié)果如下：

可得

(33)

(34)

這里，C1(0)由式(33)給出，易得出μ2、β2、T2的值。因此有了定理3。

定理3 當(dāng)τ=τ30時，由式(34)可得以下結(jié)論：

(1)μ2確定Hopf分支的方向，如果μ2>0，則Hopf分支是超臨界的，反之是次臨界的；

(2)β2確定分支周期解的穩(wěn)定性，如果β2<0，則周期解是穩(wěn)定的，反之是不穩(wěn)定的；

(3)T2確定分支周期解的周期，如果T2>0，則周期解的周期是增加的，反之是減少的。

4 數(shù)值模擬

為驗證情形三、情形五的結(jié)論(情形二與情形三類似)，考慮參數(shù)如下：

(35)

根據(jù)Descarte’s符號準(zhǔn)則，系統(tǒng)(35)有唯一平衡點E*(x*,y*)=E*(17.902 3,3.919 1)。當(dāng)τ1=τ2=0時，兩種群密度的數(shù)量變化趨于穩(wěn)定。

當(dāng)τ1>0,τ2=0時，τ20=4.67。當(dāng)τ1=4.4<τ20=4.67時，種群密度保持穩(wěn)定，見圖1；當(dāng)τ1=6>τ20=4.67時，系統(tǒng)出現(xiàn)周期解，見圖2。

當(dāng)τ1>0,τ2>0時，取τ1=3，τ*=3.849 6。當(dāng)τ2=2.0<τ*=3.849 6時，平衡點E*是穩(wěn)定的，兩種群數(shù)量最終趨向于平衡點E*,如圖3所示；當(dāng)τ2=5>τ*=3.849 6時，平衡點E*附近出現(xiàn)周期解，如圖4所示。

(a) x的變化 (b) y的變化 (c) 系統(tǒng)的變化圖1 當(dāng)τ1=4.4<τ20=4.67,τ2=0時，系統(tǒng)(35)的動力學(xué)行為

(a) x的變化 (b) y的變化 (c) 系統(tǒng)的變化圖2 當(dāng)τ1=6>τ20=4.67,τ2=0時，系統(tǒng)(35)的動力學(xué)行為

(a) x的變化 (b) y的變化 (c) 系統(tǒng)的變化圖3 當(dāng)τ1=3,τ2=2.0<τ*=3.849 6時，系統(tǒng)(35)的動力學(xué)行為

(a) x的變化 (b) y的變化 (c) 系統(tǒng)的變化圖4 當(dāng)τ1=3,τ2=5>τ*=3.849 6時，系統(tǒng)(35)的動力學(xué)行為

為驗證情形四的結(jié)論，選擇參數(shù)：r=0.4,d1=0.05,d2=0.01,p=0.1,δ=0.04,c=0.7,d=0.07,k=0.03。此時有η1<0,η2<0,η3<0,η4<0,η5>0。根據(jù)Descarte’s符號準(zhǔn)則，至少存在一個正平衡點E*。

當(dāng)τ1=τ2=τ時，ω3=0.284 5，τ30=4.15，C1(0)=-1.160 5-2.147 7i，μ2=-71.24，β2=-2.53，T2=-7.54。

當(dāng)τ=4.0<τ30=4.15時，平衡點E*是穩(wěn)定的,見圖5；當(dāng)τ=5.1>τ30=4.15時，平衡點E*失去穩(wěn)定性，見圖6。

(a) x的變化 (b) y的變化 (c) 系統(tǒng)的變化圖5 當(dāng)τ1=τ2=τ=4.0<τ30=4.15時，系統(tǒng)(35)的動力學(xué)行為

(a) x的變化 (b) y的變化 (c) 系統(tǒng)的變化圖6 當(dāng)τ1=τ2=τ=5.1>τ30=4.15時，系統(tǒng)(35)的動力學(xué)行為

表1是在情形四的條件下避難所δ的變化對Hopf分支臨界值的影響。在情形四的條件下，只改變δ的值，而其余參數(shù)不變。研究發(fā)現(xiàn)通過增大δ的值，時滯也會逐漸增大。意味著庇護(hù)所的存在，對種群數(shù)量的穩(wěn)定起到積極作用。

表1 在情形四的條件下避難所δ的變化對Hopf分支臨界值的影響

5 結(jié)論

本文研究一類具有避難所和雙時滯的捕食系統(tǒng)Hopf分支問題，分析了時滯在不同情況下系統(tǒng)正平衡點的穩(wěn)定性狀況。理論分析表明，時滯的存在是決定系統(tǒng)發(fā)生Hopf分支的充分條件。其次選擇恰當(dāng)?shù)膮?shù)，驗證了理論分析的正確性。從生物學(xué)的角度來看，在種群的初始密度固定時，當(dāng)反饋時滯很小(小于臨界值時)，系統(tǒng)將維持穩(wěn)定的狀態(tài)，兩種群的數(shù)量最終趨向一個固定的數(shù)值；當(dāng)反饋時滯(大于臨界值時)，捕食者和獵物種群之間的動態(tài)平衡是呈現(xiàn)出一種相互限制和周期性循環(huán)的形式。此外，數(shù)值分析表明，避難所系數(shù)的增大，能夠改變時滯的臨界值，從而增強系統(tǒng)的穩(wěn)定性。然而，本文并沒有考慮恐懼效應(yīng)對系統(tǒng)穩(wěn)定性的影響，這一方面的研究可以在未來的工作中進(jìn)行補充。