近年来, 恐惧效应对生物种群的影响引起众多学者的关注. Zanette等^[1]通过对鸣禽(Songbirds)的实验研究发现, 在没有直接掠杀的情况下, 仅捕食者的恐惧效应就可导致成年鸣禽每年繁殖后代的数量减少百分之四十. 实验结果还表明, 仅捕食者给食饵带来的掠食风险感知就能在一定程度上影响食饵的行为和生理特征, 并引发连锁反应, 其影响程度甚至超过直接捕杀. 随后, Clinchy和Zanette等^[2]提出了适用于野生动物的测量"恐惧神经回路"的实验室技术. 该实验研究表明, 食饵暴露在捕食者面前可以引起食饵的持续心理压力, 以至影响食饵种群规模的大小. 此外, 恐惧效应还存在于许多陆地脊椎动物种群中, 如鸟类^[3-10]、麋鹿^[11]、雪鞋兔^[12]以及儒艮^[13].

Wang和Zanette等^[14]为了考虑恐惧效应对种群的影响, 首次将恐惧效应引入到捕食-食饵模型中. 基于实验数据, 在食饵出生率项上增加恐惧效应比例系数, 建立了如下具有恐惧效应的生态动力系统:

(1.1) $ \begin{equation} \left\{ \begin{array}{l} { } \frac{{\rm d}u}{{\rm d}t} = ur_{0}f(k_{0}, v)-du-au^{2}-g(u)v, \\ { } \frac{{\rm d}v}{{\rm d}t} = v(-m+cg(u)), \end{array}\right. \end{equation} $

其中$ u, v $分别表示$ t $时刻时食饵和捕食者的种群数量或密度, $ r_{0}, d, a $分别表示食饵的出生率、死亡率和拥挤效应参数, $ m $表示捕食者的死亡率, $ c $表示捕食者与食饵之间的作用强度, $ g(u) $是功能反应函数, $ k_{0} $表示恐惧效应参数, $ f(k_{0}, v) = \frac{1}{1+k_{0}v} $反映捕食者对食饵的恐惧效应, 其中参数$ k_{0} $越大说明捕食者给食饵带来的恐惧效应越大. 作者首先讨论了恐惧效应对动力系统中种群分布的影响. 然后, 当系统(1.1)的功能反应函数为Holling Ⅱ型时, 作者又给出了系统平衡点的稳定性, Hopf分支的存在性及其分支方向的结果.

随后, 带恐惧效应的捕食-食饵系统引起学者的广泛关注. 闫建博和刘贤宁^[15]研究了具有Beddington-DeAngelis功能反应函数的捕食系统; 夏青艳和张睿^[16]研究了一类具有群体防御能力的Ⅰ类功能反应的捕食系统; Pal等人^[17]研究了一类具有Leslie-Gower型的捕食系统. 进一步地, 考虑空间因素, 戴斌祥和孙光讯^[18]研究了带恐惧效应的扩散捕食-食饵系统的Turing斑图. Wang Xiaoying和邹幸福^[19]提出并研究了一个反应扩散对流捕食- 食饵模型, 并将恐惧效应纳入模型的局部反应项中. 陈珊珊^[20]研究了如下的反应扩散系统

(1.2) $ \begin{equation} \left\{ \begin{array}{ll} { } \frac{\partial u}{\partial t} = d_{1}\Delta u+\frac{ru}{1+kv}-du-au^{2}-\frac{buv}{1+qu}, {\quad} x\in\Omega, t>0, \\ { } \frac{\partial v}{\partial t} = d_{2}\Delta v-m_{1}v-m_{2}v^{2}+\frac{cuv}{1+qu}, {\quad} x\in\Omega, t>0, \\ { } \frac{\partial u}{\partial \nu} = \frac{\partial v}{\partial \nu} = 0, {\quad} x\in\partial\Omega, t>0, \\ { } u(x, 0) = u_{0}(x)\geq(\not\equiv)0, v(x, 0) = v_{0}(x)\geq(\not\equiv)0, \end{array}\right. \end{equation} $

其中$ d_{1}, d_{2} $分别为食饵和捕食者的扩散系数, $ r, d, a, k $分别表示食饵的内禀增长率、死亡率, 拥挤效应参数和恐惧效应参数, $ m_{1}, m_{2} $分别表示捕食者的死亡率与拥挤效应参数, $ b, c $分别表示捕食者与食饵之间的作用强度, $ q $表示食饵逃避攻击的能力, 且以上所有参数均为正常数. $ (x, t)\in \Omega\times {{\Bbb R}}^{+} $, $ \Omega\subset {{\Bbb R}}^{2} $为有界区域, $ \nu $为$ \partial\Omega $上单位外法向量. 作者研究了这个具有Holling Ⅱ型功能反应函数的空间扩散系统, 利用Lyapunov函数等理论证明了在大转化率下不存在非恒定正稳态解.

数值模拟能够直观地展示捕食-食饵系统的动力学现象, 如平衡点的稳定性、Hopf分支和空间斑图等. 李海银^[21]考虑了密度制约且具有Beddington-DeAngelis功能反应函数的时滞捕食-被捕食系统, 利用数值模拟对系统稳定性变换和系统的Hopf分支进行了验证. 赵洪涌和张学兵等^[22]提出了一种带时滞的生物经济系统, 通过理论和数值模拟展示了系统具有丰富的动力学行为. Pijush等人^[23]研究了具有恐惧效应的捕食者-食饵模型, 并数值模拟了时滞使得系统动力学在稳定焦点和极限周期振荡之间多次切换, 对于更高的延迟参数值, 系统最终会进入混乱状态. 张道祥和孙光讯^[24]研究了一类带有Holling Ⅲ型功能反应和线性收获效应的时滞扩散捕食者-食饵系统的空间动力学, 利用数值模拟展示了系统具有丰富的动力学行为.

基于以上讨论, 本文将从理论和数值模拟两个方面研究如下带有恐惧效应和时滞的Holling Ⅲ型反应扩散捕食- 食饵系统

(1.3) $ \begin{equation} \left\{ \begin{array}{ll} { } \frac{\partial u(x, t)}{\partial t} = d_{1}\Delta u(x, t)+\frac{ru(x, t)}{1+k_{0} v(x, t)}- du(x, t)-au^{2}(x, t)-\frac{bu^{2}(x, t)v(x, t)}{1+qu^{2}(x, t)}, \ x\in\Omega, t>0, \\ { }\frac{\partial v(x, t)}{\partial t} = d_{2}\Delta v(x, t)-m_{1}v(x, t)-m_{2}v^{2}(x, t)+\frac{cu^{2}(x, t-\tau)v(x, t)}{1+qu^{2}(x, t-\tau)}, {\quad} x\in\Omega, t>0, \\ { } \frac{\partial u(x, t)}{\partial \nu} = \frac{\partial v(x, t)}{\partial \nu} = 0, {\quad} x \in\partial\Omega, t>0, \\ u(x, t) = u_{0}(x, t)\geq0, v(x, t) = v_{0}(x, t)\geq0, {\quad} t\in [-\tau, 0], x\in \Omega, \end{array}\right. \end{equation} $

其中$ u(x, t), v(x, t) $分别表示食饵和捕食者在$ x $处$ t $时刻的数量或密度. $ u(x, t-\tau) $表示食饵在$ x $处$ t-\tau $时刻的数量或密度, 正常数$ \tau $为食饵种群孕育后代引起的时间滞后效应参数. 正常数$ k_{0} $表示捕食者对食饵产生的恐惧效应参数. 其余参数的生物学意义与系统(1.2)中的一致.

系统(1.3)对应的常微分方程(ODE)系统为

(1.4) $ \begin{equation} \left\{ \begin{array}{l} { }\frac{{\rm d}u(t)}{{\rm d}t} = \frac{ru(t)}{1+k_{0}v(t)}-{\rm d}u(t)-au^{2}(t)-\frac{bu^{2}(t)v(t)}{1+qu^{2}(t)}, t>0, \\ { } \frac{{\rm d}v(t)}{{\rm d}t} = -m_{1}v(t)-m_{2}v^{2}(t)+\frac{cu^{2}(t-\tau)v(t)}{1+qu^{2}(t-\tau)}, t>0, \\ u(t) = u_{0}(t)\geq0, v(t) = v_{0}(t)\geq0, t\in [-\tau, 0]. \end{array}\right. \end{equation} $

本文亮点工作有: (1)首先将时滞效应引入到系统(1.3), 从而首次探讨了恐惧效应对具有Holling Ⅲ型功能反应函数的生态模型Hopf分支的影响. (2)与已有的文献相比, 本文重点展现了恐惧效应$ k_{0} $与系统稳定区间$ [0, \tau_{0}) $的非线性关系的有趣现象($ \tau_{0} $为Hopf分支临界值). 本文结构如下: 第2节中, 首先给出正平衡点存在的条件, 并通过特征方程得出正平衡点局部渐近稳定的条件. 第3节中, 通过线性稳定性分析给出了存在Hopf分支的条件, 并计算出了Hopf分支点. 第4节中, 通过数值模拟, 证实了恐惧效应影响系统正平衡点的稳定性, 并通过研究恐惧效应与系统稳定区间的关系, 发现了恐惧效应的增加使系统由稳定变为不稳定, 再变为稳定的现象. 最后给出了本文的结论.

考虑到系统(1.3)的生态学意义, 我们只研究系统(1.3)的正平衡点, 令

$ \left\{ \begin{array}{l} { } f_{11}(u, v) = \frac{ru}{1+k_{0}v}-du-au^{2}-\frac{bu^{2}v}{1+qu^{2}}, \\ { } f_{12}(u, v) = -m_{1}v-m_{2}v^{2}+\frac{cu^{2}v}{1+qu^{2}}, \end{array} \right. $

则系统的平衡点满足如下方程组

$ \left\{ \begin{array}{l} f_{11}(u, v) = 0, \\ f_{12}(u, v) = 0. \end{array} \right. $

定理2.1  若满足条件

$ H_{1}:c-m_{1}q>0, r-d-a\sqrt{\frac{m_{1}}{c-m_{1}q}}>0, $

系统(1.3)至少存在一个正平衡点$ E_{\star} = (u_{\star}, v_{\star}) $.

证  由方程$ f_{12}(u, v) = 0 $, 我们得出

$ v_{\star} = \frac{cu^{2}_{\star}}{m_{2}(1+qu^{2}_{\star})}-\frac{m_{1}}{m_{2}}. $

将其代入$ f_{11}(u, v) = 0 $, 得到

(2.1) $ \begin{equation} a_{7}u_{\star}^{7}+a_{6}u_{\star}^{6}+a_{5}u_{\star}^{5}+a_{4}u_{\star}^{4}+a_{3}u_{\star}^{3}+a_{2}u_{\star}^{2}+a_{1}u_{\star}+a_{0} = 0, \end{equation} $

其中

$ \begin{eqnarray*} &&a_{0} = rm^{2}_{2}-dm_{2}^{2}+k_{0}dm_{1}m_{2};\\ &&a_{1} = -am_{2}^{2}+ak_{0}m_{1}m_{2}-bk_{0}m_{1}^{2}+bm_{1}m_{2};\\ &&a_{2} = 3rqm_{2}^{2}-dk_{0}cm_{2}-3dqm_{2}^{2}+3dqk_{0}m_{1}m_{2};\\ &&a_{3} = -ak_{0}cm_{2}-3aqm_{2}^{2}+3aqk_{0}m_{1}m_{2}-2bqk_{0}m_{1}^{2}+2bck_{0}m_{1}-bcm_{2}+2bqm_{1}m_{2};\\ && a_{4} = 3rq^{2}m^{2}_{2}-dcqk_{0}m_{2}-3dq^{2}m^{2}_{2}+3dq^{2}k_{0}m_{1}m_{2}-dqk_{0}m_{2};\\ &&a_{5} = aqm_{2}(-2ck_{0}-3qm_{2}+3qk_{0}m_{1})-bc^{2}k_{0}+bq(-qk_{0}m_{1}^{2}+2k_{0}cm_{1}-cm_{2}+m_{1}m_{2}q);\\ &&a_{6} = rm_{2}^{2}q^{3}-dq^{2}m_{2}(k_{0}c+qm_{2}-qm_{1}k_{0});\\ &&a_{7} = -aq^{2}m^{2}_{2}(k_{0}(c-qm_{1})+qm_{2}). \end{eqnarray*} $

为了使$ v_{\star}>0 $, 即要使$ (c-m_{1}q)u_{\star}^{2}>m_{1} $, 再由$ H_{1} $知$ c-m_{1}q>0 $, 则需要$ u_{\star}>\sqrt{\frac{m_{1}}{c-m_{1}q}}: = u_{\star \min} $.

下面证明方程(2.1)至少存在一个正解$ u_{\star} $使得$ u_{\star}>u_{\star \min} $.

令函数$ F(u) = a_{7}u^{7}+a_{6}u^{6}+a_{5}u^{5}+a_{4}u^{4}+a_{3}u^{3}+a_{2}u^{2}+a_{1}u+a_{0}. $

当$ u\rightarrow +\infty $时, $ F(u) $与最高次项系数$ a_{7} $同号. 而

$ a_{7} = -aq^{2}m^{2}_{2}(k_{0}(c-qm_{1})+qm_{2})<0, $

故$ F(+\infty)<0. $

当$ u = u_{\star \min} $时, 有

$ F(u_{\star min}) = (1+q\frac{m_{1}}{c-m_{1}q})^{3}m_{2}^{2}(r-d-a\sqrt{\frac{m_{1}}{c-m_{1}q}}), $

显然在$ H_{1} $的条件下$ F(u_{\star \min})>0. $

由根的存在性定理知, 存在$ u_{\star}\in(u_{\star \min}, +\infty) $, 使得$ F(u_{\star}) = 0 $. 故当满足条件$ H_{1} $时, 系统(1.3)至少存在一个正平衡点$ E_{\star} = (u_{\star}, v_{\star}) $.

接下来的讨论, 我们均假设系统(1.3)已存在正平衡点$ E_{\star} $.

为了讨论系统(1.3)在正平衡点$ E_{\star} $处的稳定性，我们不妨先定义

$ U(x, t) = u(x, t)-u_{\star}, V(x, t) = v(x, t)-v_{\star}, $

系统(1.3)在正平衡点$ E_{\star} = (u_{\star}, v_{\star}) $处的线性化方程为

(2.2) $ \begin{equation} \left. \begin{array}{l} \left( \begin{array}{cc} { } \frac{\partial U(x, t)}{{\rm d}t}\\ { } \frac{\partial V(x, t)}{{\rm d}t} \end{array} \right) = {D}\left( \begin{array}{cc} \Delta U(x, t)\\ \Delta V(x, t) \end{array} \right) +{J_{1}}\left( \begin{array}{cc} U(x, t)\\ V(x, t) \end{array} \right) +{J_{2}}\left( \begin{array}{cc} U(x, t-\tau)\\ V(x, t-\tau) \end{array} \right), \\ \end{array} \right. \end{equation} $

其中$ D = \left( \begin{array}{cc} d_{1}&0\\ 0&d_{2} \end{array} \right), $$ J_{1} = \left( \begin{array}{cc} a_{11}&a_{12}\\ 0&a_{22} \end{array} \right) $和$ J_{2} = \left( \begin{array}{cc} 0&0\\ a_{21}&0 \end{array} \right). $

$ a_{11} = \frac{bu_{\star}v_{\star}(qu^{2}_{\star}-1)}{(1+qu^{2}_{\star})^{2}}-au_{\star}; a_{12} = -\frac{rk_{0}u_{\star}}{(1+k_{0}v_{\star})^{2}}-\frac{bu^{2}_{\star}}{1+qu^{2}_{\star}}; a_{21} = \frac{2cu_{\star}v_{\star}}{(1+qu^{2}_{\star})^{2}}; a_{22} = -m_{2}v_{\star}. $

则系统(2.2)的特征方程为

(2.3) $ \begin{equation} \lambda^{2}+\delta_{1k}\lambda+\delta_{2k}+\delta_{3}e^{-\lambda\tau} = 0, \end{equation} $

其中$ k $是波数, 且$ k\geq0 $,

$ \delta_{1k} = -a_{11}-a_{22}+d_{1}k^{2}+d_{2}k^{2}; \delta_{2k} = d_{1}d_{2}k^{4}-(a_{11}d_{2}+a_{22}d_{1})k^{2}+a_{11}a_{22}; \delta_{3} = -a_{12}a_{21}. $

我们作出如下假设:

$ H_{2}:a_{11}+a_{22}<0; $

$ H_{3}:a_{11}d_{2}+a_{22}d_{1}<0; $

$ H_{4}:a_{11}a_{22}-a_{12}a_{21}>0. $

定理2.2  对任意$ k\geq0 $, 当$ \tau = 0 $时, 若条件$ H_{2}-H_{4} $成立, 则系统(1.3)的正平衡点$ E_{\star} $是局部渐近稳定的.

证  当$ \tau = 0 $时, 特征方程(2.3)变为

(2.4) $ \begin{equation} \lambda^{2}+\delta_{1k}\lambda+\delta_{2k}+\delta_{3} = 0, \end{equation} $

我们把上述特征方程(2.4)的两个根分别设为$ \lambda_{1}, \lambda_{2} $, 对任意$ k\geq0 $, 有

$ \lambda_{1}+\lambda_{2} = -\delta_{1k}, $

$ \lambda_{1}\lambda_{2} = \delta_{2k}+\delta_{3} = d_{1}d_{2}k^4-(a_{11}d_{2}+a_{22}d_{1})k^2+a_{11}a_{22}-a_{12}a_{21}. $

当条件$ H_{2}-H_{4} $成立时, 我们可以得到$ \lambda_{1}+\lambda_{2}<0, \lambda_{1}\lambda_{2}>0. $即方程(2.4)有两个具有负实部的根(两个负实根或一对具有负实部的共轭复根), 因此系统(1.3)的正平衡点$ E_{\star} $是局部渐近稳定的.

注2.1  结合上述定理可知, 在不考虑时滞的情况下, 由于条件$ H_{4} $成立与否受$ k_{0} $的影响, 因此恐惧效应影响着系统(1.3)的正平衡点$ E_{\star} $的稳定性.

本节中, 将时滞作为分支参数, 使用线性稳定性理论分析捕食-食饵系统的Hopf分支问题. 设$ {\rm i}\omega $是方程(2.3)的纯虚根, 则有

$ -\omega^{2}+{\rm i}\delta_{1k}\omega+\delta_{2k}+\delta_{3}(\cos\omega\tau-{\rm i}\sin\omega\tau) = 0, $

分离上述方程的实部和虚部, 我们推出

(3.1) $ \begin{equation} -\omega^2+\delta_{2k} = -\delta_{3}\cos(\omega\tau), \end{equation} $

(3.2) $ \begin{equation} \delta_{1k}\omega = \delta_{3}\sin(\omega\tau), \end{equation} $

将方程(3.1), 方程(3.2)的左右两边平方再相加, 得到

(3.3) $ \begin{equation} w^4+(\delta_{1k}^2-2\delta_{2k})w^2+\delta_{2k}^{2}-\delta_{3}^{2} = 0, \end{equation} $

令$ \omega^{2} = z, $则方程(3.3)变为

(3.4) $ \begin{equation} z^{2}+(\delta_{1k}^2-2\delta_{2k})z+\delta_{2k}^{2}-\delta_{3}^{2} = 0, \end{equation} $

我们作出如下假设:

$ H_{5}:a_{11}d_{1}+a_{22}d_{2}<0; $

$ H_{6}:a_{11}a_{22}+a_{12}a_{21}>0; $

$ H_{7}:a_{11}a_{22}+a_{12}a_{21}<0. $

定理3.1  对任意$ \tau>0 $, $ k>0 $时, 若条件$ H_{3}-H_{6} $成立, 则系统(1.3)的正平衡点$ E_{\star} $是局部渐近稳定的.

证  我们设方程(3.4)有两个根$ z_{1k}, z_{2k} $, 由韦达定理, 得到

$ z_{1k}+z_{2k} = 2\delta_{2k}-\delta_{1k}^2 = -d_{1}^{2}k^{4}-d_{2}^{2}k^{4}-a_{11}^{2}-a_{22}^{2}+2(a_{11}d_{1}+a_{22}d_{2})k^{2}; $

$ z_{1k}z_{2k} = (\delta_{2k}+\delta_{3})(\delta_{2k}-\delta_{3}) = (\delta_{2k}+\delta_{3})[d_{1}d_{2}k^4-(a_{11}d_{2}+a_{22}d_{1})k^2+a_{11}a_{22}+a_{12}a_{21}]. $

当条件$ H_{3}, H_{4}, H_{6} $成立时, 有$ \delta_{2k}+\delta_{3}>0, \delta_{2k}-\delta_{3}>0 $, 即$ z_{1k}z_{2k}>0 $. 当$ H_{5} $成立时, $ z_{1k}+z_{2k}<0 $.

综上所述, 当条件$ H_{3}-H_{6} $成立时, 方程(3.4)的两个根都具有负实部, 所以系统(1.3)的正平衡点$ E_{\star} $是局部渐近稳定的.

注3.1  条件$ H_{6} $成立与否受到$ k_{0} $的影响, 因此当考虑时滞时, 恐惧效应参数$ k_{0} $影响着正平衡点$ E_{\star} $的稳定性.

定理3.2  当$ k = 0 $, 条件$ H_{4}, H_{7} $成立时, 方程(3.4)有唯一的正实根$ z_{0} $, 此时

$ z_{0} = \frac{2\delta_{20}-\delta_{10}^{2}+\sqrt{(\delta_{10}^{2}-2\delta_{20})^2-4(\delta_{20}^{2}-\delta_{3}^{2})}}{2}. $

证  当$ k = 0 $时, 方程(3.4)为

(3.5) $ \begin{equation} z^{2}+(\delta_{10}^2-2\delta_{20})z+\delta_{20}^{2}-\delta_{3}^{2} = 0, \end{equation} $

方程(3.5)的两个根$ z_{1} $, $ z_{2} $的关系为

$ z_{1}+z_{2} = -a_{11}^{2}-a_{22}^{2}<0, $

$ z_{1}z_{2} = (a_{11}a_{22}-a_{12}a_{21})(a_{11}a_{22}+a_{12}a_{21}), $

当$ H_{4}, H_{7} $成立, $ z_{1}z_{2}<0 $, 方程(3.5)有唯一的正实根$ z_{0} $.

由二次函数求根公式, 方程(3.5)的正根为

$ z_{0} = \frac{2\delta_{20}-\delta_{10}^{2}+\sqrt{(\delta_{10}^{2}-2\delta_{20})^2-4(\delta_{20}^{2}-\delta_{3}^{2})}}{2}. $

又由于$ \omega^{2} = z $, 则$ \omega_{0} = \sqrt{z_{0}}. $

定理3.3  当$ k = 0 $, $ \tau = \tau_{0}^{j} = \frac{1}{\omega_{0}}(\arccos(\frac{\omega_{0}^{2}-\delta_{20}}{\delta_{3}})+2j\pi), j = 0, 1, 2, \cdots $时, 系统(1.3)在$ \tau = \tau_{0}^{j} $处发生Hopf分支.

证  当$ k = 0 $时, 由方程(3.1), 我们可以得到

$ \cos(\omega_{0}\tau) = \frac{\omega_{0}^2-\delta_{20}}{\delta_{3}}, $

因此有

$ \tau_{0}^{j} = \frac{1}{\omega_{0}}(\arccos(\frac{\omega_{0}^2-\delta_{20}}{\delta_{3}})+2j\pi), j = 0, 1, 2, \cdots . $

易知, 数列$ \{\tau_{0}^{j}(k_{0})\}_{j = 0}^{\infty} $是单调递增的, $ \tau_{0}^{0} = \min{\{\tau_{0}^{j}, j = 0, 1, 2, \cdots \}}: = \tau_{0} $.

接下来, 我们讨论系统(1.4)的Hopf分支的横截性条件, 即式$ {\frac{{\rm dRe}\lambda(\tau)}{{\rm d}\tau}}|_{\tau = \tau_{0}^{j}} $的符号. 假设$ \lambda(\tau) = \alpha(\tau)+{\rm i}\omega(\tau) $是方程(2.3)的根, 此处$ \alpha(\tau_{0}^{j}) = 0 $, $ \omega(\tau_{0}^{j}) = \omega_{0} $. 对方程(2.3) 的两边对$ \tau $进行求导, 得到

$ (2\lambda+\delta_{10}-\delta_{3}\tau e^{-\lambda\tau})\frac{{\rm d}\lambda(\tau)}{{\rm d}\tau} = \lambda\delta_{3}e^{-\lambda\tau}, $

则有

(3.6) $ \begin{equation} (\frac{{\rm d}\lambda(\tau)}{{\rm d}\tau})^{-1} = \frac{(2\lambda+\delta_{10})e^{\lambda\tau}}{\lambda\delta_{3}}-\frac{\tau}{\lambda}, \end{equation} $

将$ \tau_{0}^{j} $代入方程(3.6)中, 可以得到

$ {\rm Re}\{(\frac{{\rm d}\lambda}{{\rm d}\tau})^{-1}|_{\tau = \tau_{0}^{j}}\} = \frac{\sqrt{(\delta_{10}^{2}-2\delta_{20})^{2}-4(\delta_{20}^{2}-\delta_{3}^{2})}}{\delta_{3}^{2}}, $

显然, $ {\rm Re}\{(\frac{{\rm d}\lambda}{{\rm d}\tau})^{-1}|_{\tau = \tau_{0}^{j}}\}>0. $

定理3.4  对于系统(1.3), 假设条件$ H_{1}-H_{4} $, $ H_{7} $成立,

$ \rm (1) $当$ \tau\in[0, \tau_{0}), $系统(1.3)的正平衡点$ E_{\star} $是局部渐近稳定的;

$ \rm (2) $当$ \tau\in(\tau_{0}, +\infty) $, 系统(1.3)的正平衡点$ E_{\star} $是不稳定的.

注3.2  由上述定理可知, 恐惧效应参数$ k_{0} $影响着系统(1.3)正平衡点的稳定区间$ [0, \tau_{0}) $.

在这一节里, 我们使用MATLAB软件进行数值模拟, 以此来验证理论推导所得到的关于系统(1.3)的结论. 在所有数值算例中, 由于食饵和捕食者的空间斑图是类似的, 所以仅给出捕食者$ v $的空间斑图. 参数取为$ r = 6, d = 0.2, a = 1, b = 5, q = 2, c = 2. 5, m_{1} = 0.4, $$ m_{2} = 1, d_{1} = 0. 01, d_{2} = 0.1. $

当系统(1.3)的恐惧效应参数$ k_{0} = 3 $时, 通过计算得到系统的正平衡点为$ E_{\star} = (1.2469, $$ 0.5458), $ Hopf分支临界点为$ \tau_{0} = 1.9491 $, 并且这些参数满足定理3.4的条件$ H_{1}-H_{4} $和$ H_{7} $. 这样由定理3.4可知, 在上述参数不变的情况下, 当$ \tau\in[0, \tau_{0}) $时, 系统的正平衡点$ E_{\star} $将是局部渐近稳定的, 而当$ \tau\in(\tau_{0}, +\infty) $时, 系统的正平衡点$ E_{\star} $将变得不稳定.

图 1和图 2显示的分别是时滞$ \tau $取不同值时捕食者和食饵随时间的密度演化图和相图. 图 1中, 子图(a)和子图(b)分别表示时滞为$ \tau = 1.4\in[0, \tau_{0}) $时的食饵和捕食者的密度演化图, 图中食饵和捕食者在$ t>500 $后的密度值几乎不变, 且数值均为平衡时的值, 即$ u = 1.2469, v = 0.5458, $这表明系统的正平衡点$ E_{\star} $是渐近稳定的. 子图(c)和子图(d)分别表示时滞为$ \tau = 2\in(\tau_{0}, +\infty) $时的食饵和捕食者的密度演化图, 图中在$ t>1500 $时，食饵和捕食者的密度均呈振幅相同的波状, 这符合Hopf分支现象, 系统正平衡点$ E_{\star} $是不稳定的. 这验证了定理3.4的结论. 图 2中子图(a) 表示时滞为$ \tau = 1.4\in[0, \tau_{0}) $时捕食者的相图, 相图表明无论时间为多少, 食饵和捕食者的密度值不变, 即其相图呈现为一个点(1.2469, 0.5458), 且该点即为正平衡点$ E_{\star} $, 这说明系统的正平衡点$ E_{\star} $是稳定的. 图 2中子图(b)为时滞$ \tau = 2\in(\tau_{0}, +\infty) $时, 呈现为极限环, 这证明$ \tau_{0} $附近会出现Hopf分支现象, 进一步证实了理论结果.

图 3和图 4展示的是捕食者的空间斑图的演化过程, 其中各子图右侧的条状图可作为密度值的参照. 图 3 ($ \tau = 1.4<\tau_{0} $)中的四个子图表示的时刻分别为$ t = 200, t = 400, $$ t = 600, t = 20000 $. 从图中可以看出, 由于初始扰动, 当$ t<400 $时出现密度分布不均的现象, 但在$ t = 600 $时系统密度分布变得均匀, 且持续到$ t = 20000 $时密度分布仍均匀, 且参照密度分布参照条可知, 捕食者的密度约为0.54582, 与正平衡点$ v_{\star} = 0.5458 $相同, 进一步证实了系统的正平衡点$ E_{\star} $呈稳态; 图 4表示的是时滞$ \tau = 2\in(\tau_{0}, +\infty) $的情况, 各子图分别为(a) $ t = 400, $ (b) $ t = 2000, $ (c) $ t = 4000, $ (d) $ t = 20000 $时刻下的空间斑图. 从图中可以看到螺旋波的产生过程, 这时正平衡点$ E_{\star} $是不稳定的. 图 3和图 4显示的结果不仅证实了定理3.4的理论结论, 而且展示了时滞的改变会影响种群密度的空间分布.

在本节的参数条件下, 我们可以根据不同的恐惧效应参数$ k_{0} $的取值, 求出不同的Hopf分支临界值$ \tau_{0} $, 从而对应不同的稳定区间$ [0, \tau_{0}) $. 通过计算得到表 1, 从中可得出如下结论: (1)在$ k_{0}<10 $时, 随着恐惧效应参数$ k_{0} $的增加, Hopf分支临界值$ \tau_{0} $不断减小, 而当$ k_{0}>10 $时, Hopf分支临界值$ \tau_{0} $不断增大; (2)随着恐惧效应参数$ k_{0} $的持续增大, 正平衡点$ E_{\star} $的值$ u_{\star}, v_{\star} $均不断减小.

图 5直观地展现了表格1中的结论(1), 即恐惧效应参数$ k_{0} $与Hopf分支的临界值$ \tau_{0} $呈非线性关系. 对于每一个取定的$ k_{0} $, 与其对应的$ \tau_{0} $决定了系统正平衡点$ E_{\star} $的稳定区间$ [0, \tau_{0}) $的大小, 故可看出随着恐惧效应的增大, 稳定区间先缩小后增大. 图 6中给出了$ t = 20000 $时不同恐惧效应参数值下捕食者的斑图, 其中子图(a) 中恐惧效应参数为$ k_{0} = 10 $, 子图(b)中$ k_{0} = 80 $, 其余参数全部与图 3 ($ k_{0} = 3 $)相同. 我们发现恐惧效应$ k_{0} = 10 $时的空间斑图 6(a)呈现为螺旋波, 恐惧效应$ k_{0} = 80 $时的空间斑图 6(b)则呈现为均匀状态. 结合图 3的结果, 我们得到在其他参数不变的情况下, 随着恐惧效应的增加, 系统将由稳定状态($ k_{0} = 3 $) 变为不稳定状态($ k_{0} = 10 $), 之后再变为稳定状态($ k_{0} = 80 $), 这进一步展示了恐惧效应与系统正平衡点的稳定区间的非线性关系的有趣现象.

但恐惧效应的不断增加是否真的对捕食-食饵系统的长期发展有利?接下来, 我们通过图 7展现了表格1的结论(2), 即恐惧效应$ k_{0} $与捕食者和食饵平衡状态下的值$ u_{\star}, v_{\star} $的关系, 随着恐惧效应$ k_{0} $的持续增加, 平衡点的值$ u_{\star}, v_{\star} $不断下降. 当恐惧效应足够大时, 平衡状态下捕食者$ v_{\star} $已经相当小. 而种群密度过低时可能会面临灭绝的危险(如不可预测的自然灾害, 传染病流行等), 这毋庸置疑对种群的生存是十分不利的!

基于恐惧效应对食饵的出生率具有重要的影响, 本文通过理论推导和数值模拟两方面研究了恐惧效应对具有时滞的反应扩散食饵-捕食者系统的Hopf分支和稳定区间的影响. 结果表明:

(1) 不论有无时滞, 恐惧效应均会影响系统正平衡点$ E_{\star} $的稳定性;

(2) 以时滞为参数时, 恐惧效应将影响系统的Hopf分支临界值$ \tau_{0} $, 在一定参数条件下, 当恐惧效应持续增加时, 稳定区间$ [0, \tau_{0}) $先缩小后增大, 恐惧效应与稳定区间的大小呈现非线性关系;

(3) 数值模拟结果显示, 系统具有复杂的动力学行为, 如出现了极限环, 螺旋波等现象.