張道祥, 孫光訊, 徐明麗, 陳金瓊, 周 文

(安徽師范大學(xué) 數(shù)學(xué)與統(tǒng)計學(xué)院, 安徽 蕪湖 241002)

0 引 言

文獻(xiàn)[1]用反應(yīng)擴(kuò)散模型說明了實際應(yīng)用中的典型斑圖現(xiàn)象. 目前, 關(guān)于時空斑圖的研究已取得許多成果[2-7]： 文獻(xiàn)[2]研究了非線性捕食者收獲效應(yīng)下時滯對Turing斑圖和螺旋波斑圖形成的影響; 文獻(xiàn)[3]討論了時滯對傳染病模型中Turing斑圖生成的影響; 文獻(xiàn)[4]研究了時滯反應(yīng)擴(kuò)散系統(tǒng)中螺旋波的形成和Hopf分支產(chǎn)生周期解的穩(wěn)定性及其方向. 文獻(xiàn)[8]研究了下列帶有Holling-Ⅱ型功能反應(yīng)和雙曲死亡率的捕食者-食餌模型:

(1)

其中:U(t)和V(t)分別表示在t時刻食餌和捕食者的密度;r,K,a1,n,a2,δ均為正常數(shù), 分別表示食餌的內(nèi)稟增長率、 食餌的環(huán)境容納量、 捕食者對食餌的捕食率、 捕食者對食餌的最大攝入比率、 捕食者的增長率以及死亡率;,η是正常數(shù), 分別表示在浮游生物死亡的背景下由于水和自我保護(hù)導(dǎo)致的光衰減系數(shù). 文獻(xiàn)[8]討論了常微分系統(tǒng)的局部穩(wěn)定性和自擴(kuò)散效應(yīng)下Hopf分支的方向以及產(chǎn)生周期解的穩(wěn)定性. 本文考慮具有非線性食餌收獲效應(yīng)的捕食者-食餌的反應(yīng)擴(kuò)散系統(tǒng)[9]:

(2)

在自然界中, 時滯會導(dǎo)致系統(tǒng)的周期解由穩(wěn)定變?yōu)椴环€(wěn)定[10]; 在神經(jīng)網(wǎng)絡(luò)中, 時滯會導(dǎo)致網(wǎng)絡(luò)產(chǎn)生振蕩和不穩(wěn)定現(xiàn)象, 從而破壞網(wǎng)絡(luò)的原有性能[11]; 在化學(xué)中, 時滯影響控制系統(tǒng)的穩(wěn)定性[12]等.因此, 研究具有時滯效應(yīng)影響的數(shù)學(xué)模型具有一定的應(yīng)用價值. 為簡化模型, 本文做如下變換:

考慮時滯效應(yīng), 引入上述變換可將系統(tǒng)(2)轉(zhuǎn)化為如下時滯偏微分系統(tǒng):

(3)

其中: (X,t)∈Ω×+,Ω?2為邊界光滑的有界區(qū)域;ν為?Ω上單位外法向量;τ是正常數(shù), 表示食餌種群孕育后代引起的時間滯后效應(yīng);u(X,t-τ)表示食餌在X處t時刻的密度.

對于系統(tǒng)(3), 文獻(xiàn)[9]研究了在自擴(kuò)散效應(yīng)下時滯對Hopf分支產(chǎn)生周期解的穩(wěn)定性及其方向的影響; 文獻(xiàn)[13]討論了當(dāng)h=0,d12=d21=0時, 時滯效應(yīng)下系統(tǒng)周期解的存在性. 基于生態(tài)學(xué)中物種空間分布的重要性, 本文考慮系統(tǒng)(3)中時滯效應(yīng)和擴(kuò)散效應(yīng)對具有非線性食餌收獲效應(yīng)系統(tǒng)空間斑圖生成的影響.

1 穩(wěn)定性與分支分析

1.1 線性穩(wěn)定性分析

首先, 考慮系統(tǒng)(3)對應(yīng)的常微分系統(tǒng):

(4)

基于生態(tài)學(xué)意義, 本文只考慮系統(tǒng)(3)和系統(tǒng)(4)的正平衡點(diǎn). 通過求解易知, 當(dāng)h<ρ<β成立時, 系統(tǒng)(3)和系統(tǒng)(4)存在唯一的正平衡點(diǎn)E*=(u*,v*), 其中v*=u/(βγ)且u*滿足如下三次方程:

A0u3+A1u2+A2u+A3=0,

(5)

其中:

A0=βγ;A1=β2γ+βγρ+s-βγ;A2=β2γρ+βγh+ρs-β2γ-βγρ;A3=β2γ(h-ρ).

