随着现代工业和农业的迅速发展,大量的毒物和污染物进入了全球生态系统.因此,环境污染是最重要的一个社会生态问题.环境中各种有毒物质的存在对无保护种群的生存有着极大的威胁,有毒物质的无控制输入将会导致种群走向灭绝,所以环境污染下种群的持续生存问题,成为了生物数学研究的又一个热门课题.目前,已有许多学者对此做了大量的研究,得到了许多成果^[1-6].同时,我们注意到,在现实环境中,存在许多随机或偶然的因素影响着生物种群的变化^[7-11],比如:地震,洪水,海啸等会对种群产生瞬时的影响,而环境噪声会不断在不同程度上影响着种群的增长率、环境容量、竞争系数和系统的其他参数.因此,研究随机环境下具有脉冲效应的种群动力学行为,得到广大学者的青睐.

张树文等人^[12]研究了一个具有脉冲输入毒素的单种群随机干扰模型,并假设种群个体死亡后变成毒素^[13],得到了边界周期解以概率1全局吸引的充分条件,证明了系统均值有界且获得了种群灭绝与平均持续生存的阈值.而实际生态系统中,几乎没有独立的单种群存在,而是多种群共同存在的,考虑捕食-食饵关系的两种群生态模型更符合实际,因此,本文在种群个体死亡后变成毒素的假设下,考虑如下具有毒素脉冲输入和干扰的随机捕食-食饵模型

(1.1) $ \begin{equation} \left\{ \begin{array}{lll} \left. \begin{array}{lll} {\rm d}x(t)=x(t)(r-k_{1}c(t)-ax(t)-\rho_{1}y(t)){\rm d}t+\sigma_{1}x(t){\rm d}B_{1}(t), \\ {\rm d}y(t)=y(t)(-d-k_{2}c(t)+\rho_{2}x(t)){\rm d}t+\sigma_{2}y(t){\rm d}B_{2}(t), \\ {\rm d}c(t)=(-hc(t)+b_{1}k_{1}c(t)x(t)+b_{2}k_{2}c(t)y(t)){\rm d}t, \\ \end{array} \right\}t\neq nT, \\ \left. \begin{array}{lll} \Delta x(t)=0, \\ \Delta y(t)=0, \\ \Delta c(t)=P, \end{array} \right\}t=nT, \end{array}\right. \end{equation} $

其中${x}(t), {y}(t), c(t)$分别代表$t$时刻食饵,捕食者,种群所处环境的毒素密度. $r, a$分别表示食饵种群$x(t)$的内禀增长率和密度制约系数, $d$表示捕食者$y(t)$的死亡率, $h$表示毒素$c(t)$的降解率, $\rho_{1}$是捕食率, $\rho_{2}\ (\rho_{2}<\rho_{1})$是食饵转化率, $k_{1}, k_{2}$代表种群因毒素死亡率, $b_{1}, b_{2}$代表毒素转化率, $\Delta q(t)=q(t^+)-q(t), (q=x, y, c)$, $T$为脉冲周期, $P$代表脉冲时刻输入的毒素浓度的值, $\sigma_{i}\ (i=1, 2)$是白噪声强度, $B_{1}(t), B_{2}(t)$是定义在完备概率空间上相互独立的标准布朗运动,并假设以上所有参数均为正数.

在本文中,我们总假设$(\Omega, {\cal F}, \{{\cal F}_{t}\}_{t\geq0}, P)$为完备的概率空间,其中$\{{\cal F}_{t}\}_{t\geq0}$是$\Omega$上的一个$\sigma$ -代数且满足通常条件(即右连续, ${\cal F}_{0}$包含所有的零测集).并记

${\mathbb{R}}_{+}^{3}=\{(x_{1}, x_{2}, x_{3})\in {\mathbb{R}}^{3}:x_{i}>0, i=1, 2, 3\}.$

给出本文所需要的定义、引理、定理和一些记号如下:

$\langle f(t)\rangle_{t}=\lim\limits_{t\rightarrow+\infty}\frac{1}{t}\int_{0}^{t}f(s){\rm d}s, \quad\langle f(t)\rangle^{*}=\limsup\limits_{t\rightarrow+\infty}\frac{1}{t}\int_{0}^{t}f(s){\rm d}s, \quad\langle f(t)\rangle_{*}=\liminf\limits_{{t\rightarrow+\infty}}\frac{1}{t}\int_{0}^{t}f(s){\rm d}s.$

设$x(t)=(x_{1}(t), x_{2}(t), \cdots, x_{n}(t))(t\geq0)$是随机微分方程

${\rm d}x(t)=f(x(t), t){\rm d}t+g(x(t), t){\rm d}B(t)$

的解,其中$f\in {\mathit{£}}^{1}({\mathbb{R}}^{n}\times {\mathbb{R}}_{+}, {\mathbb{R}}^{n})$, $g\in{\mathit{£}}^{2}({\mathbb{R}}^{n}\times {\mathbb{R}}_{+}, {\mathbb{R}}^{n\times m})$, $B(t)$是$m$维布朗运动.

定义2.1^[14]   (1)若种群$x(t)$满足$\lim\limits_{t\rightarrow+\infty}x(t)=0$ a.s.,则称种群$x(t)$是灭绝的;

(2)若种群$x(t)$满足$\limsup\limits_{t\rightarrow+\infty}\frac{1}{t}\int_{0}^{t}x(s){\rm d}s=0$ a.s.,则称种群$x(t)$是非平均持续生存的;

(3)若存在常数$M>0$,使得种群$x(t)$满足$\liminf\limits_{t\rightarrow+\infty}\frac{1}{t}\int_{0}^{t}x(s){\rm d}s\geq M$ a.s.,则称种群$x(t)$是强平均持续生存的;

(4)若存在常数$N>0$,使得种群$x(t)$满足$\limsup\limits_{t\rightarrow+\infty}\frac{1}{t}\int_{0}^{t}x(s){\rm d}s\leq N<+\infty$ a.s.,则称种群$x(t)$是弱平均持续生存的.

