Relaxation oscillation and canard explosion in a slow–fast predator–prey model with Beddington–DeAngelis functional response

Abstract

In this paper, we consider a predator–prey model with Beddington–DeAngelis functional response. Considering the predator’s rate of growth and death is much lower than that of prey’s, the model becomes a slow–fast system that mathematically leads to a singular perturbation problem. Using geometric singular perturbation theory due to Fenichel and blow up technique, we have investigated the system and obtained very rich and complicated dynamical phenomena including the existence of relaxation oscillation, canard cycles near the Hopf bifurcation point and the interesting phenomenon of canard explosion.

Appendix A

The First Lyapunov Coefficient Adopting the procedure discussed in [50], we now compute the first Lyapunov coefficient in order to study the stability of Hopf-bifurcating periodic solution. Translating the equilibrium point $(X_{2\epsilon },Y_{2\epsilon })$ of the system (49) into the origin by using the change of variables $u=X_2-X_{2\epsilon }, v=Y_2-Y_{2\epsilon }$ and expanding the reduced system by Taylor’s series about $(u, v) = (0, 0)$, we obtain

$$\begin{aligned} \frac{d\mathrm{u}}{\mathrm{d}t}&=\alpha _1 u+\alpha _2v+\alpha _3 u^2+\alpha _4v^2+\alpha _5 uv+\alpha _6 u^3, \end{aligned}$$

(69)

$$\begin{aligned} \frac{\mathrm{d}v}{\mathrm{d}t}&=\beta _1 u+\beta _2 v+\beta _3 u^2, \end{aligned}$$

(70)

where $\alpha _1=-a_5\sqrt{\epsilon }+{\mathcal {O}}(2)$, $\alpha _2=-b_1+{\mathcal {O}}(2)$, $\alpha _3=-b_2+{\mathcal {O}}(2)$, $\alpha _4=0$, $\alpha _5=-a_2\sqrt{\epsilon }$, $\alpha _6=a_3\sqrt{\epsilon }$, $\beta _1=b_3+{\mathcal {O}}(2)$, $\beta _2=a_5\sqrt{\epsilon }$, $\beta _3=a_4\sqrt{\epsilon }$.

Here, all the coefficients are evaluated at $\lambda =\lambda _2^*$ or $\left( \lambda =\lambda _H(\sqrt{\epsilon })=-\frac{a_5b_3}{2b_2x_*y_*}\epsilon +{\mathcal {O}}(\epsilon ^{3/2})\right) .$

Hence, it follows that at $\lambda =\lambda _H(\sqrt{\epsilon })$, $\alpha _1+\beta _2=0$ and also $\alpha '_1(\lambda _H)+\beta '_2(\lambda _H)=\frac{2b_2}{b_3}x_*y_*\ne 0$. We assume that $\alpha _1(\lambda _H)\beta _2(\lambda _H)-\alpha _2(\lambda _H)\beta _1(\lambda _H)=b_1b_3+{\mathcal {O}}(2)>0$.

Now following [50], the first Lyapunov coefficient $l_1$, for a general system of the form

$$\begin{aligned} {\dot{x}}&=ax+by+{\mathcal {O}}(2)=\rho (x,y) \end{aligned}$$

(71)

$$\begin{aligned} {\dot{y}}&=cx+dy+{\mathcal {O}}(2)=\zeta (x,y) \end{aligned}$$

(72)

where the above-mentioned conditions are satisfied is given by

$$\begin{aligned} \begin{aligned} l_1&=\frac{-3\pi }{2b{\omega _0}^\frac{3}{2}}\Big \{\Big [ac\left( \rho _{xy}^2+\rho _{xy}\zeta _{yy}+\rho _{yy}\zeta _{xy}\right) \\&\quad +ab\left( \zeta _{xy}^2+\rho _{xx}\zeta _{xy}+\rho _{xy}\zeta _{yy}\right) \\&\quad +\,c^2\left( \rho _{xy}\rho _{yy}+2\rho _{yy}\zeta _{yy}\right) -2ac\left( \zeta _{yy}^2-\rho _{xx}\rho _{yy}\right) \\&\quad -2ab\left( \rho _{xx}^2 -\zeta _{xx}\zeta _{yy}\right) -\,b^2\left( \zeta _{xy}\zeta _{xx}+2\rho _{xx}\zeta _{xx}\right) \\&\quad \left. +\left( bc-2a^2\right) \left( \zeta _{xy}\zeta _{yy}-\rho _{xy}\rho _{xx}\right) \right] \\&\quad -\left( a^2+bc\right) \left[ 3\left( c\zeta _{yyy}-b\rho _{xxx}\right) \right. \\&\quad +2a\left( \rho _{xxy}+\zeta _{xyy}\right) +\left( c\rho _{xyy}-b\zeta _{xxy}\right) \Big ]\Big \} \end{aligned} \end{aligned}$$

(73)

evaluated at the critical parameter value where the above conditions are satisfied, $\pm i\omega _0$ are being the eigenvalues of the Jacobian matrix. Hence, in our case following the above result, on simplification, the first Lyapunov coefficient $l_1$ is given by

$$\begin{aligned} l_1=\frac{3\pi }{(b_1b_3)^\frac{3}{2}}A\sqrt{\epsilon }+{\mathcal {O}}(\epsilon ) \end{aligned}$$

(74)

where $A=a_2b_2b_3+4a_4b_1b_2-4a_5b_2^2+9a_3b_1b_3$.