其次, 考慮系統(tǒng)(3)的穩(wěn)定性. 由文獻(xiàn)[14-15]知, 可假設(shè)τ是小量, 將u(X,t-τ)寫成如下形式:

(6)

將式(6)代入系統(tǒng)(3)中, 可得

(7)

其中R(u)=u/(β+u). 將R(u(X,t)-τ(?u(X,t)/?t))在(u,v)處Taylor展開, 且不考慮高階非線性部分, 則系統(tǒng)(7)轉(zhuǎn)化為

(8)

其中Ru(u)=dR(u)/du. 為簡便, 記

(9)

結(jié)合式(8)和式(9), 可得如下系統(tǒng):

(10)

下面對系統(tǒng)(10)的平衡點(diǎn)(u*,v*)進(jìn)行線性穩(wěn)定性分析. 在系統(tǒng)(10)的平衡點(diǎn)(u*,v*)處做如下微擾:

u(X,t)=u*+δu(X,t),v(X,t)=v*+δv(X,t).

(11)

對系統(tǒng)(10)在平衡點(diǎn)(u*,v*)處進(jìn)行線性化, 可得

(12)

其中:

在Neumann邊界條件下, 系統(tǒng)(10)的解有如下形式:

δu=δu*eλtei(kxx+kyy),δv=δv*eλtei(kxx+kyy),

(13)

(14)

(15)

將式(13)代入系統(tǒng)(10)中, 可得系統(tǒng)(10)在平衡點(diǎn)(u*,v*)處的特征方程:

λ2-T(k)λ+D(k)=0,

(16)

其中:

1.2 Hopf分支

若存在kH>0, 使得T(kH)=0,D(kH)>0,T(k)≠0,D(k)≠0, ?k≠kH, 則系統(tǒng)(10)發(fā)生Hopf分支. 此時矩陣A的特征值為純虛數(shù), 即tr(A)=T(kH)=0, 其中tr(A)表示A的跡.

假設(shè):

證明: 存在kj∈={0,1,2,…}, 使得當(dāng)k=kj=kH時, tr(Aj)=0. 因此, 存在一個τ的臨界值τ(i), 對于i≥0, 有

1)q(τ)≠0, 對任意的τ∈I均成立;

綜上所述, 當(dāng)τ=τc時, 系統(tǒng)(10)在平衡點(diǎn)(u*,v*)處發(fā)生了Hopf分支.

1.3 Turing分支

假設(shè):

由于時滯和擴(kuò)散可導(dǎo)致系統(tǒng)(10)的Turing不穩(wěn)定. 下面考慮系統(tǒng)(10)發(fā)生Turing不穩(wěn)定的條件, 即Turing分支的條件. 由Turing斑圖形成的必要條件知, 對于系統(tǒng)(10)的常微分形式下平衡點(diǎn)(u*,v*)是穩(wěn)定的, 而對于系統(tǒng)(10)平衡點(diǎn)(u*,v*)是不穩(wěn)定的. 由上述穩(wěn)定性分析和Turing不穩(wěn)定的條件可知, 對于系統(tǒng)(10)的常微分方程在平衡點(diǎn)(u*,v*)的Jacobi矩陣J, 有tr(J)<0和det(J)>0, 且對任意的k>0, 總有tr(A)<0成立. 若要系統(tǒng)(10)在平衡點(diǎn)(u*,v*)處是不穩(wěn)定的, 則存在某個波數(shù)k>0, 使得det(A)<0. 由此可得如下結(jié)論:

記h為分支參數(shù). 通過計算可知Turing分支參數(shù)h的臨界值hT滿足如下Turing分支曲線:

圖1 系統(tǒng)(3)的分支Fig.1 Bifurcation of system (3)

通過上述線性穩(wěn)定性理論和分支理論[16-17], 可得h-τ關(guān)系如圖1所示.圖1中紫色曲線為正解性曲線, 由Turing分支曲線(藍(lán)線)和Hopf分支曲線(紅線), 可將整個區(qū)域分成4個參數(shù)空間, 其中: 區(qū)域Ⅰ表示Turing空間, 即在常微分系統(tǒng)中, 平衡點(diǎn)是穩(wěn)定的, 而在偏微分系統(tǒng)中, 平衡點(diǎn)是不穩(wěn)定的; 區(qū)域Ⅱ表示Turing-Hopf空間, 即在常微分系統(tǒng)和偏微分系統(tǒng)下, 系統(tǒng)的平衡點(diǎn)均不穩(wěn)定; 區(qū)域Ⅲ表示由Hopf分支引起的不穩(wěn)定, 即Hopf空間; 區(qū)域Ⅳ表示在常微分系統(tǒng)和偏微分系統(tǒng)下, 平衡點(diǎn)均是穩(wěn)定的.