定义2.2^[14]  设$X^{*}(t)=(x^{*}(t), y^{*}(t), c^{*}(t))$是系统(1.1)满足初始条件$x^{*}(0)>0, $$ y^{*}(0)>0, $$ c^{*}(0)>0$的解,对系统(1.1)满足初始条件$x(0)>0, y(0)>0, c(0)>0$的任意解$X(t)=(x(t), y(t), c(t))$,如果

$\lim\limits_{t\rightarrow+\infty}|x^{*}(t)-x(t)|=0, \quad \lim\limits_{t\rightarrow+\infty}|y^{*}(t)-y(t)|=0, \quad \lim\limits_{t\rightarrow+\infty}|c^{*}(t)-c(t)|=0$

均成立,则称$X^{*}(t)=(x^{*}(t), y^{*}(t), c^{*}(t))$是全局吸引的.

考虑如下非自治随机单种群Logstic模型

(2.1)$\begin{equation} {\rm d}N(t)=N(t)[(a(t)-b(t)N(t)){\rm d}t+\alpha(t){\rm d}B(t)], \ \ t\geq0, \end{equation} $

其中$a(t), b(t), \alpha(t) $在$[0, +\infty)$上是连续有界函数,且$a(t)>0, ~b(t)>0$.

引理2.1^[15]  系统(2.1)存在以$N(0)=N_{0}>0$为初始值的唯一全局连续正解

$N(t)=\frac{\exp\{\int_{0}^{t}[a(s)-\frac{\alpha^{2}(s)}{2}]{\rm d}s+\int_{0}^{t}\alpha(s){\rm d}B(s)\}}{\frac{1}{N_{0}}+\int_{0}^{t}b(s)\exp\{\int_{0}^{s}[a(\tau)-\frac{\alpha_{\tau}^{2}}{2}]{\rm d}\tau+\int_{0}^{s}\alpha(\tau){\rm d}B (\tau)\}}{\rm d}s, \ \ t\geq0.$

引理2.2^[16]  系统