2 數(shù)值模擬

下面用MATLAB軟件對系統(tǒng)(3)的空間斑圖進(jìn)行數(shù)值模擬. 所有數(shù)值模擬均在離散格子200×200或400×400內(nèi)進(jìn)行, 空間間隔Δh=1, 時間間隔Δt=0.1; 所有數(shù)值模擬均采用齊次Neumann邊界條件. 由于食餌和捕食者的空間分布類似, 所以數(shù)值模擬中僅考慮食餌空間斑圖的形成. 其中一些參數(shù)取值如下:

β=0.42,ρ=0.28,γ=1.2,s=0.8,d12=0.05.

(20)

2.1 區(qū)域Ⅰ中Turing斑圖的形成

考慮在圖1區(qū)域Ⅰ中Turing斑圖的形成以及時滯τ對Turing斑圖形成的影響. 初始條件為種群的初始隨機(jī)分布, 其他參數(shù)取值為:α=1.5,β=0.8,d11=1,d21=11,d22=11.5.圖2為在不同的食餌捕獲效應(yīng)參數(shù)h下系統(tǒng)(3)的空間Turing斑圖, 其中(A),(B),(C)和(D),(E),(F)分別表示當(dāng)τ=0,1.5時的情形.圖2(A)和(D)是當(dāng)h=0.261時相同迭代步數(shù)下的Turing斑圖, 對比兩圖可見, 斑圖結(jié)構(gòu)一致, 均為點(diǎn)狀, 但圖2(D)相對圖2(A)點(diǎn)狀斑圖形成更快. 當(dāng)h值增至0.275和0.296時, 同理可見圖2(B)和(E)、
圖2(C)和(F)均為點(diǎn)狀和條狀共存的斑圖, 只是圖2(B)和(E)中點(diǎn)狀多于條狀,圖2(C)和(F)中條狀多于點(diǎn)狀, 且斑圖的結(jié)構(gòu)未發(fā)生改變, 與文獻(xiàn)[2-3]中結(jié)論一致.

2.2 區(qū)域Ⅱ中的空間斑圖

圖3為在常微分系統(tǒng)和偏微分系統(tǒng)下平衡點(diǎn)均不穩(wěn)定時, 食餌u隨時間演化的空間斑圖. 初始條件及各參數(shù)取值為:u(X,t0)=u*-0.000 2×(x-100),v(X,t0)=v*-0.000 3×(y-100),t0∈[-τ,0][7],α=1.5,h=0.1,τc=1.458 5,τ=1.65>τc,d11=0.1,d21=0.8,d22=0.95.圖3(A)為當(dāng)t=0時的斑圖, 隨著時間的演化, 出現(xiàn)規(guī)則的條狀斑圖和點(diǎn)狀斑圖共存的情形(圖3(B)～(E)), 不同的是圖3(B)～(D)中條狀優(yōu)于點(diǎn)狀, 而圖3(E)中點(diǎn)狀優(yōu)于條狀, 最后點(diǎn)狀斑圖占據(jù)整個區(qū)域(圖3(F)), 系統(tǒng)的動力學(xué)行為不再發(fā)生改變.

圖3 當(dāng)h=0.1, τ=1.65時, 食餌u隨時間演化的空間斑圖Fig.3 Spatial patterns of time evolution of prey u when h=0.1, τ=1.65

2.3 區(qū)域Ⅲ中螺旋波斑圖

考慮圖1的區(qū)域Ⅲ內(nèi)螺旋波的形成. 初始條件[4]為u(X,t0)=u*-0.000 2×(x-200),v(X,t0)=v*-0.000 3×(y-200),t0∈[-τ,0].

2.3.1 時滯對螺旋波斑圖形成的影響 圖4為系統(tǒng)(3)中食餌u的螺旋波斑圖隨時間的演化情況. 各參數(shù)取值為:α=2.15,h=0.05,d11=0.2,d21=0.03,d22=0.6, 其他參數(shù)由式(20)確定, 此時τc=1.524 9.圖4(A)為當(dāng)t=0時的斑圖, 起初螺旋波斑圖是規(guī)則的(圖4(B)～(D)), 隨著時間的演化, 由螺旋波中心開始, 其規(guī)則結(jié)構(gòu)依次消失, 逐漸演化成不規(guī)則的螺旋波斑圖(圖4(E)～(G)), 最后不規(guī)則的螺旋波斑圖占據(jù)整個區(qū)域(圖4(H)). 由圖4可見, 不規(guī)則的螺旋波斑圖由規(guī)則的螺旋波斑圖隨時間的演化而逐漸形成.