$\left\{ \begin{array}{lll} {\rm d}W(t)=-aW(t){\rm d}t, \quad t\neq nT, \\ \Delta W(t)=P, \quad t= nT\end{array}\right.$

存在一个正周期解

$W^{*}(t)=\frac{ P\exp[-a(t-nT)]}{1-\exp(-aT)}, \quad t\in(nT, (n+1)T], ~n\in Z_{+}, $

其中$W^{*}(0^{+})=\frac{ P}{1-\exp(-aT)}$,且对此系统的所有的解$W(t)$均有当$t\rightarrow+\infty$时, $W(t)\rightarrow W^{*}(t)$成立.

引理2.3^[11]  设$x(t)\in C(\Omega\times[0, +\infty), {\mathbb{R}}^+)$, $B_{i}(t)\ (1\leq i\leq n)$为标准的布朗运动,则有

(1)若存在正数$T, \mu_{0}$使得

$\ln x(t)\leq\mu t-\mu_{0}\int_{0}^{t}x(s){\rm d}s+\sum\limits_{i=1}^{n}\beta_{i}B_{i}(t), \ t\geq T$

成立,其中$\beta_{i}\ (1\leq i\leq n)$是常数,那么

$\left\{ \begin{array}{l} \langle x(t)\rangle^{*}\leq\mu/\mu_{0}\quad & a.s.,\quad \mu>0,\\ \lim\limits_{t\rightarrow+\infty}x(t)=0\quad & a.s.,\quad \mu<0.\end{array} \right.$

(2)若存在正数$T, \mu_{0}$和$\mu>0$使得

$\ln x(t)\geq\mu t-\mu_{0}\int_{0}^{t}x(s){\rm d}s+\sum\limits_{i=1}^{n}\beta_{i}B_{i}(t), \ t\geq T$

成立,其中$\beta_{i}\ (1\leq i\leq n)$是常数,那么

$\langle x(t)\rangle_{*}\geq\mu/\mu_{0}\quad a.s..$

定理2.1^[14](伊藤公式)  设$x(t)(t\geq0)$是伊藤过程,其随机微分为

${\rm d}x(t)=f(t){\rm d}t+g(t){\rm d}B(t), $

其中$f\in{\mathit{£}}^{1}({\mathbb{R}}_{+}, {\mathbb{R}}^{n})$, $g\in{\mathit{£}}^{2}({\mathbb{R}}_{+}, {\mathbb{R}}^{n\times m})$.若$V(x(t), t)\in C^{2, 1}({\mathbb{R}}^{n}\times{\mathbb{R}}_{+};{\mathbb{R}})$,则$V(x(t), t)$仍然是伊藤过程,具有如下随机微分

${\rm d}V(x(t), t)=V_{t}(x(t), t){\rm d}t+V_{x}(x(t), t){\rm d}x(t)+\frac{1}{2}{\rm d}x^{\intercal}(t)V_{xx}(x(t), t){\rm d}x(t).$

定理3.1  对于任意给定的初值$(x(0^{+}), y(0^{+}), c(0^{+})) \in{\mathbb{R}}_{+}^{3}$,系统$(1.1)$存在唯一全局正解$(x(t), y(t), c(t))$, $t\geq0$,且此解以概率1停留在${\mathbb{R}}_{+}^{3}$中.

证  由系统的第一个方程知

${\rm d}x(t)\leq x(t)(r-ax(t)){\rm d}t+\sigma_{1}x(t){\rm d}B_{1}(t), $

构造比较系统

(3.1)$\begin{equation}\left\{ \begin{array}{lll} {\rm d}u_{1}(t)=u_{1}(t)(r-au_{1}(t)){\rm d}t+\sigma_{1}u_{1}(t){\rm d}B_{1}(t), \\ u_{1}(0^{+})=x(0^{+})>0, \end{array}\right.\end{equation}$

由引理2.1知系统(3.1)存在唯一全局连续正解如下

$u_{1}^{*}(t)=\frac{\exp\{(r-\frac{1}{2}\sigma_{1}^{2})t+\sigma_{1}B_{1}(t)\}}{\frac{1}{u_{1}(0^{+})}+a\int_{0}^{t}\exp\{(r-\frac{1}{2}\sigma_{1}^{2})s+\sigma_{1}B_{1}(s)\}{\rm d}s}, $

再由随机微分方程的比较定理得

(3.2)$\begin{equation}0\leq x(t)\leq u_{1}^{*}(t).\end{equation}$

由系统(1.1)的第二个方程和(3.2)式知

$\begin{array}{lll}{\rm d}y(t)&\leq y(t)(-d+\rho_{2}x(t)){\rm d}t+\sigma_{2}y(t){\rm d}B_{2}(t)\\&\leq y(t)(-d+\rho_{2}u_{1}^{*}(t)){\rm d}t+\sigma_{2}y(t){\rm d}B_{2}(t), \end{array}$

构造比较系统

(3.3)$\begin{equation}\left\{ \begin{array}{lll} {\rm d}u_{2}(t)=u_{2}(t)(-d+\rho_{2}u_{1}^{*}(t)){\rm d}t+\sigma_{2}u_{2}(t){\rm d}B_{2}(t), \\ u_{2}(0^{+})=y(0^{+})>0, \end{array}\right.\end{equation}$

易知系统(3.3)存在唯一的全局连续正解

$u_{2}^{*}(t)=u_{2}(0^{+})\exp\bigg\{\int_{0}^{t}[-d+\rho_{2}u_{1}^{*}(s)]{\rm d}s+\sigma_{2}B_{2}(t)\bigg\}, $

根据随机微分方程的比较定理可知

(3.4)$\begin{equation}0\leq y(t)\leq u_{2}^{*}(t).\end{equation}$

最后,由系统(1.1)的第三个方程, (3.2)和(3.4)式知

${\rm d}c(t)\leq c(t)(-h+b_{1}k_{1}u_{1}^{*}(t)+b_{2}k_{2}u_{2}^{*}(t)){\rm d}t, $

构造如下比较系统

(3.5)$\begin{equation}\left\{ \begin{array}{lll} {\rm d}u_{3}(t)=u_{3}(t)(-h+b_{1}k_{1}u_{1}^{*}(t)+b_{2}k_{2}u_{2}^{*}(t)){\rm d}t, \quad t\neq nT, \\ \Delta u_{3}(t)=P, \quad t= nT, \\ u_{3}(0^{+})=c(0^{+})>0, \end{array}\right.\end{equation}$

根据引理2.2,可得上述脉冲微分方程的全局渐近稳定周期解为

$u_{3}^{*}(t)=\frac{P\exp\{\int_{nT}^{t}[-h+b_{1}k_{1}u_{1}^{*}(s)+b_{2}k_{2}u_{2}^{*}(s)]{\rm d}s\}}{1-\exp\{\int_{nT}^{(n+1)T}[-h+b_{1}k_{1}u_{1}^{*}(s)+b_{2}k_{2}u_{2}^{*}(s)]{\rm d}s\}}, $

由脉冲微分方程的比较定理得

(3.6)$\begin{equation}0\leq c(t)\leq u_{3}^{*}(t).\end{equation}$

因此,由(3.2), (3.4)和(3.6)式可知,对任意的初始值$(x(0^{+}), y(0^{+}), c(0^{+}))$,当$t\geq0$时,系统(1.1)存在唯一的全局正解$(x(t), y(t), c(t))$.

定理3.2  如果$X(t)=(x(t), y(t), c(t))$是系统$(1.1)$满足初值$(x(0^{+}), y(0^{+}), c(0^{+}))\in{\mathbb{R}}_{+}^{3}$的解,则存在一个正数$M$,使得$Ex(t)\leq M, Ey(t)\leq M, Ez(t)\leq M$.

证  定义$V(x(t), y(t), c(t))=x(t)+y(t)+c(t)$,对$V(X)$ (其中$X=(x(t), y(t), c(t))$)运用伊藤公式有

(3.7)$\begin{eqnarray}{\rm d}V(X)&=&{\rm d}x(t)+{\rm d}y(t)+{\rm d}c(t) \\&=&(rx(t)-k_{1}c(t)x(t)-ax^{2}(t)-\rho_{1}x(t)y(t)-dy(t)-k_{2}c(t)y(t)+\rho_{2}x(t)y(t) \\&&-hc(t)+b_{1}k_{1}c(t)x(t)+b_{2}k_{2}c(t)y(t)){\rm d}t+\sigma_{1}x(t){\rm d}B_{1}(t)+\sigma_{2}y(t){\rm d}B_{2}(t) \\&\leq& (rx(t)-ax^{2}(t)-dy(t)-hc(t)){\rm d}t+\sigma_{1}x(t){\rm d}B_{1}(t)+\sigma_{2}y(t){\rm d}B_{2}(t), \end{eqnarray}$

当$t=nT$,

(3.8)$\begin{equation}\Delta V(t)=P, \end{equation}$

$\forall \ t\in(nT, (n+1)T]$,将(3.7)式从$nT$到t积分得

$\begin{eqnarray*}V(t)&\leq &V(nT^{+})+\int_{nT}^{t}[rx(s)-ax^{2}(s)-dy(s)-hc(s)]{\rm d}s\\&&+\int_{nT}^{t}\sigma_{1}x(s){\rm d}B_{1}(s)+\int_{nT}^{t}\sigma_{2}y(s){\rm d}B_{2}(s).\end{eqnarray*}$

对上式两边取期望有

$EV(t)\leq EV(nT^{+})+\int_{nT}^{t}E[rx(s)-ax^{2}(s)-dy(s)-hc(s)]{\rm d}s, $

易知

$\begin{eqnarray*}\frac{{\rm d}EV(t)}{{\rm d}t}&\leq& rEx(t)-aEx^{2}(t)-dEy(t)-hEc(t)\\&\leq &-a(Ex(t))^{2}+2rEx(t)-rEx(t)-dEy(t)-hEc(t), \end{eqnarray*}$

由于$-a(Ex(t))^{2}+2rEx(t)$的最大值为$\frac{r^{2}}{a}$,取$m=\min\{r, d, h\}$,有

(3.9)$\begin{equation} \frac{{\rm d}EV(t)}{{\rm d}t}\leq \frac{r^{2}}{a}-mEV(t), \end{equation}$

另一方面,由(3.8)式可知

(3.10)$\begin{equation}\Delta EV(t)=P, \end{equation}$

结合(3.9)和(3.10)式可知,若考虑如下脉冲微分方程

(3.11)$\begin{equation}\left\{ \begin{array}{lll} {\rm d}Q(t)=\Big[\frac{r^{2}}{a}-mQ(t)\Big]{\rm d}t, & \quad t\neq nT, \\[2mm] \Delta Q(t)=P, & \quad t=nT, \end{array}\right.\end{equation}$

则系统(3.11)的解为

$Q(t)=\frac{r^{2}}{am}+(Q(0^{+})-\frac{r^{2}}{am})\exp[-m(t-nT)], \quad t\in(nT, (n+1)T], $

其中

$Q(0^{+})=\frac{r^{2}}{am}+\frac{P}{1-\exp(-mT)}, $

由脉冲微分方程比较定理得

$\lim\limits_{t\rightarrow\infty}EV(t)\leq\lim\limits_{t\rightarrow\infty}EQ(t)\leq \frac{r^{2}}{am}+\frac{P}{1-\exp(-mT)}, $

即系统是均值有界的.

为进一步研究系统,令$x(t)=0, y(t)=0$,考虑如下系统

(3.12)$\begin{equation}\left\{ \begin{array}{lll} {\rm d}c(t)=-hc(t){\rm d}t, &\quad t\neq nT, \\ \Delta c(t)=P, &\quad t= nT, \end{array}\right.\end{equation}$

由引理2.2知系统(3.12)存在唯一正周期解$c^{\intercal}(t)=\frac{ P\exp[-h(t-nT)]}{1-\exp(-hT)}, ~t\in(nT, (n+1)T], n\in Z_{+}$,且当$t\rightarrow+\infty$时,有$c(t)\rightarrow c^{\intercal}(t)$.

显然,在任意脉冲区间$(nT, (n+1)T], n\in Z_{+}$上,系统存在一个边界周期解$(0, 0, c^{\intercal}(t))$.证毕.

定理3.3  如果系统(1.1)满足

$2r>\sigma_{1}^{2} \ \mbox{且} \ P>\frac{hT(r-\frac{\sigma_{1}^{2}}{2})}{k_{1}}, $

则系统(1.1)的边界周期解$(0, 0, c^{\intercal}(t))$以概率1全局吸引.

证  由系统(1.1)的第三个式子可知

$\left\{ \begin{array}{lll} {\rm d}c(t)\geq-hc(t){\rm d}t, & \quad t\neq nT, \\ \Delta c(t)=P, &\quad t= nT, \end{array}\right.$

因此考虑如下脉冲微分方程

(3.13)$\begin{equation}\left\{ \begin{array}{lll} {\rm d}v(t)=-hv(t){\rm d}t, \quad & t\neq nT, \\ \Delta v(t)=P, \quad & t= nT, \\ v(0^{+})=c(0^{+})>0, \end{array}\right.\end{equation}$

由脉冲微分方程比较定理得$c(t)\geq v(t)$,又由系统(3.12)可知,当$t\rightarrow+\infty$时, $v(t)\rightarrow c^{\intercal}(t)$.则当$t$充分大时,对任意的$\varepsilon_{1}>0$,有$c(t)>c^{\intercal}(t)-\varepsilon_{1}>0$.

即,对任意小的$\varepsilon_{1}>0$,存在$T_{1}>0$,使得当$t>T_{1}$时,有$c(t)>c^{\intercal}(t)-\varepsilon_{1}>0$,使得不等式

${\rm d}x(t)\leq x(t)(r-k_{1}(c^{\intercal}(t)-\varepsilon_{1})-ax(t)-\rho_{1}y(t)){\rm d}t+\sigma_{1}x(t){\rm d}B_{1}(t)$

成立.

构造比较系统

(3.14)$\begin{equation}\left\{ \begin{array}{lll} {\rm d}w(t)=w(t)(r-k_{1}(c^{\intercal}(t)-\varepsilon_{1})-aw(t)-\rho_{1}y(t)){\rm d}t+\sigma_{1}w(t){\rm d}B_{1}(t), \\ w(0)=x(0)>0, \end{array}\right.\end{equation}$

利用伊藤公式,沿着系统(3.14)对$\ln w(t)$求随机微分可得

$ \begin{eqnarray*}{\rm d}\ln w(t)&=&\frac{1}{w(t)}{\rm d}w(t)-\frac{1}{2w^{2}(t)}({\rm d}w(t))^{2}\\ &=&\Big[r-\frac{\sigma_{1}^{2}}{2}-k_{1}(c^{\intercal}(t)-\varepsilon_{1})-aw(t)-\rho_{1}y(t)\Big]{\rm d}t+\sigma_{1}{\rm d}B_{1}(t), \end{eqnarray*}$

两边从$0$到$t$积分并除以$t$得

(3.15)$\begin{eqnarray}\frac{\ln w(t)}{t}&=&\frac{\ln w(0)}{t}+\frac{1}{t}\int_{0}^{t}\Big[r-\frac{\sigma_{1}^{2}}{2}-k_{1}(c^{\intercal}(s)-\varepsilon_{1})-aw(s)-\rho_{1}y(s)\Big]{\rm d}s+\frac{\sigma_{1}B_{1}(t)}{t} \\&\leq& \frac{\ln w(0)}{t}+\frac{1}{t}\int_{0}^{t}\Big[r-\frac{\sigma_{1}^{2}}{2}-k_{1}(c^{\intercal}(s)-\varepsilon_{1})\Big]{\rm d}s+\frac{\sigma_{1}B_{1}(t)}{t} \\&=&\frac{\ln w(0)}{t}+\Big(r-\frac{\sigma_{1}^{2}}{2}\Big)-\frac{1}{t}\int_{0}^{t}[k_{1}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s+\frac{\sigma_{1}B_{1}(t)}{t}, \end{eqnarray}$

由于$c^{\intercal}(s)$是以$T$为周期的周期解,则有

$\frac{1}{t}\int_{0}^{t}[k_{1}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s\geq \frac{1}{(n+1)T}\int_{0}^{nT}[k_{1}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s, $

$\frac{1}{t}\int_{0}^{t}[k_{1}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s\leq \frac{1}{nT}\int_{0}^{(n+1)T}[k_{1}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s, $

当$t\rightarrow+\infty$时,即当$n\rightarrow+\infty$时有

$\frac{1}{(n+1)T}\int_{0}^{nT}[k_{1}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s \rightarrow \frac{1}{T}\int_{0}^{T}[k_{1}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s, $

$\frac{1}{nT}\int_{0}^{(n+1)T}[k_{1}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s \rightarrow \frac{1}{T}\int_{0}^{T}[k_{1}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s, $

因此当$t\rightarrow+\infty$时,有

(3.16)$\begin{equation} \frac{1}{t}\int_{0}^{t}[k_{1}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s \rightarrow \frac{1}{T}\int_{0}^{T}[k_{1}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s = \frac{Pk_1}{hT}.\end{equation}$

由于$\varepsilon_{1}$足够小,则当$t\rightarrow+\infty$时,由大数定律并结合(3.16)式,对(3.15)式取上确界极限有

$\limsup\limits_{t\rightarrow+\infty}\frac{1}{t}\ln w(t)\leq r-\frac{\sigma_{1}^{2}}{2}-\frac{1}{T}\int_{0}^{T}k_{1}c^{\intercal}(s){\rm d}s=r-\frac{Pk_{1}}{hT}-\frac{\sigma_{1}^{2}}{2}<0, $

所以当$t\rightarrow+\infty$时, $w(t)\rightarrow0$,即当$t\rightarrow+\infty$时, $x(t)\rightarrow0$.

则对任意的$0<\varepsilon_{2}< \frac{1}{\rho_{2}}(d+\frac{\sigma_{2}^{2}}{2}+\frac{k_{2}P}{hT})$,存在$T_{2}>T_{1}>0$,当$t>T_{2}$时,有$0<x(t)<\varepsilon_{2}$.因此有不等式

$\begin{array}{lll}{\rm d}y(t)\leq y(t)[-d-k_{2}(c^{\intercal}(t)-\varepsilon_{1})+\rho_{2}\varepsilon_{2}]{\rm d}t+\sigma_{2}y(t){\rm d}B_{2}(t).\\\end{array}$

构造比较系统

(3.17)$\begin{equation}\left\{ \begin{array}{lll} {\rm d}z(t)=z(t)[-d-k_{2}(c^{\intercal}(t)-\varepsilon_{1})+\rho_{2}\varepsilon_{2}]{\rm d}t+\sigma_{2}z(t){\rm d}B_{2}(t), \\ z(0)=y(0)>0, \end{array}\right.\end{equation}$

利用伊藤公式,沿着系统(3.17)对$\ln z(t)$求随机微分可得

$ \begin{eqnarray*}{\rm d}\ln z(t)&=&\frac{1}{z(t)}{\rm d}z(t)-\frac{1}{2z^{2}(t)}({\rm d}z(t))^{2}\\ &=&\Big[-d-\frac{\sigma_{2}^{2}}{2}-k_{2}(c^{\intercal}(t)-\varepsilon_{1})+\rho_{2}\varepsilon_{2}\Big]{\rm d}t+\sigma_{2}{\rm d}B_{2}(t), \end{eqnarray*} $

两边从$0$到$t$积分并除以$t$得

(3.18)$\begin{equation}\frac{\ln z(t)}{t}=\frac{\ln z(0)}{t}+\Big(-d-\frac{\sigma_{2}^{2}}{2}+\rho_{2}\varepsilon_{2}\Big)-\frac{1}{t}\int_{0}^{t}[k_{2}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s+\frac{\sigma_{1}B_{2}(t)}{t}, \end{equation}$

由于$c^{\intercal}(s)$是以$T$为周期的周期解,同理,当$t\rightarrow+\infty$时,有

(3.19)$\begin{equation} \frac{1}{t}\int_{0}^{t}[k_{2}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s \rightarrow \frac{1}{T}\int_{0}^{T}[k_{2}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s = \frac{k_{2}P}{hT}.\end{equation}$

由于$\varepsilon_{1}$足够小,则当$t\rightarrow+\infty$时,由大数定律并结合(3.19)式,对(3.18)式取上确界极限有

$\limsup\limits_{t\rightarrow+\infty}\frac{1}{t}\ln z(t)= -d-\frac{\sigma_{2}^{2}}{2}+\rho_{2}\varepsilon_{2}-\frac{k_{2}P}{hT}.$

又因为$0<\varepsilon_{2}< \frac{1}{\rho_{2}}(d+\frac{\sigma_{2}^{2}}{2}+\frac{k_{2}P}{hT})$,则$-d-\frac{\sigma_{2}^{2}}{2}+\rho_{2}\varepsilon_{2}-\frac{k_{2}P}{hT}<0$,所以有$\limsup\limits_{t\rightarrow+\infty}\frac{1}{t}\ln z(t)<0$.即当$t\rightarrow+\infty$时,有$z(t)\rightarrow0$.则当$t\rightarrow+\infty$时,有$y(t)\rightarrow0$.

同理可证,当$t\rightarrow+\infty$时,有$c(t)\rightarrow c^{\intercal}(t)$.

因此,当$2r>\sigma_{1}^{2}$且$P>\frac{hT(r-\frac{\sigma_{1}^{2}}{2})}{k_{1}}$成立时,系统(1.1)的边界周期解$(0, 0, c^{\intercal}(t))$是全局吸引的.

定理3.4  当$2r>\sigma_{1}^{2}$且$P=\frac{hT(r-\frac{\sigma_{1}^{2}}{2})}{k_{1}}$时,系统(1.1)的解满足

$\limsup\limits_{t\rightarrow+\infty}\frac{1}{t}\int_{0}^{t}x(s){\rm d}s\leq\frac{r-\frac{Pk_{1}}{hT}-\frac{\sigma_{1}^{2}}{2}}{a}=0 \ \ {\rm a.s., }$

则种群$x(t)$是非平均持续生存的,种群$y(t)$是灭绝的.

证  结合系统(3.15)和脉冲微分方程比较定理可得: $c(t)\geq c^{\intercal}(t)-\varepsilon_{1}$.记$\delta^{*}=r-\frac{Pk_{1}}{hT}-\frac{\sigma_{1}^{2}}{2}$,根据定理3.3的证明过程可知,对任意给定的$\varepsilon_{3}>0$,存在$T_{3}>T_{2}>0$,使得当$t>T_{3}$时有

$\frac{\ln w(0)}{t}\leq\frac{\varepsilon_{3}}{3}, \quad (r-\frac{\sigma_{1}^{2}}{2})-\frac{1}{t}\int_{0}^{t}[k_{1}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s\leq\delta^{*}+\frac{\varepsilon_{3}}{3}, \quad \frac{\sigma_{1}B_{1}(t)}{t}\leq\frac{\varepsilon_{3}}{3}, $

将上式代入

$\frac{\ln w(t)}{t}\leq\frac{\ln w(0)}{t}+(r-\frac{\sigma_{1}^{2}}{2})-\frac{1}{t}\int_{0}^{t}[k_{1}(c^{\intercal}(s)-\varepsilon_{1})]{\rm d}s-\frac{a}{t}\int_{0}^{t}w(s){\rm d}s+\frac{\sigma_{1}B_{1}(t)}{t}, $

有

(3.20)$\begin{equation}\ln w(t)\leq\lambda t-a\int_{0}^{t}w(s){\rm d}s, \quad t\geq T_{3}, \end{equation}$

其中$\lambda=\delta^{*}+\varepsilon_{3}$,因此当$t\geq T_{3}$时,根据比较定理有

(3.21)$\begin{equation}x(t)\leq w(t), \quad \ln w(t)\leq\lambda t-a\int_{0}^{t}x(s){\rm d}s, \end{equation}$

由(3.20)和(3.21)式可知

(3.22)$\begin{equation}\ln x(t)\leq\lambda t-a\int_{0}^{t}x(s){\rm d}s.\end{equation}$

令$h(t)=\int_{0}^{t}x(s){\rm d}s$,则(3.22)式可化为

$\frac{{\rm d}h(t)}{{\rm d}t}\exp(ah(t))\leq\exp(\lambda t), \quad t\geq T_{3}, $

上式两端同时从$T_{3}$到$t$积分得

$\frac{1}{a}[\exp(ah(t))-\exp(ah(T_{3}))]\leq\frac{1}{\lambda}[\exp(\lambda t)-\exp(\lambda T_{3})], $

整理得

$\exp(ah(t))\leq\exp(ah(T_{3}))+\frac{a}{\lambda}[\exp(\lambda t)-\exp(\lambda T_{3})], $

上式两端同时取对数有

$h(t)\leq\frac{1}{a}\ln\Big\{\exp(ah(T_{3}))+\frac{a}{\lambda}[\exp(\lambda t)-\exp(\lambda T_{3})] \Big\}, $

两端同时除以$t$,并取上确界极限有

$\limsup\limits_{t\rightarrow+\infty}\frac{1}{t}\int_{0}^{t}x(s){\rm d}s\leq\frac{1}{a}\limsup\limits_{t\rightarrow+\infty}\frac{1}{t}\ln\Big\{ \exp(ah(T_{3}))+\frac{a}{\lambda}[\exp(\lambda t)-\exp(\lambda T_{3})] \Big\}, $

由洛必达法则可得

$\limsup\limits_{t\rightarrow+\infty}\frac{1}{t}\int_{0}^{t}x(s){\rm d}s\leq\frac{\lambda}{a}, $

其中$\lambda=\delta^{*}+\varepsilon_{3}$,由$\varepsilon_{3}$的任意性可知

$\limsup\limits_{t\rightarrow+\infty}\frac{1}{t}\int_{0}^{t}x(s){\rm d}s\leq\frac{\delta^{*}}{a}=\frac{r-\frac{Pk_{1}}{hT}-\frac{\sigma_{1}^{2}}{2}}{a} \ \ {\rm a.s., }$

则当$r=\frac{Pk_{1}}{hT}+\frac{\sigma_{1}^{2}}{2}$时, $\limsup\limits_{t\rightarrow+\infty}\frac{1}{t}\int_{0}^{t}x(s){\rm d}s=0$ a.s.,即种群$x(t)$是非平均持续生存的.

同理,沿着系统的解对$\ln y(t)$求随机微分可得

$ \begin{eqnarray*}{\rm d}\ln y(t)&=&[-d-k_{2}c(t)+\rho_{2}x(t)]{\rm d}t+\sigma_{2}{\rm d}B_{2}(t)\\&\leq&[-d-k_{2}(c^{\intercal}(t)-\varepsilon_{1})+\rho_{2}x(t)]{\rm d}t+\sigma_{2}{\rm d}B_{2}(t) \end{eqnarray*}$

两端同时从$0$到$t$取积分并除以$t$得

$\frac{\ln y(t)}{t}\leq\frac{\ln y(0)}{t}+\frac{1}{t}\int_{0}^{t}\Big[-d-\frac{\sigma_{2}^{2}}{2}-k_{2}(c^{\intercal}(t)-\varepsilon_{1})+\rho_{2}x(t)\Big]{\rm d}t+\frac{\sigma_{2}B_{2}(t)}{t}, $

由$\varepsilon_{1}$的任意性及大数定律,对上式两端取极限有

$\limsup\limits_{t\rightarrow+\infty}\frac{1}{t}\ln y(t)\leq-d-\frac{\sigma_{2}^{2}}{2}-\frac{1}{t}\int_{0}^{t}k_{2}(c^{\intercal}(s)){\rm d}s<0, $

则$\lim\limits_{t\rightarrow+\infty}y(t)=0$,即种群$y(t)$是灭绝的.

定理3.5  当$2r>\sigma_{1}^{2}$且$P<\frac{[\rho_{2}(r-\frac{\sigma_{1}^{2}}{2})-a(d+\frac{\sigma_{2}^{2}}{2})]hT}{ak_{2}+\rho_{2}k_{1}} < \frac{hT(r-\frac{\sigma_{1}^{2}}{2})}{k_{1}}$时,有

$\limsup\limits_{t\rightarrow+\infty}\frac{1}{t}\int_{0}^{t}x(s){\rm d}s\leq\frac{r-\frac{Pk_{1}}{hT}-\frac{\sigma_{1}^{2}}{2}}{a} \ a.s., \ \ \\liminf\limits_{t\rightarrow+\infty}\frac{1}{t}\int_{0}^{t}y(s){\rm d}s\geq \frac{\lambda a}{\rho_{1}\rho_{2}} \ a.s., $

其中$\lambda=\frac{\rho_{2}(r-\frac{Pk_{1}}{hT}-\frac{\sigma_{1}^{2}}{2})}{a}-d-\frac{\sigma_{2}^{2}}{2}-\frac{Pk_{2}}{hT}>0$,即种群$x(t)$是弱平均持续生存的,种群$y(t)$是强平均持续生存的.

证  考虑系统

(3.23)$\begin{equation}{\rm d}\Psi(t)= \Psi(t)(r-k_{1}(c^{\intercal}(t)-\varepsilon_{1})-a\Psi(t)){\rm d}t+\sigma_{1}\Psi(t){\rm d}B_{1}(t), \Psi(0)=x(0)>0, \end{equation}$

根据文献[17]及任意小的$\varepsilon_{1}$,可得

(3.24)$\begin{equation}\left\{ \begin{array}{lll} \lim\limits_{t\rightarrow\infty}\Psi(t)=0 \ a.s., & \mbox{若} \ P>\frac{hT(r-\frac{\sigma_{1}^{2}}{2})}{k_{1}}, \\[3mm] \lim\limits_{t\rightarrow\infty}\frac{1}{t}\int_{0}^{t}\Psi(s){\rm d}s=\frac{r-\frac{Pk_{1}}{hT}-\frac{\sigma_{1}^{2}}{2}}{a} \ a.s., & \mbox{若} \ P<\frac{hT(r-\frac{\sigma_{1}^{2}}{2})}{k_{1}}, \end{array}\right.\end{equation}$

由微分方程的比较定理,可知$x(t)\leq \Psi(t), \ t\geq0 $,则有$\limsup\limits_{t\rightarrow+\infty}\frac{1}{t}\int_{0}^{t}x(s) {\rm d}s\leq\frac{r-\frac{Pk_{1}}{hT}-\frac{\sigma_{1}^{2}}{2}}{a}$ a.s.,即当$P<\frac{hT(r-\frac{\sigma_{1}^{2}}{2})}{k_{1}}$时,种群$x(t)$是弱平均持续生存的.

根据系统(1.1)的第一个等式和(3.23)式,对$\ln x(t), \ln\Psi (t)$求随机微分可得

(3.25)$\begin{equation}\left\{ \begin{array}{lll} \frac{\ln x(t)-\ln x(0)}{t}=r-\frac{\sigma_{1}^{2}}{2}-\frac{a}{t}\int_{0}^{t}x(s){\rm d}s-\frac{\rho_{1}}{t}\int_{0}^{t}y(s){\rm d}s-\frac{k_{1}}{t}\int_{0}^{t}c(s){\rm d}s+\frac{\sigma_{1}B_{1}(t)}{t}, \\\frac{\ln \Psi(t)-\ln \Psi(0)}{t}=r-\frac{\sigma_{1}^{2}}{2}-\frac{a}{t}\int_{0}^{t}\Psi(s){\rm d}s-\frac{k_{1}}{t}\int_{0}^{t}c(s){\rm d}s+\frac{\sigma_{1}B_{1}(t)}{t}, \end{array}\right. \end{equation}$

由(3.25)式可得

$0\geq\frac{\ln x(t)-\ln \Psi(t)}{t}=\frac{a}{t}\int_{0}^{t}(\Psi(s)-x(s)){\rm d}s-\frac{\rho_{1}}{t}\int_{0}^{t}y(s){\rm d}s, $

即

(3.26)$\begin{equation}\frac{1}{t}\int_{0}^{t}(\Psi(s)-x(s)){\rm d}s\leq\frac{\rho_{1}}{at}\int_{0}^{t}y(s){\rm d}s, \end{equation}$

结合(3.26)式对系统(1.1)的第二个式子运用伊藤公式可得

$ \begin{eqnarray*}\frac{\ln y(t)-\ln y(0)}{t}&=&-d-\frac{\sigma_{2}^{2}}{2}-\frac{k_{2}}{t}\int_{0}^{t}c(s){\rm d}s+\frac{\rho_{2}}{t}\int_{0}^{t}\Psi(s){\rm d}s\\&&-\frac{\rho_{2}}{t}\int_{0}^{t}(\Psi(s)-x(s)){\rm d}s+\frac{\sigma_{2}B_{2}(t)}{t}\\&\geq &-d-\frac{\sigma_{2}^{2}}{2}-\frac{k_{2}}{t}\int_{0}^{t}c(s){\rm d}s+\frac{\rho_{2}}{t}\int_{0}^{t}\Psi(s){\rm d}s\\&&-\frac{\rho_{1}\rho_{2}}{at}\int_{0}^{t}y(s){\rm d}s+\frac{\sigma_{2}B_{2}(t)}{t}, \end{eqnarray*} $

由(3.24)式可知,对任意的$\varepsilon>0$,存在$T>0$,使得当$t>T$时,有

$\frac{1}{t}\int_{0}^{t}\Psi(s){\rm d}s>\frac{r-\frac{Pk_{1}}{hT}-\frac{\sigma_{1}^{2}}{2}}{a}-\frac{\varepsilon}{\rho_{2}}.$

则对任意的$t>T$,有

$ \begin{eqnarray*}\frac{\ln y(t)-\ln y(0)}{t}&\geq&-d-\frac{\sigma_{2}^{2}}{2}-\frac{Pk_{2}}{hT}+\frac{\rho_{2}(r-\frac{Pk_{1}}{hT}-\frac{\sigma_{1}^{2}}{2})}{a}-\varepsilon-\frac{\rho_{1}\rho_{2}}{at}\int_{0}^{t}y(s){\rm d}s+\frac{\sigma_{2}B_{2}(t)}{t}\\&=&\lambda-\varepsilon-\frac{\rho_{1}\rho_{2}}{at}\int_{0}^{t}y(s){\rm d}s+\frac{\sigma_{2}B_{2}(t)}{t}, \end{eqnarray*}$

其中$\lambda=-d-\frac{\sigma_{2}^{2}}{2}-\frac{Pk_{2}}{hT}+\frac{\rho_{2}(r-\frac{Pk_{1}}{hT}-\frac{\sigma_{1}^{2}}{2})}{a}>0$,即

$\ln y(t)\geq (\lambda-\varepsilon)t-\frac{\rho_{1}\rho_{2}}{a}\int_{0}^{t}y(s){\rm d}s+\sigma_{2}B_{2}(t), $

由引理2.3和$\varepsilon$的任意性知

$\liminf\limits_{t\rightarrow+\infty}\frac{1}{t}\int_{0}^{t}y(s){\rm d}s\geq \frac{\lambda a}{\rho_{1}\rho_{2}} \ \ {\rm a.s..}$

证毕.

令周期数$T=10$,取参数$r=0.1, k_1=0.2, k_2=0.1, a=0.05, \rho_1=0.5, \rho_2=0.3, $$ d=0.3, h=0.5, b_1=0.1, b_2=0.1, \sigma_{1}=0.4, \sigma_{2}=0.3$.并取初始值为$(x_{1}(0), y_{1}(0), c_{1}(0))=(1, 0.3, 2)$和$(x_{2}(0), y_{2}(0), c_{2}(0))=(1, 0.3, 3)$.

(1)令脉冲值$P=1$,则满足条件$0.2=2r>\sigma_{1}^{2}=0.16$且$1=P>\frac{hT(r-\frac{\sigma_{1}^{2}}{2})}{k_{1}}=0.5$,由定理3.3知,系统(1.1)满足初始值下的解$(x_{1}(t), y_{1}(t), c_{1}(t)), (x_{2}(t), y_{2}(t), c_{2}(t))$全局吸引到$(0, 0, c^{\intercal}(t))$,见图 1.

(2)令脉冲值$P=5$,且保持初始值及其他参数不变,则满足条件$0.2=2r>\sigma_{1}^{2}=0.16$且$5=P>\frac{hT(r-\frac{\sigma_{1}^{2}}{2})} {k_{1}}=0.5$,由定理3.3知,系统(1.1)满足初始值下的解$(x_{1}(t), y_{1}(t), c_{1}(t)), $$ (x_{2}(t), y_{2}(t), c_{2}(t))$全局吸引到$(0, 0, c^{\intercal}(t))$,见图 2.

(3)令脉冲值$P=10$,且保持初始值及其他参数不变,则满足条件$0.2=2r>\sigma_{1}^{2}=0.16$且$10=P>\frac{hT(r-\frac{\sigma_{1} ^{2}}{2})}{k_{1}}=0.5$,由定理3.3知,系统(1.1)满足初始值下的解$(x_{1}(t), y_{1}(t), c_{1}(t)), $$(x_{2}(t), y_{2}(t), c_{2}(t))$全局吸引到$(0, 0, c^{\intercal}(t))$,见图 3.

由图 1, 图 2, 图 3可知,其他参数一定时,周期投放毒素脉冲值越大,种群灭绝的速度越快.

(4)令脉冲值$P=0.1$, $\sigma_{1}=0.2$,且保持初始值及其他参数不变,则满足条件$0.2=2r>\sigma_{1}^{2}=0.04$,且$0.1=P<\frac{[\rho_{2}(r-\frac{\sigma_{1}^{2}}{2})-a(d+\frac{\sigma_{2}^{2}}{2})]hT}{ak_{2}+\rho_{2}k_{1}}=0.9038<2.375=\frac{hT(r-\frac{\sigma_{1}^{2}}{2})}{k_{1}}$,由定理3.5知,种群$x(t)$是弱平均持续生存的,种群$y(t)$是强平均持续生存的,见图 4.

本文通过研究周期脉冲毒素输入的随机捕食-食饵系统,证明了系统的均值有界性,全局正解的存在唯一性,并得到种群灭绝与持续生存的阈值$P^{*}=\frac{hT(r-\frac{\sigma_{1}^{2}}{2})}{k_{1}}$,使得

(1)当$2r>\sigma_{1}^{2}$,且$P>P^{*}$时,系统(1.1)存在边界周期解$(0, 0, c^{\intercal}(t))$且以概率1全局吸引;

(2)当$2r>\sigma_{1}^{2}$,且$P=P^{*}$时,系统(1.1)中种群$x(t)$是非平均持续生存的,种群$y(t)$是灭绝的;

(3)当$2r>\sigma_{1}^{2}$,且$P<\frac{[\rho_{2}(r-\frac{\sigma_{1}^{2}}{2})-a(d+\frac{\sigma_{2}^{2}}{2})]hT}{ak_{2}+\rho_{2}k_{1}} < P^{*}$时,种群$x(t)$是弱平均持续生存的,种群$y(t)$是强平均持续生存的.

由理论结果和数值模拟可知,当脉冲毒素输入超过一定阈值$P^{*}$时,会导致种群灭绝,并且周期投放毒素脉冲值越大,种群将越快的趋于灭绝;而当脉冲输入毒素和白噪声足够小时,两种群都将持续生存.