圖4 當(dāng)τ=1.9時, 食餌u隨時間演化的空間螺旋波斑圖Fig.4 Spatial spiral patterns of time evolution of prey u when τ=1.9

圖5 u(t)的時間演化曲線Fig.5 Time evolution curve of u(t)

圖5為食餌u隨時間演化的密度曲線. 由圖5可見, 隨著時間的變化, 食餌u的密度呈振蕩的變化趨勢. 因此, 在足夠長的時間內(nèi), 食餌u的空間斑圖向不規(guī)則的趨勢變化.圖6為在不同時滯τ下食餌u的空間螺旋波斑圖. 各參數(shù)取值為:α=2.15,h=0.05,d11=0.2,d21=0.4,d22=0.6, 其他參數(shù)由式(20)確定. 由圖6可見: 當(dāng)τ=2.1時, 規(guī)則的螺旋波斑圖占據(jù)整個區(qū)域(圖6(A)); 隨著τ的增大, 中間部分仍是規(guī)則的, 而邊緣出現(xiàn)不規(guī)則的螺旋波斑圖; 且隨著時滯的增大, 中間部分規(guī)則的螺旋波變得越來越細(xì), 外圍部分越來越混亂(圖6(B)～(D)).

2.3.2 非線性收獲效應(yīng)對螺旋波斑圖形成的影響

圖7為在同一時刻對系統(tǒng)(3)隨h變化的螺旋波斑圖. 各參數(shù)取值為:α=2.15,d11=0.2,d21=0.15,d22=0.85,τ=1.9, 其他參數(shù)由式(20)確定. 由圖7可見: 當(dāng)h=0.05時, 系統(tǒng)(3)產(chǎn)生的螺旋波斑圖是規(guī)則的(圖7(A)); 隨著h增至0.08和0.11, 開始出現(xiàn)不規(guī)則的螺旋波(圖7(B),(C)), 且由邊緣處依次向中間延伸破碎; 當(dāng)h增至0.16時, 完全不規(guī)則的螺旋波斑圖占據(jù)整個區(qū)域(圖7(D)); 隨著收獲效應(yīng)參數(shù)h的增大, 螺旋波的破碎進(jìn)程加快. 因此在實際應(yīng)用中, 人類適當(dāng)?shù)氖斋@可對物種進(jìn)行保護(hù)[18].

圖6 當(dāng)時間t=10 000時, 食餌u隨τ變化的空間螺旋波斑圖Fig.6 Spatial spiral patterns of prey u with different τ when t=10 000

圖7 當(dāng)時間t=10 000時, 食餌u隨h變化的空間螺旋波斑圖Fig.7 Spatial spiral patterns of prey u with different h when t=10 000

圖8 在時間t=10 000時, 食餌u隨d21變化的空間螺旋波斑圖Fig.8 Spatial spiral patterns of prey u with different d21 when t=10 000

2.3.3 交叉擴(kuò)散效應(yīng)對螺旋波斑圖形成的影響 圖8為在同一時刻系統(tǒng)(3)隨d21變化的螺旋波斑圖. 各參數(shù)取值為:α=2.15,h=0.05,d11=0.2,d22=0.6,τ=1.9, 其他參數(shù)由式(20)確定. 由圖8可見: 當(dāng)d21=0.42時, 系統(tǒng)(3)出現(xiàn)了規(guī)則的螺旋波(圖8(A)); 隨著d21減至0.29和0.12, 螺旋波斑圖逐漸變得不規(guī)則(圖8(B),(C)), 且由邊緣處依次向中間延伸破碎; 當(dāng)d21減至0.03時, 系統(tǒng)(3)出現(xiàn)了完全不規(guī)則的螺旋波斑圖(圖8(D)); 隨著擴(kuò)散系數(shù)d21的減小, 螺旋波的破碎進(jìn)程加快.

綜上, 基于時滯對生態(tài)系統(tǒng)動力學(xué)行為的影響, 本文從理論和數(shù)值兩方面研究了一類具有時滯反應(yīng)擴(kuò)散的非線性收獲效應(yīng)的捕食者-食餌模型的空間斑圖生成問題, 可得如下結(jié)論: 1) 數(shù)值結(jié)果表明, 時滯未改變系統(tǒng)(3)的Turing斑圖結(jié)構(gòu), 可能會影響Turing斑圖生成的進(jìn)程; 2) 理論與數(shù)值結(jié)果表明, 時滯和非線性收獲效應(yīng)以及交叉擴(kuò)散效應(yīng)均影響了螺旋波斑圖的生成和破碎; 3) 理論與數(shù)值結(jié)果表明, 系統(tǒng)(3)具有豐富的動力學(xué)行為.