# Dynamics of a Predator-Prey System with a Mate-Finding Allee Effect

Authors: Ruiwen Wu; Xiuxiang Liu
Journal: Abstract and Applied Analysis
 (2014)

---

## Abstract

We consider a ratio-dependent predator-prey system with a mate-finding Allee effect on prey. The stability properties of the equilibria and a complete bifurcation analysis, including the existence of a saddle-node, a Hopf bifurcation, and, a Bogdanov-Takens bifurcations, have been proved theoretically and numerically. The blow-up method has been applied to investigate the structure of a neighborhood of the origin. Our mathematical results show the mate-finding Allee effect can reduce the complexity of system behaviors by making the complicated equilibrium less complicated, and it can be a destabilizing force as well, which makes the system has a high possibility of being threatened with extinction in ecology.


---

## Body

## 1. Introduction

Just as pointed out by Berryman in [1] that “this dynamical relationship between predators and their prey has long been and will continue to be one of the dominant themes in both ecology and mathematical ecology due to its universal existence and importance”, both ecologists and mathematicians are interested in the dynamical analysis of predator-prey models. It is well known that the classical Gause type predator-prey system is written by
(1)
d
N
d
τ
=
G
(
N
)
N
-
f
(
N
,
P
)
P
,
d
P
d
τ
=
e
f
(
N
,
P
)
-
m
P
,
where N
(
τ
) and P
(
τ
) are the densities of prey and predator, respectively, e represents the trophic efficiency or the conversion efficiency of predator, and the parameter m characterizes the predator natural mortality rate. The functional response f
(
N
,
P
), which quantifies the amount of prey consumed per predator per unit time and plays an important role in predator-prey dynamics, is conventionally modeled as prey-dependent, where the prey consumption rate by an average predator is only a function of prey density alone; that is, f
(
N
,
P
)
=
f
(
N
). Different prey-dependent response types (e.g., the mass-action approach in Lotka-Volterra model and Holling types I–III) have been used to model the predator-prey interactions and get success in describing some ecological communities. When the spatial structure of one or both of the interacting populations is involved, it would be more plausible to take the predator-dependent functional form, where both predator and prey densities affect the response, for example, Hassell and Varley function response [2] and Beddington-DeAngelis functional response functional [3]. In this paper, we are interested in one important form of predator-dependent functional response, called ratio-dependent response; that is, the functional response depends on the term N
/
P and the corresponding Monod-Holling hyperbolic form is
(2)
f
(
N
P
)
=
a
N
N
+
a
h
P
,
where a is the maximum prey consumption rate and h is the predator handling time, which is proposed by Arditi and Ginzburg [4] and studied by many authors; see, for example, Berezovskaya et al. [5, 6], Kuang et al. [7–10], and Zhang et al. [11, 12] in which we know this ratio-dependent functional response provides more reasonable explanations and accurate predictions when community-level situations of food chains and food webs are considered.Another functionG
(
N
) in (1) is the per capita growth rate of prey population in the absence of predators. The classical view of popular dynamics is that the major ecological force at work is the release from the constraints of intraspecific competition when a population is small or at low density. The fewer we are, the more we all have and the better welfare be there. In general, G
(
N
) takes the form
(3)
G
(
N
)
=
r
(
1
-
N
K
)
,
where r and K are the intrinsic growth rate and the environment carrying capacity of prey, respectively. However, the ecologist Allee gathered sufficient experimental and observational data to conclude that the evolution of social structures was not only driven by competition, but that cooperation was another, if not the most, fundamental principle in animal species [13]. The dynamical consequences of this importance of animal aggregations directly led to what Odum called in 1953 “the Allee principle,” now known as the Allee effect [14].There are two conventional types of Allee effects: component Allee effect and demographic Allee effect. The positive relationship between any component of individual fitness and population size can be regarded as a component Allee effect, while a demographic Allee effect is linked to the level of overall individual fitness [15]. Moreover, a demographic Allee effect is strong if there exists an Allee threshold, below which the per capita growth rate becomes negative. The classical and simplest types of demographic Allee effects are considered commonly, which takes the mathematical form of
(4)
G
(
N
)
=
r
(
N
-
A
)
(
1
-
N
K
)
,
where A denotes the Allee threshold if A
>
0. Much work has been done on strong Allee effects; see, Hilker et al. [16, 17], Shi et al. [18, 19], and González-Olivares et al. [20, 21] and their citations. A wide range of mechanisms which may result in Allee effects are considered and discovered, from fertilization efficiency in sessile organisms to pollen limitation in plants and to cooperative hunting in animals. So far, one of the most important and probably the most studied mechanisms for Allee effects is called made-finding; that is, individuals in a population fail to find a suitable mate at low density, thus resulting in fewer reproductive outputs and examples include Glanville fritillary butterfly [15], sheep ticks, and whales [15, 22], and such mechanism may cause a mate-finding (component) Allee effect potentially. Mathematically, a mate-finding process can be modeled by a female mating rate which has positive dependence on density, denoted by M
(
N
). The hyperbolic function
(5)
M
(
N
)
=
N
N
+
θ
,
where θ scales the mate-finding Allee effects, is widely used [15, 22–24].In our paper, we aim to study a predator-prey system which is subject to the ratio-dependent functional response and a mate-finding Allee effect, and we believe it will give people a better understanding of Allee effects in predator-prey systems. More precisely, we consider the following system:(6)
d
N
d
τ
=
b
N
N
θ
+
N
-
d
N
(
1
+
N
K
)
-
a
N
P
N
+
a
h
P
,
d
P
d
τ
=
e
a
N
P
N
+
a
h
P
-
m
P
,
with the initial conditions N
(
0
)
>
0, P
(
0
)
>
0. By introducing the dimensionless variables given by t
=
d
τ, x
=
N
/
K, and y
=
a
h
P
/
K, then the system (6) becomes
(7)
d
x
d
t
=
g
x
x
δ
+
x
-
x
(
1
+
x
)
-
α
x
y
x
+
y
:
=
F
1
(
x
,
y
)
,
d
y
d
t
=
α
1
x
y
x
+
y
-
β
y
:
=
F
2
(
x
,
y
)
,
where g
=
b
/
d, α
=
1
/
d
h, α
1
=
a
e
/
d, β
=
m
/
d, and δ
=
θ
/
K (relative strength of the mate-finding Allee effect). Now we have rescaled the system from (6) to (7) by reducing the number of parameters from seven to five. For the sake of consistency, we make the same dimensionless transformation of g
=
b
/
d as that in the system studied by Pavlová et al. [25]. Note that b
>
d
>
0, so(
A
0
)
g
=
b
/
d
>
1.We assume that (
A
0
) holds in our paper and all the other parameters are positive; that is, δ
,
α
,
α
1
,
β
>
0.We organize the paper as follows. In Section2, we show the system is dissipative and the stability and bifurcation analysis of complicated equilibrium E
0, predator-free equilibria, and positive equilibria are given. We discuss our findings and summarize our conclusions in Section 3, and numerical simulations are also carried out to support our findings.


## 2. Stability and Bifurcation Analysis

### 2.1. Dissipativetity

Note thatx
y
/
(
x
+
y
) is not defined at (
0,0
), and we can redefine the derivative as follows:
(8)
F
1
(
x
,
y
)
=
F
2
(
x
,
y
)
=
0
,
(
x
,
y
)
=
(
0,0
)
.
Finally we define F
∶
=
(
F
1
F
2
). Clearly, with the extended definition, there globally exists a unique solution of system (7) for any given nonnegative initial condition. Obviously, R
+
2
=
{
(
x
,
y
)
∣
x
,
y
⩾
0
} is an invariant set. The following lemma shows that system (7) is dissipative.Lemma 1.
Let(
x
(
t
)
,
y
(
t
)
) be a solution of system (7); then, one has
(9)
limsup
⁡
t
→
+
∞
(
x
(
t
)
+
α
α
1
y
(
t
)
)
⩽
(
g
+
β
-
1
)
2
4
β
.Proof.
LetV
(
t
)
=
x
(
t
)
+
(
α
/
α
1
)
y
(
t
). Differentiating V, one has
(10)
V
′
(
t
)
=
g
x
2
x
+
δ
-
x
2
+
(
β
-
1
)
x
-
β
V
(
t
)
<
-
x
2
+
(
g
+
β
-
1
)
x
-
β
V
(
t
)
⩽
(
g
+
β
-
1
)
2
4
-
β
V
(
t
)
.
Thus, we have limsup
⁡
t
→
+
∞
⁡
V
(
t
)
⩽
(
g
+
β
-
1
)
2
/
4
β and system (7) is dissipative. This completes the proof.


### 2.2. Complicated EquilibriumE
0
(
0,0
) and Its Stability

Since the Jacobian matrix cannot be evaluated atE
0
(
0,0
) (noting that F
1 and F
2 are not differentiable at (
0,0
)), the classical stability anlysis methods are not applied. We apply the algorithm presented in [6] to investigate the structure of a neighborhood of E
0
(
0,0
) and show that system (7) has a stable E
0
(
0,0
) for all system parameters. Such result is obtained numerically in Figures 1, 2, 4, and 6(a).Phase portraits of system (7) if α
1
<
β and for g
=
1.44, δ
1
=
0.04, α
=
0.4, α
1
=
0.1, and β
=
0.3.
(a)
δ
=
0.042, δ
>
δ
1
(b)
δ
=
0.04, δ
=
δ
1
(c)
δ
=
0.032, δ
<
δ
1Phase portraits of system (7) if α
>
α
1 and for g
=
1.44, α
=
0.8, α
1
=
0.4, and β
=
0.3.
(a)
A saddleE
1
* and a stable E
2
* (δ
=
0.008)
(b)
A saddleE
1
* and a stable E
2
* surrounded by a homoclinic loop (δ
=
0.0092006)
(c)
A saddleE
1
* and an unstable E
2
* (δ
=
0.0099)For the sake of simplicity in calculation, we focus on the following system which is equivalent to system (7):
(11)
d
x
d
t
=
-
x
3
+
(
g
-
1
-
δ
)
x
2
-
δ
x
-
α
x
2
y
+
α
δ
x
y
x
+
y
,
d
y
d
t
=
α
1
x
2
y
+
α
1
δ
x
y
x
+
y
-
β
y
(
x
+
δ
)
.Definition 2.
A vector fieldW
(
x
,
y
)
=
P
(
x
,
y
)
(
∂
/
∂
x
)
+
Q
(
x
,
y
)
(
∂
/
∂
y
) is nondegenerate if it satisfies the following.(A1)
PolynomialsP
n
(
x
,
y
) and Q
n
(
x
,
y
) have no common factors of the form μ
x
+
ν
y, where P
(
x
,
y
)
=
P
n
(
x
,
y
)
+
P
~
(
x
,
y
), Q
(
x
,
y
)
=
Q
n
(
x
,
y
)
+
Q
~
(
x
,
y
), and at least one of the constants μ, ν is nonzero.
(A2)
PolynomialL
(
x
,
y
) has no factors of the form (
μ
x
+
ν
y
)
k, where k
>
1, and L
(
x
,
y
)
=
x
Q
n
(
x
,
y
)
-
y
P
n
(
x
,
y
).By the time scale changed
t
→
d
t
/
(
x
+
y
), the system (11) takes the form
(12)
d
x
d
t
=
-
x
4
-
x
3
y
+
(
g
-
1
-
δ
)
x
3
+
(
g
-
1
-
δ
-
α
)
x
2
y
-
δ
x
2
-
(
α
+
1
)
δ
x
y
:
=
P
,
d
y
d
t
=
(
α
1
-
β
)
x
2
y
-
β
x
y
2
+
(
α
1
-
β
)
δ
x
y
-
β
δ
y
2
:
=
Q
.
Let
(13)
P
2
(
x
,
y
)
=
-
δ
x
2
-
(
α
+
1
)
δ
x
y
,
Q
2
(
x
,
y
)
=
(
α
1
-
β
)
δ
x
y
-
β
δ
y
2
,
P
~
(
x
,
y
)
=
-
x
4
+
x
3
y
+
(
g
-
1
-
δ
)
x
3
+
(
g
-
1
-
δ
-
α
)
x
2
y
,
Q
~
(
x
,
y
)
=
(
α
1
-
β
)
x
2
y
-
β
x
y
2
,
and then one has
(14)
P
(
x
,
y
)
=
P
2
(
x
,
y
)
+
P
~
(
x
,
y
)
,
P
~
(
x
,
y
)
=
o
(
|
x
,
y
|
3
)
,
Q
(
x
,
y
)
=
Q
2
(
x
,
y
)
+
Q
~
(
x
,
y
)
,
Q
~
(
x
,
y
)
=
o
(
|
x
,
y
|
3
)
,
L
(
x
,
y
)
=
x
Q
2
(
x
,
y
)
-
y
P
2
(
x
,
y
)
L
i
(
x
,
y
)
=
δ
x
y
[
(
α
1
-
β
+
1
)
x
+
(
α
-
β
+
1
)
y
]
.
It is easy to see that the vector field W
(
x
,
y
) of (12) is nondegenerate. After applying the blow-up transformations (
x
,
y
)
→
(
x
,
u
) with
(15)
u
=
y
x
,
x
≠
0
,
and the time change d
t
→
x
d
t, one has
(16)
d
x
d
t
=
x
P
2
(
1
,
u
)
+
G
1
(
x
,
u
)
,
d
u
d
t
=
L
1
(
u
)
+
G
2
(
x
,
u
)
,
where
(17)
L
1
(
u
)
=
L
(
1
,
u
)
=
δ
u
[
(
α
1
-
β
+
1
)
+
(
α
-
β
+
1
)
u
]
,
G
1
(
x
,
u
)
=
P
~
(
x
,
u
x
)
x
,
G
2
(
x
,
u
)
=
(
Q
~
(
x
,
u
x
)
-
P
~
(
x
,
u
x
)
u
)
x
2
.
Thus, if (
α
-
β
+
1
)
(
α
1
-
β
+
1
)
>
0, L
1
(
u
) has only one nonnegative root u
1
=
0, and consequently, system (16) has one equilibrium O
1
(
0
,
u
1
). L
1
(
u
) has two nonnegative roots: u
1
=
0 and u
2
=
-
(
α
1
-
β
+
1
)
/
(
α
-
β
+
1
) if (
α
-
β
+
1
)
(
α
1
-
β
+
1
)
<
0, and thus system (16) has two equilibria O
1
(
0
,
u
1
) and O
2
(
0
,
u
2
) on u-axis. We state the following lemma.Lemma 3.
We have the following results.(1)
Ifα
-
β
+
1
<
0 and α
1
-
β
+
1
<
0, O
1
(
0
,
u
1
) is a single stable node.
(2)
Ifα
-
β
+
1
>
0 and α
1
-
β
+
1
>
0, O
1
(
0
,
u
1
) is a single saddle.
(3)
Ifα
-
β
+
1
>
0 and α
1
-
β
+
1
<
0, O
1
(
0
,
u
1
) is a stable node and O
2
(
0
,
u
2
) is a saddle.
(4)
Ifα
-
β
+
1
<
0 and α
1
-
β
+
1
>
0, O
1
(
0
,
u
1
) is a saddle and O
2
(
0
,
u
2
) is a stable node.Proof.
We only verify (3) since the other three results can be shown in a similar way. It is easy to see that(18)
P
2
(
1
,
u
1
)
=
-
δ
<
0
,
L
1
′
(
u
1
)
=
δ
(
α
1
-
β
+
1
)
,
P
2
(
1
,
u
2
)
=
δ
(
α
α
1
-
α
β
+
α
1
)
α
-
β
+
1
,
L
1
′
(
u
2
)
=
-
δ
(
α
1
-
β
+
1
)
.
If α
-
β
+
1
>
0 and α
1
-
β
+
1
<
0, that is, α
1
<
β
-
1
<
α, then L
1
′
(
u
1
)
<
0, L
1
′
(
u
2
)
>
0, and
(19)
P
2
(
1
,
u
2
)
=
δ
(
α
α
1
-
α
β
+
α
1
)
α
-
β
+
1
<
δ
(
-
α
+
α
1
)
α
-
β
+
1
<
0
.
Based on [6, Proposition 3], we can conclude that that O
1 is a stable node and O
2 is a saddle, respectively.We obtain the following system in a similar way by making transformations(
x
,
y
)
→
(
v
,
y
) with
(20)
v
=
x
y
,
y
≠
0
,
and d
t
→
y
d
t; then, we obtain
(21)
d
v
d
t
=
L
2
(
v
)
+
J
1
(
v
,
y
)
,
d
y
d
t
=
y
Q
2
(
v
,
1
)
+
J
2
(
v
,
y
)
,
where
(22)
L
2
(
v
)
=
-
L
(
v
,
1
)
=
-
δ
v
[
(
α
1
-
β
+
1
)
v
+
(
α
-
β
+
1
)
]
,
J
1
(
v
,
y
)
=
(
P
~
(
v
y
,
y
)
-
Q
~
(
v
y
,
y
)
v
)
y
2
,
J
2
(
v
,
y
)
=
Q
~
(
v
y
,
y
)
y
.
Thus, if (
α
-
β
+
1
)
(
α
1
-
β
+
1
)
>
0, system (21) has only one nonnegative equilibrium O
1
*(0,0). The system (21) has two equilibria O
1
*
(
0,0
) and O
2
*
(
v
2
,
0
) on the v-axis if (
α
-
β
+
1
)
(
α
1
-
β
+
1
)
<
0. However, O
2
* corresponds to O
2; thus, we do not need to study it again. O
1
* is a new equilibrium and it does not exist in the (
x
,
u
)-plane. Note that L
2
′
(
0
)
=
-
δ
(
α
+
1
-
β
) and Q
2
(
0,1
)
=
-
β
δ
<
0. Again, based on [6, Proposition 3], then comes the following conclusion.Lemma 4.
We say that(1)
ifα
-
β
+
1
<
0, O
1
*
(
0,0
) is a saddle;
(2)
ifα
-
β
+
1
>
0, O
1
*
(
0,0
) is a stable node.According to Lemmas3 and 4, it then follows from [6] that the following result holds for the origin E
0
(
0,0
) of system (7), which is equivalent to system (11).Theorem 5.
If(
α
-
β
+
1
)
(
α
1
-
β
+
1
)
≠
0, then there is at least one attracting parabolic sector in the neighborhood of (
0,0
) in the first quadrant, and no hyperbolic or elliptic sector can be found; that is, the complicated equilibrium E
0
(
0,0
) is stable.Remark 6.
If(
α
-
β
+
1
)
(
α
1
-
β
+
1
)
>
0, the neighborhood of E
0
(
0,0
) in the first quadrant has a parabolic (attracting) sector, but if (
α
-
β
+
1
)
(
α
1
-
β
+
1
)
<
0, there are two parabolic (attracting) sectors; in fact, we can regard those two sectors as a whole attracting one. If (
α
-
β
+
1
)
(
α
1
-
β
+
1
)
=
0, O
i, i
=
1,2 are also complicated equilibria and needed to apply the blow-up transformation again. Now we can see that the relative strength of mate-finding Allee effect δ greatly changes the structure of the neighborhood of E
0
(
0,0
) in our system (7), because P
2
(
1
,
u
1
)
=
-
δ and Q
2
(
0,1
)
=
-
β
δ always stay negative and O
1
(
0,0
) and O
1
*
(
0,0
) cannot be unstable; thus, the diversity of phase portraits around E
0
(
0,0
) of ratio-dependent prey-predator systems decreases. More importantly, the system becomes fragile and even has a higher risk of extinction due to mate-finding Allee effect.


### 2.3. Existence and Stability of Predator-Free Equilibria

In order to analyze the predator-free equilibria of the system (7), we set
(23)
Δ
1
=
(
g
-
δ
-
1
)
2
-
4
δ
.
We can notice that, for y
=
0, there exist two boundary equilibria E
10
=
(
x
1
,
0
) and E
20
=
(
x
2
,
0
) if 0
<
δ
<
g
-
1 and Δ
1
>
0; that is,(
B
1
)
0
<
δ
<
δ
1,whereδ
1
:
=
(
g
-
1
)
2. Let x
i, i
=
1,2 be the roots of the quadratic equation x
2
+
(
1
+
δ
-
g
)
x
+
δ
=
0; that is,
(24)
x
1
=
1
2
(
g
-
δ
-
1
-
(
g
-
δ
-
1
)
2
-
4
δ
)
,
x
2
=
1
2
(
g
-
δ
-
1
+
(
g
-
δ
-
1
)
2
-
4
δ
)
,
and they appear in pairs. Besides, there is no equilibrium on y-axis.The Jacobian matrix evaluated atE
i
0, i
=
1,2 is given by
(25)
J
i
0
=
(
x
i
(
x
i
+
δ
)
2
(
g
δ
-
(
x
i
+
δ
)
2
)
-
α
0
α
1
-
β
)
,
and one eigenvalue is λ
1
(
E
i
0
)
=
(
x
i
/
(
x
i
+
δ
)
2
)
(
g
δ
-
(
x
i
+
δ
)
2
), i
=
1,2, and the second eigenvalue is λ
2
(
E
i
0
)
=
α
1
-
β, i
=
1,2. Note that the sign of λ
1
(
E
i
0
)
,
i
=
1,2, depends on ψ
i
:
=
g
δ
-
δ
-
x
i, i
=
1,2. Using the expression of x
i, i
=
1,2, from (24), we can have
(26)
ψ
1
=
1
2
(
2
g
δ
-
δ
-
g
+
1
+
(
g
-
δ
-
1
)
2
-
4
δ
)
,
ψ
2
=
1
2
(
2
g
δ
-
δ
-
g
+
1
-
(
g
-
δ
-
1
)
2
-
4
δ
)
.
Under the restriction (
B
1
), it follows that
(27)
ψ
1
>
δ
(
g
-
δ
-
1
)
>
0
,
ψ
2
=
1
2
(
-
(
g
-
δ
)
2
+
1
-
Δ
1
)
<
-
Δ
1
2
<
0
.
Hence λ
1
(
E
10
) is always positive, while λ
1
(
E
20
) is negative. By using the facts above, we have the following theorem.Theorem 7.
Assume (B
1) holds. Then(1)
ifα
1
<
β, E
10 is a saddle and E
20 is a stable node (Figure 1(c)),
(2)
ifα
1
>
β, E
10 is an unstable node and E
20 is a saddle (Figure 2).Remark 8.
Ifδ
=
δ
1, system (7) has a unique predator-free boundary equilibrium point E
*
0
=
(
x
*
0
,
0
) (also known as the saddle-node equilibrium), where
(28)
x
*
0
=
1
2
(
g
-
δ
1
-
1
)
=
δ
1
.
Figure 1(b) shows the instantaneous equilibrium E
*
0 in the case of α
1
<
β and the behaviors of trajectories are divided by the stable manifold W
s
(
E
*
0
) of E
*
0 (the green curve).


### 2.4. Positive Equilibria and Their Stability and Bifurcation Results

#### 2.4.1. Existence and Stability of Positive Equilibria

The interior (positive) equilibria can be evaluated by the intersections of the zero isoclines(29)
g
x
x
+
δ
-
(
1
+
x
)
=
α
y
x
+
y
,
α
1
x
x
+
y
=
β
,
in the first quadrant. The predator zero isocline, with the slope (
α
1
-
β
)
/
β, is a straight line passing through the origin, and it lies in the first quadrant if the following restriction holds(
C
1
)
α
1
>
β.Besides, we can define(30)
σ
:
=
1
+
α
-
α
β
α
1
.
Note that σ
>
1 under the restriction (
C
1
). And for x
2
+
(
δ
+
σ
-
g
)
x
+
δ
σ
=
0, we also set
(31)
Δ
2
=
(
g
-
δ
-
σ
)
2
-
4
δ
σ
.
If 0
<
δ
<
g
-
σ and Δ
2
>
0, that is,(
B
2
)
0
<
δ
<
δ
2,whereδ
2
:
=
(
g
-
σ
)
2, there are two positive equilibria E
1
*
=
(
x
1
*
,
y
1
*
) and E
2
*
=
(
x
2
*
,
y
2
*
), which also appear simultaneously, with the coordinates
(32)
x
1
*
=
1
2
(
g
-
δ
-
σ
-
(
g
-
δ
-
σ
)
2
-
4
δ
σ
)
,
x
2
*
=
1
2
(
g
-
δ
-
σ
+
(
g
-
δ
-
σ
)
2
-
4
δ
σ
)
,
(33)
y
i
*
=
α
1
-
β
β
x
i
*
,
i
=
1,2
.
Note that 0
<
δ
2
<
δ
1, it suggests that if (
B
2
) holds, then (
B
1
) holds; that is, if there exists any interior equilibrium in the first quadrant, then two predator-free equilibria will be detected. The Jacobian matrix J
i
*, i
=
1,2 evaluated at an interior equilibrium is given by
(34)
J
i
*
=
(
g
δ
x
i
*
(
x
i
*
+
δ
)
2
-
x
i
*
+
α
x
i
*
y
i
*
(
x
i
*
+
y
i
*
)
2
-
α
x
i
*
2
(
x
i
*
+
y
i
*
)
2
α
1
y
i
*
2
(
x
i
*
+
y
i
*
)
2
-
α
1
x
i
*
y
i
*
(
x
i
*
+
y
i
*
)
2
)
,
det
⁡
(
J
i
*
)
=
α
1
x
i
*
2
y
i
*
[
(
x
i
*
+
δ
)
2
-
g
δ
]
(
x
i
*
+
y
i
*
)
2
(
x
i
*
+
δ
)
2
.
Similarly, we have det
⁡
(
J
1
*
)
<
0 and det
⁡
(
J
2
*
)
>
0. Moreover,
(35)
tr
⁡
(
J
2
*
)
=
g
δ
x
2
*
(
x
2
*
+
δ
)
2
-
x
2
*
+
β
(
α
-
α
1
)
(
α
1
-
β
)
α
1
2
,
and it is interesting to find that tr
⁡
(
J
2
*
) is an increasing function of δ. In fact, we can have
(36)
d
x
2
d
δ
=
-
1
2
(
1
+
g
-
δ
+
σ
(
g
-
δ
+
σ
)
2
-
4
δ
σ
)
<
0
,
and we can define x
2
′
:
=
d
x
2
/
d
δ. Then, it follows that
(37)
dtr
(
J
2
*
)
d
δ
=
x
2
′
[
g
δ
(
x
2
+
δ
)
2
-
1
]
+
x
2
g
(
x
2
+
δ
)
3
[
x
2
-
δ
-
2
δ
x
2
′
]
>
0
,
since
(38)
-
δ
-
2
δ
x
2
′
=
-
δ
+
δ
(
1
+
g
-
δ
+
σ
(
g
-
δ
+
σ
)
2
-
4
δ
σ
)
=
δ
(
g
-
δ
+
σ
)
(
g
-
δ
+
σ
)
2
-
4
δ
σ
>
0
.
Moreover, we can have the following conclusion.Theorem 9.
Assume that(
B
2
) and (
C
1
) hold. Then one has the following.(1)
E
1
* is a saddle.
(2)
Ifα
⩽
α
1, E
2
* is stable.
(3)
Ifα
>
α
1, E
2
* can be stable or unstable depending on δ (Figure 2).


#### 2.4.2. Hopf Bifurcation

Here assume that(
B
2
) and (
C
1
) hold, and we consider δ as the bifurcation parameter to discuss the Hopf bifurcation. It is easy to know that the Hopf bifurcation can occur only at E
2
* as E
1
* is always a saddle point. Then we try to show the existence of δ
=
δ
h such that tr
⁡
(
J
2
*
)
=
0 and det
⁡
(
J
2
*
)
>
0. From the expression of tr
⁡
(
J
2
*
) and some numerical attempts, we find that other parameters g, α, α
1, and β still have an important impact on the qualitative properties of the Hopf bifurcation though only δ is considered as the bifurcation parameter. The following assumption is used to guarantee the exsitence of the Hopf bifurcation with bifurcation parameter δ:(
C
2
)
α
>
α
1, (
1
-
(
β
/
α
1
)
)
(
α
-
β
+
(
α
β
/
α
1
)
)
<
g
-
1.Assumption (
C
2
) implies that tr
⁡
(
J
2
*
)
|
δ
=
0
<
0 and tr
⁡
(
J
2
*
)
|
δ
=
δ
2
>
0, and we have known that tr
⁡
(
J
2
*
) is a continuous increasing function of δ in [
0
,
+
∞
); thus, there exists δ
=
δ
h such that tr
⁡
(
J
2
*
)
=
0. Furthermore, the transversality condition (
d
/
d
δ
)
tr
⁡
(
J
2
*
)
|
δ
=
δ
h
>
0 is satisfied. Using the facts above, system (7) can undergo a Hopf bifurcation at E
2
* for δ
=
δ
h if (
C
2
) holds.Now we try to discuss the stability of the limit cycle of system (7) as a Hopf bifurcation occurs by computing the first Lyapunov coefficient l [26] at the E
2
*. Transformations u
1
=
x
-
x
2
* and u
2
=
y
-
y
2
* are used to translate the equilibrium E
2
*
(
x
2
*
,
y
2
*
) of the system (7) to the origin. Then we get
(39)
u
˙
1
=
a
10
u
1
+
a
01
u
2
+
a
20
u
1
2
+
a
11
u
1
u
2
+
a
02
u
2
2
+
a
03
u
1
3
+
a
21
u
1
2
u
2
+
a
12
u
1
u
2
2
+
a
03
u
2
3
+
R
1
(
u
1
,
u
2
)
,
u
˙
2
=
b
10
u
1
+
b
01
u
2
+
b
20
u
1
2
+
b
11
u
1
u
2
+
b
02
u
2
2
+
b
03
u
1
3
+
b
21
u
1
2
u
2
+
b
12
u
1
u
2
2
+
b
03
u
2
3
+
R
2
(
u
1
,
u
2
)
,
where
(40)
a
10
=
x
2
*
g
δ
(
x
2
*
+
δ
)
2
+
α
x
2
*
y
2
*
(
x
2
*
+
y
2
*
)
2
-
x
2
*
,
a
01
=
-
α
x
2
*
2
(
x
2
*
+
y
2
*
)
2
,
a
20
=
g
δ
2
(
x
2
*
+
δ
)
3
+
α
y
2
*
2
(
x
2
*
+
y
2
*
)
3
-
1
,
a
11
=
-
2
α
x
2
*
y
2
*
(
x
2
*
+
y
2
*
)
3
,
a
02
=
-
α
x
2
*
y
2
*
(
x
2
*
+
y
2
*
)
3
,
a
30
=
-
g
δ
2
(
x
2
*
+
δ
)
4
-
α
y
2
*
2
(
x
2
*
+
y
2
*
)
4
,
a
21
=
2
α
x
2
*
y
2
*
-
α
y
2
*
2
(
x
2
*
+
y
2
*
)
4
,
a
12
=
2
α
x
2
*
y
2
*
-
α
x
2
*
2
(
x
2
*
+
y
2
*
)
4
,
a
03
=
-
α
x
2
*
2
(
x
2
*
+
y
2
*
)
4
,
b
10
=
α
1
y
2
*
2
(
x
2
*
+
y
2
*
)
2
,
b
01
=
-
α
1
x
2
*
y
2
*
(
x
2
*
+
y
2
*
)
2
,
b
20
=
-
α
1
y
2
*
2
(
x
2
*
+
y
2
*
)
3
,
b
11
=
2
α
1
x
2
*
y
2
*
(
x
2
*
+
y
2
*
)
3
,
b
02
=
-
α
1
x
2
*
2
(
x
2
*
+
y
2
*
)
3
,
b
30
=
α
1
y
2
*
2
(
x
2
*
+
y
2
*
)
4
,
b
21
=
-
2
α
1
x
2
*
y
2
*
+
α
1
y
2
*
2
(
x
2
*
+
y
2
*
)
4
,
b
12
=
-
2
α
1
x
2
*
y
2
*
+
α
1
x
2
*
2
(
x
2
*
+
y
2
*
)
4
,
b
03
=
α
1
x
2
*
2
(
x
2
*
+
y
2
*
)
4
,
and R
i
(
x
,
y
), i
=
1,2 are power series in the powers of u
1
j
u
2
k satisfying j
+
k
⩾
4. The Lyapunov number l (as defined in [26]) is given by
(41)
l
=
-
3
π
2
a
01
Δ
3
/
2
×
(
[
a
10
b
10
(
a
11
2
+
a
11
b
02
+
a
02
b
11
)
ffffiffii
+
a
10
a
01
(
b
11
2
+
a
20
b
11
+
a
11
b
02
)
ffffiffii
+
b
10
2
(
a
11
a
02
+
2
a
02
b
02
)
-
2
a
10
b
10
(
b
02
2
-
a
20
a
02
)
fffifffii
-
2
a
10
a
01
(
a
20
2
-
b
20
b
02
)
-
a
01
2
(
2
a
20
b
20
+
b
11
b
20
)
fffffifii
+
(
a
01
b
10
-
2
a
10
2
)
(
b
11
b
02
-
a
11
a
20
)
]
ffiffifif
-
(
a
10
2
+
a
01
b
10
)
[

3
(
b
10
b
03
-
a
01
a
30
)
fffffffiffffffffffffffffffffiffff
+
2
a
10
(
a
21
+
b
12
)
ffffffffifffffffffffffffffffiffff
+
(
b
10
a
12
-
a
01
b
21
)
]
a
11
2
)
,
where Δ
=
a
10
b
01
-
a
01
b
10. Since the expression Lyapunov number for l is very complex, we fail to discuss the sign of l precisely though a subcritical Hopf bifurcation has been found numerically with the Lyapunov number l
=
5.545875
>
0 (Figure 3). Obviously, the assumption (
C
2
) is satisfied in Figure 3. We need to require δ
<
δ
1 (green line) and δ
<
δ
2 (blue line) to guarantee the existence of E
2
*; that is, (
B
1
) and (
B
2
) hold (Figure 3(a)). tr
⁡
(
J
2
*
) is a continuous increasing function of δ. If 0
<
δ
<
δ
h, system (7) has a stable E
2
* (Figure 2(a)). As δ increases, a homoclinic loop is created by joining the stable and unstable manifolds of the saddle E
1
* at some δ
~ (δ
~
<
δ
h) with a stable E
2
* inside (Figure 2(b)). Then as δ gets larger, an unstable limit cycle appears and coexists with a stable E
2
* (Figure 4). However, the limit cycle begins to shrink with the increasing δ. The shrinking unstable limit cycle (Figures 3(b) and 3(c)) disappears as δ passes through δ
h
(
δ
h is also the root of equation tr
⁡
(
J
2
*
)
=
0, δ
h
=
0.009724
). Then δ
h
<
δ
<
δ
2, E
2
* becomes unstable (Figure 2(c)).Bifurcation structure ofE
2
* for Hopf bifurcation and phase portraits of system (7) for g
=
1.44, α
=
0.8, α
1
=
0.4, and β
=
0.3.
(a)
Bifurcation structure
(b)
δ
=
0.009724
(
δ
=
δ
h
), the unstable limit cycle when a subcritical Hopf bifurcation occurs
(c)
δ
=
0.009724
(
δ
=
δ
h
), the trajectory around the unstable limit cycle when a subcritical Hopf bifurcation occursAn unstable limit cycle coexists with a stableE
2
*, for g
=
1.44, δ
=
0.0095, α
=
0.84, α
1
=
0.4, and β
=
0.3.
(a)
Inside the limit cycle
(b)
Outside the limit cycle


#### 2.4.3. Saddle-Node Bifurcation

We know that wheng
*, δ
*, α
*, α
1
*, and β
* satisfy (
C
1
) and the following equation
(42)
δ
-
(
g
-
σ
)
2
=
0
,
where σ
=
1
+
α
-
α
β
/
α
1, then there is only one interior equilibrium E
* of system (7), whose coordinates are given by
(43)
x
*
=
1
2
(
g
*
-
δ
*
-
σ
*
)
=
δ
*
σ
*
,
y
*
=
α
1
*
-
β
*
2
β
*
(
g
*
-
δ
*
-
σ
*
)
=
α
1
*
-
β
*
β
*
δ
*
σ
*
.
And the Jacobian matrix evaluated at E
* is
(44)
J
*
=
(
α
*
x
*
y
*
(
x
*
+
y
*
)
2
-
α
*
x
*
2
(
x
*
+
y
*
)
2
α
1
*
y
*
2
(
x
*
+
y
*
)
2
-
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
)
.
Next we choose the relative strength of mate-finding Allee effect δ as the bifurcation parameter; then, we have the following theorem.Theorem 10.
Assume (C
1) holds. For δ
=
δ
* and α
1
*
≠
α
*, a saddle-node bifurcation occurs at the unique positive equilibrium E
* of system (7).Proof.
First, we havedet
⁡
(
J
*
)
=
0; therefore, J
* has an eigenvalue λ
*
=
0, and if α
1
*
≠
α
*, then λ
*
=
0 is simple. Let W
1 and W
2 be the eigenvectors corresponding to the λ
*
=
0 for J
* and J
*
T, respectively. Then we obtain
(45)
W
1
=
(
1
α
1
*
-
β
*
β
*
)
:
=
(
w
11
w
12
)
,
W
2
=
(
1
-
α
*
β
*
α
1
*
(
α
1
*
-
β
*
)
)
.
From the expressions for W
1 and W
2, we get
(46)
W
2
T
F
δ
(
(
x
*
,
y
*
)
;
δ
*
)
=
-
g
*
x
*
2
(
x
*
+
δ
*
)
2
<
0
(
≠
0
)
,
W
2
T
D
2
F
(
W
1
,
W
1
)
=
-
2
x
*
x
*
+
δ
*
<
0
(
≠
0
)
,
where
(47)
F
δ
=
∂
F
∂
δ
,
D
2
F
(
W
1
,
W
1
)
=
∂
2
F
∂
x
2
w
11
w
11
+
∂
2
F
∂
x
∂
y
w
11
w
12
+
∂
2
F
∂
y
∂
x
w
12
w
11
+
∂
2
F
∂
y
2
w
12
w
12
.
Thus by Sotomayor’s theorem [26], system (7) undergoes a saddle-node bifurcation at E
* as δ passes δ
=
δ
* if α
1
*
≠
α
*.A numerical example has been carried out to demonstrate a saddle-node bifurcation of system (7) (Figure 5). The green curves represent the smaller unstable equilibria. The stability of the larger equilibria changes from stable (blue curves) to unstable (black curves) as δ passes through δ
h
(
δ
h
=
0.009724
).A saddle-node bifurcation of system (7) at δ
=
0.010932 (mark as SN), for g
=
1.44, α
=
0.8, α
1
=
0.4, and β
=
0.3. A Hopf bifurcation also occurs.
(a)
Prey population scenario
(b)
Predator population scenarioBifurcation structure of system (7) in (
δ
,
α
1
) with a Bogdanov-Takens bifurcation and phase portraits for g
=
1.44, δ
=
0.022859, α
=
0.4, α
1
=
0.4, and β
=
0.3.
(a)
Bifurcation structure
(b)
A cuspE
*


#### 2.4.4. Bogdanov-Takens Bifurcation

Let us focus onE
* again. After knowing that a saddle-node bifurcation can occur at E
* if δ
=
δ
* and α
1
*
≠
α
*, we go on to consider the case of δ
=
δ
2, det
⁡
(
J
*
)
=
0, and tr
⁡
(
J
*
)
=
0. If det
⁡
(
J
*
)
=
0 and tr
⁡
(
J
*
)
=
0, then the Jacobian matrix of E
* has λ
*
1
=
λ
*
2
=
0. When det
⁡
(
J
*
)
=
0, it is easily shown that if α
1
*
=
α
* then tr
⁡
(
J
*
)
=
0, and hence we choose the relative mate-finding strength δ and the predator growing ability α
1 as two bifurcation parameters. Then the following statement holds.Theorem 11.
Let (C
1) hold. The system (7) undergoes a Bogdanov-Takens bifurcation around the equilibrium point E
* when δ
=
δ
* and α
1
*
=
α
*.Proof.
Following the ideas in [27, 28], we consider the neighborhood of (
δ
*
,
α
1
*
); that is, δ
=
δ
*
+
r
1 and α
1
=
α
1
*
+
r
2, where r
i, i
=
1,2, are sufficient small, and system (7) becomes
(48)
d
x
d
t
=
g
*
x
2
δ
*
+
r
1
+
x
-
d
*
x
(
1
+
x
)
-
α
1
*
x
y
x
+
y
,
d
y
d
t
=
(
α
1
*
+
r
2
)
x
y
x
+
y
-
β
*
y
.
We make z
1
=
x
-
x
*, z
2
=
y
-
y
* for translating E
*
=
(
x
*
,
y
*
) to the origin (
0,0
), and we obtain
(49)
z
˙
1
=
[
-
g
*
x
*
2
(
x
*
+
δ
*
)
2
r
1
+
ξ
1
(
r
1
)
]
+
[
-
2
g
*
δ
*
x
*
(
x
*
+
δ
*
)
3
r
1
+
ξ
2
(
r
1
)
+
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
]
z
1
-
α
1
*
x
*
2
(
x
*
+
y
*
)
2
z
2
+
[
g
*
δ
*
(
2
x
*
-
δ
*
)
(
x
*
+
δ
*
)
4
r
1
+
ξ
3
(
r
1
)
-
α
1
*
x
*
2
(
x
*
+
y
*
)
2
z
2
ii
+
i
+
α
1
*
y
*
2
(
x
*
+
y
*
)
3
-
x
*
x
*
+
δ
*
]
z
1
2
-
2
α
1
*
x
*
y
*
(
x
*
+
y
*
)
3
z
1
z
2
+
α
1
*
x
*
2
(
x
*
+
y
*
)
3
z
2
2
+
T
1
(
z
1
,
z
2
)
:
=
Q
1
(
z
1
,
z
2
,
r
1
,
r
2
)
,
z
˙
2
=
x
*
y
*
x
*
+
y
*
r
2
+
[
α
1
*
y
*
2
(
x
*
+
y
*
)
2
+
y
*
2
r
2
(
x
*
+
y
*
)
2
]
z
1
+
[
-
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
+
x
*
2
r
2
(
x
*
+
y
*
)
2
]
z
2
-
[
α
1
*
y
*
2
(
x
*
+
y
*
)
3
+
y
*
2
(
x
*
+
y
*
)
3
r
2
]
z
1
2
+
[
2
α
1
*
x
*
y
*
(
x
*
+
y
*
)
3
+
2
x
*
y
*
(
x
*
+
y
*
)
3
r
2
]
z
1
z
2
-
[
α
1
*
x
*
2
(
x
*
+
y
*
)
3
+
x
*
2
(
x
*
+
y
*
)
3
r
2
]
z
2
2
+
T
2
(
z
1
,
z
2
)
:
=
Q
2
(
z
1
,
z
2
,
r
1
,
r
2
)
,
where ξ
i, i
=
1,2
,
3 are polynomial functions of r
1 with r
1
p satisfying p
⩾
2, and T
1
(
z
1
,
z
2
) and T
2
(
z
1
,
z
2
) are C
∞ functions z
1
j
z
2
k satisfying j
+
k
⩾
3.
Then we apply the affine transformation(50)
v
1
=
z
1
,
v
2
=
[
-
2
g
*
δ
*
x
*
(
x
*
+
δ
*
)
3
r
1
+
ξ
2
(
r
1
)
+
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
]
z
1
-
α
1
*
x
*
2
(
x
*
+
y
*
)
2
z
2
,
to obtain
(51)
v
˙
1
=
v
2
+
[
-
x
*
x
*
+
δ
*
+
ξ
~
1
(
r
1
)
]
v
1
2
+
ξ
~
2
(
r
1
)
v
1
v
2
+
x
*
+
y
*
α
1
*
x
*
2
v
2
2
+
T
~
1
(
v
1
,
v
2
,
r
1
,
r
2
)
,
v
˙
2
=
η
~
1
(
r
1
,
r
2
)
v
1
+
η
~
2
(
r
1
,
r
2
)
v
2
+
[
-
α
1
*
x
*
2
y
*
(
x
*
+
δ
*
)
(
x
*
+
y
*
)
2
+
η
~
3
(
r
1
,
r
2
)
]
v
1
2
+
η
~
4
(
r
1
,
r
2
)
v
1
v
2
+
[
1
x
*
+
η
~
5
(
r
1
,
r
2
)
]
v
2
2
+
T
~
2
(
v
1
,
v
2
,
r
1
,
r
2
)
,
where ξ
~
i
(
r
1
)
,
i
=
1,2
,
3 and η
~
i
(
r
1
,
r
2
)
,
i
=
1,2
,
3,4
,
5 are polynomial functions with r
1
p satisfying p
⩾
1 and r
1
p
r
2
q satisfying p
+
q
⩾
1, respectively. T
~
i, i
=
1,2, are C
∞ functions in their variables with v
1
j
v
2
k satisfying j
+
k
⩾
3.
Furthermore, ifr
1
=
r
2
=
0, we can see that(BT0)
(52)
J
(
*
,
0,0
)
=
(
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
-
α
1
*
x
*
2
(
x
*
+
y
*
)
2
α
1
*
y
*
2
(
x
*
+
y
*
)
2
-
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
)
≠
θ
2
×
2
,
(BT1)
[
-
(
x
*
/
(
x
*
+
δ
*
)
)
+
ξ
~
1
(
r
1
)
+
η
~
4
(
r
1
,
r
2
)
]
|
r
1
=
0
,
r
2
=
0
=
-
x
*
/
(
x
*
+
δ
*
)
<
0,
(BT2)
[
-
(
α
1
*
x
*
2
y
*
/
(
(
x
*
+
δ
*
)
(
x
*
+
y
*
)
2
)
)
+
η
~
3
(
r
1
,
r
2
)
]
|
r
1
=
0
,
r
2
=
0
=
-
α
1
*
x
*
2
y
*
/
(
(
x
*
+
δ
*
)
(
x
*
+
y
*
)
2
)
<
0,
(BT3)
the map(53)
(
(
z
1
z
2
)
,
(
r
1
r
2
)
)
⟶
(
(
Q
1
Q
2
)
,
tr
⁡
(
∂
Q
i
∂
z
i
)
,
det
⁡
⁡
(
∂
Q
i
∂
z
i
)
)
,
i
=
1,2
is regular at(54)
(
z
1
z
2
)
=
(
0
0
)
,
(
r
1
r
2
)
=
(
0
0
)
.According to [27, 28], the system (7) undergoes a Bogdanov-Takens bifurcation around the equilibrium point E
* when δ
=
δ
* and α
1
*
=
α
*. This completes the proof.Remark 12.
Furthermore, the quantityω [27, 28] which determines the structure of the bifurcation is given by
(55)
ω
=
sign
⁡
[
α
1
*
x
*
3
y
*
(
x
*
+
δ
*
)
2
(
x
*
+
y
*
)
2
]
=
1
.The bifurcation diagram of system (7) in (
δ
,
α
1
) is presented in Figure 6(a). There are no positive equilibria of system (7) in Domain 1. In the domain bounded by Limit point cycle (green) and Hopf curve (blue), there exist a saddle E
1
* and an unstable E
2
*. In the domain bounded by Hopf curve (blue) and Homoclinic curve (dot and magenta), there is a saddle E
1
* and a stable E
2
* surrounded by an unstable limit cycle. In Domain 2, system (7) has a saddle E
1
* and a stable E
2
*. At BT point, E
* is the unique positive equilibrium of system (7) in the first quadrant, and the local dynamics near E
* is characterized by two trajectories (green curves, Figure 6(b)), which meet at E
*, but with opposite time orientation, and also behave like a cusp-like configuration.


## 2.1. Dissipativetity

Note thatx
y
/
(
x
+
y
) is not defined at (
0,0
), and we can redefine the derivative as follows:
(8)
F
1
(
x
,
y
)
=
F
2
(
x
,
y
)
=
0
,
(
x
,
y
)
=
(
0,0
)
.
Finally we define F
∶
=
(
F
1
F
2
). Clearly, with the extended definition, there globally exists a unique solution of system (7) for any given nonnegative initial condition. Obviously, R
+
2
=
{
(
x
,
y
)
∣
x
,
y
⩾
0
} is an invariant set. The following lemma shows that system (7) is dissipative.Lemma 1.
Let(
x
(
t
)
,
y
(
t
)
) be a solution of system (7); then, one has
(9)
limsup
⁡
t
→
+
∞
(
x
(
t
)
+
α
α
1
y
(
t
)
)
⩽
(
g
+
β
-
1
)
2
4
β
.Proof.
LetV
(
t
)
=
x
(
t
)
+
(
α
/
α
1
)
y
(
t
). Differentiating V, one has
(10)
V
′
(
t
)
=
g
x
2
x
+
δ
-
x
2
+
(
β
-
1
)
x
-
β
V
(
t
)
<
-
x
2
+
(
g
+
β
-
1
)
x
-
β
V
(
t
)
⩽
(
g
+
β
-
1
)
2
4
-
β
V
(
t
)
.
Thus, we have limsup
⁡
t
→
+
∞
⁡
V
(
t
)
⩽
(
g
+
β
-
1
)
2
/
4
β and system (7) is dissipative. This completes the proof.


## 2.2. Complicated EquilibriumE
0
(
0,0
) and Its Stability

Since the Jacobian matrix cannot be evaluated atE
0
(
0,0
) (noting that F
1 and F
2 are not differentiable at (
0,0
)), the classical stability anlysis methods are not applied. We apply the algorithm presented in [6] to investigate the structure of a neighborhood of E
0
(
0,0
) and show that system (7) has a stable E
0
(
0,0
) for all system parameters. Such result is obtained numerically in Figures 1, 2, 4, and 6(a).Phase portraits of system (7) if α
1
<
β and for g
=
1.44, δ
1
=
0.04, α
=
0.4, α
1
=
0.1, and β
=
0.3.
(a)
δ
=
0.042, δ
>
δ
1
(b)
δ
=
0.04, δ
=
δ
1
(c)
δ
=
0.032, δ
<
δ
1Phase portraits of system (7) if α
>
α
1 and for g
=
1.44, α
=
0.8, α
1
=
0.4, and β
=
0.3.
(a)
A saddleE
1
* and a stable E
2
* (δ
=
0.008)
(b)
A saddleE
1
* and a stable E
2
* surrounded by a homoclinic loop (δ
=
0.0092006)
(c)
A saddleE
1
* and an unstable E
2
* (δ
=
0.0099)For the sake of simplicity in calculation, we focus on the following system which is equivalent to system (7):
(11)
d
x
d
t
=
-
x
3
+
(
g
-
1
-
δ
)
x
2
-
δ
x
-
α
x
2
y
+
α
δ
x
y
x
+
y
,
d
y
d
t
=
α
1
x
2
y
+
α
1
δ
x
y
x
+
y
-
β
y
(
x
+
δ
)
.Definition 2.
A vector fieldW
(
x
,
y
)
=
P
(
x
,
y
)
(
∂
/
∂
x
)
+
Q
(
x
,
y
)
(
∂
/
∂
y
) is nondegenerate if it satisfies the following.(A1)
PolynomialsP
n
(
x
,
y
) and Q
n
(
x
,
y
) have no common factors of the form μ
x
+
ν
y, where P
(
x
,
y
)
=
P
n
(
x
,
y
)
+
P
~
(
x
,
y
), Q
(
x
,
y
)
=
Q
n
(
x
,
y
)
+
Q
~
(
x
,
y
), and at least one of the constants μ, ν is nonzero.
(A2)
PolynomialL
(
x
,
y
) has no factors of the form (
μ
x
+
ν
y
)
k, where k
>
1, and L
(
x
,
y
)
=
x
Q
n
(
x
,
y
)
-
y
P
n
(
x
,
y
).By the time scale changed
t
→
d
t
/
(
x
+
y
), the system (11) takes the form
(12)
d
x
d
t
=
-
x
4
-
x
3
y
+
(
g
-
1
-
δ
)
x
3
+
(
g
-
1
-
δ
-
α
)
x
2
y
-
δ
x
2
-
(
α
+
1
)
δ
x
y
:
=
P
,
d
y
d
t
=
(
α
1
-
β
)
x
2
y
-
β
x
y
2
+
(
α
1
-
β
)
δ
x
y
-
β
δ
y
2
:
=
Q
.
Let
(13)
P
2
(
x
,
y
)
=
-
δ
x
2
-
(
α
+
1
)
δ
x
y
,
Q
2
(
x
,
y
)
=
(
α
1
-
β
)
δ
x
y
-
β
δ
y
2
,
P
~
(
x
,
y
)
=
-
x
4
+
x
3
y
+
(
g
-
1
-
δ
)
x
3
+
(
g
-
1
-
δ
-
α
)
x
2
y
,
Q
~
(
x
,
y
)
=
(
α
1
-
β
)
x
2
y
-
β
x
y
2
,
and then one has
(14)
P
(
x
,
y
)
=
P
2
(
x
,
y
)
+
P
~
(
x
,
y
)
,
P
~
(
x
,
y
)
=
o
(
|
x
,
y
|
3
)
,
Q
(
x
,
y
)
=
Q
2
(
x
,
y
)
+
Q
~
(
x
,
y
)
,
Q
~
(
x
,
y
)
=
o
(
|
x
,
y
|
3
)
,
L
(
x
,
y
)
=
x
Q
2
(
x
,
y
)
-
y
P
2
(
x
,
y
)
L
i
(
x
,
y
)
=
δ
x
y
[
(
α
1
-
β
+
1
)
x
+
(
α
-
β
+
1
)
y
]
.
It is easy to see that the vector field W
(
x
,
y
) of (12) is nondegenerate. After applying the blow-up transformations (
x
,
y
)
→
(
x
,
u
) with
(15)
u
=
y
x
,
x
≠
0
,
and the time change d
t
→
x
d
t, one has
(16)
d
x
d
t
=
x
P
2
(
1
,
u
)
+
G
1
(
x
,
u
)
,
d
u
d
t
=
L
1
(
u
)
+
G
2
(
x
,
u
)
,
where
(17)
L
1
(
u
)
=
L
(
1
,
u
)
=
δ
u
[
(
α
1
-
β
+
1
)
+
(
α
-
β
+
1
)
u
]
,
G
1
(
x
,
u
)
=
P
~
(
x
,
u
x
)
x
,
G
2
(
x
,
u
)
=
(
Q
~
(
x
,
u
x
)
-
P
~
(
x
,
u
x
)
u
)
x
2
.
Thus, if (
α
-
β
+
1
)
(
α
1
-
β
+
1
)
>
0, L
1
(
u
) has only one nonnegative root u
1
=
0, and consequently, system (16) has one equilibrium O
1
(
0
,
u
1
). L
1
(
u
) has two nonnegative roots: u
1
=
0 and u
2
=
-
(
α
1
-
β
+
1
)
/
(
α
-
β
+
1
) if (
α
-
β
+
1
)
(
α
1
-
β
+
1
)
<
0, and thus system (16) has two equilibria O
1
(
0
,
u
1
) and O
2
(
0
,
u
2
) on u-axis. We state the following lemma.Lemma 3.
We have the following results.(1)
Ifα
-
β
+
1
<
0 and α
1
-
β
+
1
<
0, O
1
(
0
,
u
1
) is a single stable node.
(2)
Ifα
-
β
+
1
>
0 and α
1
-
β
+
1
>
0, O
1
(
0
,
u
1
) is a single saddle.
(3)
Ifα
-
β
+
1
>
0 and α
1
-
β
+
1
<
0, O
1
(
0
,
u
1
) is a stable node and O
2
(
0
,
u
2
) is a saddle.
(4)
Ifα
-
β
+
1
<
0 and α
1
-
β
+
1
>
0, O
1
(
0
,
u
1
) is a saddle and O
2
(
0
,
u
2
) is a stable node.Proof.
We only verify (3) since the other three results can be shown in a similar way. It is easy to see that(18)
P
2
(
1
,
u
1
)
=
-
δ
<
0
,
L
1
′
(
u
1
)
=
δ
(
α
1
-
β
+
1
)
,
P
2
(
1
,
u
2
)
=
δ
(
α
α
1
-
α
β
+
α
1
)
α
-
β
+
1
,
L
1
′
(
u
2
)
=
-
δ
(
α
1
-
β
+
1
)
.
If α
-
β
+
1
>
0 and α
1
-
β
+
1
<
0, that is, α
1
<
β
-
1
<
α, then L
1
′
(
u
1
)
<
0, L
1
′
(
u
2
)
>
0, and
(19)
P
2
(
1
,
u
2
)
=
δ
(
α
α
1
-
α
β
+
α
1
)
α
-
β
+
1
<
δ
(
-
α
+
α
1
)
α
-
β
+
1
<
0
.
Based on [6, Proposition 3], we can conclude that that O
1 is a stable node and O
2 is a saddle, respectively.We obtain the following system in a similar way by making transformations(
x
,
y
)
→
(
v
,
y
) with
(20)
v
=
x
y
,
y
≠
0
,
and d
t
→
y
d
t; then, we obtain
(21)
d
v
d
t
=
L
2
(
v
)
+
J
1
(
v
,
y
)
,
d
y
d
t
=
y
Q
2
(
v
,
1
)
+
J
2
(
v
,
y
)
,
where
(22)
L
2
(
v
)
=
-
L
(
v
,
1
)
=
-
δ
v
[
(
α
1
-
β
+
1
)
v
+
(
α
-
β
+
1
)
]
,
J
1
(
v
,
y
)
=
(
P
~
(
v
y
,
y
)
-
Q
~
(
v
y
,
y
)
v
)
y
2
,
J
2
(
v
,
y
)
=
Q
~
(
v
y
,
y
)
y
.
Thus, if (
α
-
β
+
1
)
(
α
1
-
β
+
1
)
>
0, system (21) has only one nonnegative equilibrium O
1
*(0,0). The system (21) has two equilibria O
1
*
(
0,0
) and O
2
*
(
v
2
,
0
) on the v-axis if (
α
-
β
+
1
)
(
α
1
-
β
+
1
)
<
0. However, O
2
* corresponds to O
2; thus, we do not need to study it again. O
1
* is a new equilibrium and it does not exist in the (
x
,
u
)-plane. Note that L
2
′
(
0
)
=
-
δ
(
α
+
1
-
β
) and Q
2
(
0,1
)
=
-
β
δ
<
0. Again, based on [6, Proposition 3], then comes the following conclusion.Lemma 4.
We say that(1)
ifα
-
β
+
1
<
0, O
1
*
(
0,0
) is a saddle;
(2)
ifα
-
β
+
1
>
0, O
1
*
(
0,0
) is a stable node.According to Lemmas3 and 4, it then follows from [6] that the following result holds for the origin E
0
(
0,0
) of system (7), which is equivalent to system (11).Theorem 5.
If(
α
-
β
+
1
)
(
α
1
-
β
+
1
)
≠
0, then there is at least one attracting parabolic sector in the neighborhood of (
0,0
) in the first quadrant, and no hyperbolic or elliptic sector can be found; that is, the complicated equilibrium E
0
(
0,0
) is stable.Remark 6.
If(
α
-
β
+
1
)
(
α
1
-
β
+
1
)
>
0, the neighborhood of E
0
(
0,0
) in the first quadrant has a parabolic (attracting) sector, but if (
α
-
β
+
1
)
(
α
1
-
β
+
1
)
<
0, there are two parabolic (attracting) sectors; in fact, we can regard those two sectors as a whole attracting one. If (
α
-
β
+
1
)
(
α
1
-
β
+
1
)
=
0, O
i, i
=
1,2 are also complicated equilibria and needed to apply the blow-up transformation again. Now we can see that the relative strength of mate-finding Allee effect δ greatly changes the structure of the neighborhood of E
0
(
0,0
) in our system (7), because P
2
(
1
,
u
1
)
=
-
δ and Q
2
(
0,1
)
=
-
β
δ always stay negative and O
1
(
0,0
) and O
1
*
(
0,0
) cannot be unstable; thus, the diversity of phase portraits around E
0
(
0,0
) of ratio-dependent prey-predator systems decreases. More importantly, the system becomes fragile and even has a higher risk of extinction due to mate-finding Allee effect.


## 2.3. Existence and Stability of Predator-Free Equilibria

In order to analyze the predator-free equilibria of the system (7), we set
(23)
Δ
1
=
(
g
-
δ
-
1
)
2
-
4
δ
.
We can notice that, for y
=
0, there exist two boundary equilibria E
10
=
(
x
1
,
0
) and E
20
=
(
x
2
,
0
) if 0
<
δ
<
g
-
1 and Δ
1
>
0; that is,(
B
1
)
0
<
δ
<
δ
1,whereδ
1
:
=
(
g
-
1
)
2. Let x
i, i
=
1,2 be the roots of the quadratic equation x
2
+
(
1
+
δ
-
g
)
x
+
δ
=
0; that is,
(24)
x
1
=
1
2
(
g
-
δ
-
1
-
(
g
-
δ
-
1
)
2
-
4
δ
)
,
x
2
=
1
2
(
g
-
δ
-
1
+
(
g
-
δ
-
1
)
2
-
4
δ
)
,
and they appear in pairs. Besides, there is no equilibrium on y-axis.The Jacobian matrix evaluated atE
i
0, i
=
1,2 is given by
(25)
J
i
0
=
(
x
i
(
x
i
+
δ
)
2
(
g
δ
-
(
x
i
+
δ
)
2
)
-
α
0
α
1
-
β
)
,
and one eigenvalue is λ
1
(
E
i
0
)
=
(
x
i
/
(
x
i
+
δ
)
2
)
(
g
δ
-
(
x
i
+
δ
)
2
), i
=
1,2, and the second eigenvalue is λ
2
(
E
i
0
)
=
α
1
-
β, i
=
1,2. Note that the sign of λ
1
(
E
i
0
)
,
i
=
1,2, depends on ψ
i
:
=
g
δ
-
δ
-
x
i, i
=
1,2. Using the expression of x
i, i
=
1,2, from (24), we can have
(26)
ψ
1
=
1
2
(
2
g
δ
-
δ
-
g
+
1
+
(
g
-
δ
-
1
)
2
-
4
δ
)
,
ψ
2
=
1
2
(
2
g
δ
-
δ
-
g
+
1
-
(
g
-
δ
-
1
)
2
-
4
δ
)
.
Under the restriction (
B
1
), it follows that
(27)
ψ
1
>
δ
(
g
-
δ
-
1
)
>
0
,
ψ
2
=
1
2
(
-
(
g
-
δ
)
2
+
1
-
Δ
1
)
<
-
Δ
1
2
<
0
.
Hence λ
1
(
E
10
) is always positive, while λ
1
(
E
20
) is negative. By using the facts above, we have the following theorem.Theorem 7.
Assume (B
1) holds. Then(1)
ifα
1
<
β, E
10 is a saddle and E
20 is a stable node (Figure 1(c)),
(2)
ifα
1
>
β, E
10 is an unstable node and E
20 is a saddle (Figure 2).Remark 8.
Ifδ
=
δ
1, system (7) has a unique predator-free boundary equilibrium point E
*
0
=
(
x
*
0
,
0
) (also known as the saddle-node equilibrium), where
(28)
x
*
0
=
1
2
(
g
-
δ
1
-
1
)
=
δ
1
.
Figure 1(b) shows the instantaneous equilibrium E
*
0 in the case of α
1
<
β and the behaviors of trajectories are divided by the stable manifold W
s
(
E
*
0
) of E
*
0 (the green curve).


## 2.4. Positive Equilibria and Their Stability and Bifurcation Results

### 2.4.1. Existence and Stability of Positive Equilibria

The interior (positive) equilibria can be evaluated by the intersections of the zero isoclines(29)
g
x
x
+
δ
-
(
1
+
x
)
=
α
y
x
+
y
,
α
1
x
x
+
y
=
β
,
in the first quadrant. The predator zero isocline, with the slope (
α
1
-
β
)
/
β, is a straight line passing through the origin, and it lies in the first quadrant if the following restriction holds(
C
1
)
α
1
>
β.Besides, we can define(30)
σ
:
=
1
+
α
-
α
β
α
1
.
Note that σ
>
1 under the restriction (
C
1
). And for x
2
+
(
δ
+
σ
-
g
)
x
+
δ
σ
=
0, we also set
(31)
Δ
2
=
(
g
-
δ
-
σ
)
2
-
4
δ
σ
.
If 0
<
δ
<
g
-
σ and Δ
2
>
0, that is,(
B
2
)
0
<
δ
<
δ
2,whereδ
2
:
=
(
g
-
σ
)
2, there are two positive equilibria E
1
*
=
(
x
1
*
,
y
1
*
) and E
2
*
=
(
x
2
*
,
y
2
*
), which also appear simultaneously, with the coordinates
(32)
x
1
*
=
1
2
(
g
-
δ
-
σ
-
(
g
-
δ
-
σ
)
2
-
4
δ
σ
)
,
x
2
*
=
1
2
(
g
-
δ
-
σ
+
(
g
-
δ
-
σ
)
2
-
4
δ
σ
)
,
(33)
y
i
*
=
α
1
-
β
β
x
i
*
,
i
=
1,2
.
Note that 0
<
δ
2
<
δ
1, it suggests that if (
B
2
) holds, then (
B
1
) holds; that is, if there exists any interior equilibrium in the first quadrant, then two predator-free equilibria will be detected. The Jacobian matrix J
i
*, i
=
1,2 evaluated at an interior equilibrium is given by
(34)
J
i
*
=
(
g
δ
x
i
*
(
x
i
*
+
δ
)
2
-
x
i
*
+
α
x
i
*
y
i
*
(
x
i
*
+
y
i
*
)
2
-
α
x
i
*
2
(
x
i
*
+
y
i
*
)
2
α
1
y
i
*
2
(
x
i
*
+
y
i
*
)
2
-
α
1
x
i
*
y
i
*
(
x
i
*
+
y
i
*
)
2
)
,
det
⁡
(
J
i
*
)
=
α
1
x
i
*
2
y
i
*
[
(
x
i
*
+
δ
)
2
-
g
δ
]
(
x
i
*
+
y
i
*
)
2
(
x
i
*
+
δ
)
2
.
Similarly, we have det
⁡
(
J
1
*
)
<
0 and det
⁡
(
J
2
*
)
>
0. Moreover,
(35)
tr
⁡
(
J
2
*
)
=
g
δ
x
2
*
(
x
2
*
+
δ
)
2
-
x
2
*
+
β
(
α
-
α
1
)
(
α
1
-
β
)
α
1
2
,
and it is interesting to find that tr
⁡
(
J
2
*
) is an increasing function of δ. In fact, we can have
(36)
d
x
2
d
δ
=
-
1
2
(
1
+
g
-
δ
+
σ
(
g
-
δ
+
σ
)
2
-
4
δ
σ
)
<
0
,
and we can define x
2
′
:
=
d
x
2
/
d
δ. Then, it follows that
(37)
dtr
(
J
2
*
)
d
δ
=
x
2
′
[
g
δ
(
x
2
+
δ
)
2
-
1
]
+
x
2
g
(
x
2
+
δ
)
3
[
x
2
-
δ
-
2
δ
x
2
′
]
>
0
,
since
(38)
-
δ
-
2
δ
x
2
′
=
-
δ
+
δ
(
1
+
g
-
δ
+
σ
(
g
-
δ
+
σ
)
2
-
4
δ
σ
)
=
δ
(
g
-
δ
+
σ
)
(
g
-
δ
+
σ
)
2
-
4
δ
σ
>
0
.
Moreover, we can have the following conclusion.Theorem 9.
Assume that(
B
2
) and (
C
1
) hold. Then one has the following.(1)
E
1
* is a saddle.
(2)
Ifα
⩽
α
1, E
2
* is stable.
(3)
Ifα
>
α
1, E
2
* can be stable or unstable depending on δ (Figure 2).


### 2.4.2. Hopf Bifurcation

Here assume that(
B
2
) and (
C
1
) hold, and we consider δ as the bifurcation parameter to discuss the Hopf bifurcation. It is easy to know that the Hopf bifurcation can occur only at E
2
* as E
1
* is always a saddle point. Then we try to show the existence of δ
=
δ
h such that tr
⁡
(
J
2
*
)
=
0 and det
⁡
(
J
2
*
)
>
0. From the expression of tr
⁡
(
J
2
*
) and some numerical attempts, we find that other parameters g, α, α
1, and β still have an important impact on the qualitative properties of the Hopf bifurcation though only δ is considered as the bifurcation parameter. The following assumption is used to guarantee the exsitence of the Hopf bifurcation with bifurcation parameter δ:(
C
2
)
α
>
α
1, (
1
-
(
β
/
α
1
)
)
(
α
-
β
+
(
α
β
/
α
1
)
)
<
g
-
1.Assumption (
C
2
) implies that tr
⁡
(
J
2
*
)
|
δ
=
0
<
0 and tr
⁡
(
J
2
*
)
|
δ
=
δ
2
>
0, and we have known that tr
⁡
(
J
2
*
) is a continuous increasing function of δ in [
0
,
+
∞
); thus, there exists δ
=
δ
h such that tr
⁡
(
J
2
*
)
=
0. Furthermore, the transversality condition (
d
/
d
δ
)
tr
⁡
(
J
2
*
)
|
δ
=
δ
h
>
0 is satisfied. Using the facts above, system (7) can undergo a Hopf bifurcation at E
2
* for δ
=
δ
h if (
C
2
) holds.Now we try to discuss the stability of the limit cycle of system (7) as a Hopf bifurcation occurs by computing the first Lyapunov coefficient l [26] at the E
2
*. Transformations u
1
=
x
-
x
2
* and u
2
=
y
-
y
2
* are used to translate the equilibrium E
2
*
(
x
2
*
,
y
2
*
) of the system (7) to the origin. Then we get
(39)
u
˙
1
=
a
10
u
1
+
a
01
u
2
+
a
20
u
1
2
+
a
11
u
1
u
2
+
a
02
u
2
2
+
a
03
u
1
3
+
a
21
u
1
2
u
2
+
a
12
u
1
u
2
2
+
a
03
u
2
3
+
R
1
(
u
1
,
u
2
)
,
u
˙
2
=
b
10
u
1
+
b
01
u
2
+
b
20
u
1
2
+
b
11
u
1
u
2
+
b
02
u
2
2
+
b
03
u
1
3
+
b
21
u
1
2
u
2
+
b
12
u
1
u
2
2
+
b
03
u
2
3
+
R
2
(
u
1
,
u
2
)
,
where
(40)
a
10
=
x
2
*
g
δ
(
x
2
*
+
δ
)
2
+
α
x
2
*
y
2
*
(
x
2
*
+
y
2
*
)
2
-
x
2
*
,
a
01
=
-
α
x
2
*
2
(
x
2
*
+
y
2
*
)
2
,
a
20
=
g
δ
2
(
x
2
*
+
δ
)
3
+
α
y
2
*
2
(
x
2
*
+
y
2
*
)
3
-
1
,
a
11
=
-
2
α
x
2
*
y
2
*
(
x
2
*
+
y
2
*
)
3
,
a
02
=
-
α
x
2
*
y
2
*
(
x
2
*
+
y
2
*
)
3
,
a
30
=
-
g
δ
2
(
x
2
*
+
δ
)
4
-
α
y
2
*
2
(
x
2
*
+
y
2
*
)
4
,
a
21
=
2
α
x
2
*
y
2
*
-
α
y
2
*
2
(
x
2
*
+
y
2
*
)
4
,
a
12
=
2
α
x
2
*
y
2
*
-
α
x
2
*
2
(
x
2
*
+
y
2
*
)
4
,
a
03
=
-
α
x
2
*
2
(
x
2
*
+
y
2
*
)
4
,
b
10
=
α
1
y
2
*
2
(
x
2
*
+
y
2
*
)
2
,
b
01
=
-
α
1
x
2
*
y
2
*
(
x
2
*
+
y
2
*
)
2
,
b
20
=
-
α
1
y
2
*
2
(
x
2
*
+
y
2
*
)
3
,
b
11
=
2
α
1
x
2
*
y
2
*
(
x
2
*
+
y
2
*
)
3
,
b
02
=
-
α
1
x
2
*
2
(
x
2
*
+
y
2
*
)
3
,
b
30
=
α
1
y
2
*
2
(
x
2
*
+
y
2
*
)
4
,
b
21
=
-
2
α
1
x
2
*
y
2
*
+
α
1
y
2
*
2
(
x
2
*
+
y
2
*
)
4
,
b
12
=
-
2
α
1
x
2
*
y
2
*
+
α
1
x
2
*
2
(
x
2
*
+
y
2
*
)
4
,
b
03
=
α
1
x
2
*
2
(
x
2
*
+
y
2
*
)
4
,
and R
i
(
x
,
y
), i
=
1,2 are power series in the powers of u
1
j
u
2
k satisfying j
+
k
⩾
4. The Lyapunov number l (as defined in [26]) is given by
(41)
l
=
-
3
π
2
a
01
Δ
3
/
2
×
(
[
a
10
b
10
(
a
11
2
+
a
11
b
02
+
a
02
b
11
)
ffffiffii
+
a
10
a
01
(
b
11
2
+
a
20
b
11
+
a
11
b
02
)
ffffiffii
+
b
10
2
(
a
11
a
02
+
2
a
02
b
02
)
-
2
a
10
b
10
(
b
02
2
-
a
20
a
02
)
fffifffii
-
2
a
10
a
01
(
a
20
2
-
b
20
b
02
)
-
a
01
2
(
2
a
20
b
20
+
b
11
b
20
)
fffffifii
+
(
a
01
b
10
-
2
a
10
2
)
(
b
11
b
02
-
a
11
a
20
)
]
ffiffifif
-
(
a
10
2
+
a
01
b
10
)
[

3
(
b
10
b
03
-
a
01
a
30
)
fffffffiffffffffffffffffffffiffff
+
2
a
10
(
a
21
+
b
12
)
ffffffffifffffffffffffffffffiffff
+
(
b
10
a
12
-
a
01
b
21
)
]
a
11
2
)
,
where Δ
=
a
10
b
01
-
a
01
b
10. Since the expression Lyapunov number for l is very complex, we fail to discuss the sign of l precisely though a subcritical Hopf bifurcation has been found numerically with the Lyapunov number l
=
5.545875
>
0 (Figure 3). Obviously, the assumption (
C
2
) is satisfied in Figure 3. We need to require δ
<
δ
1 (green line) and δ
<
δ
2 (blue line) to guarantee the existence of E
2
*; that is, (
B
1
) and (
B
2
) hold (Figure 3(a)). tr
⁡
(
J
2
*
) is a continuous increasing function of δ. If 0
<
δ
<
δ
h, system (7) has a stable E
2
* (Figure 2(a)). As δ increases, a homoclinic loop is created by joining the stable and unstable manifolds of the saddle E
1
* at some δ
~ (δ
~
<
δ
h) with a stable E
2
* inside (Figure 2(b)). Then as δ gets larger, an unstable limit cycle appears and coexists with a stable E
2
* (Figure 4). However, the limit cycle begins to shrink with the increasing δ. The shrinking unstable limit cycle (Figures 3(b) and 3(c)) disappears as δ passes through δ
h
(
δ
h is also the root of equation tr
⁡
(
J
2
*
)
=
0, δ
h
=
0.009724
). Then δ
h
<
δ
<
δ
2, E
2
* becomes unstable (Figure 2(c)).Bifurcation structure ofE
2
* for Hopf bifurcation and phase portraits of system (7) for g
=
1.44, α
=
0.8, α
1
=
0.4, and β
=
0.3.
(a)
Bifurcation structure
(b)
δ
=
0.009724
(
δ
=
δ
h
), the unstable limit cycle when a subcritical Hopf bifurcation occurs
(c)
δ
=
0.009724
(
δ
=
δ
h
), the trajectory around the unstable limit cycle when a subcritical Hopf bifurcation occursAn unstable limit cycle coexists with a stableE
2
*, for g
=
1.44, δ
=
0.0095, α
=
0.84, α
1
=
0.4, and β
=
0.3.
(a)
Inside the limit cycle
(b)
Outside the limit cycle


### 2.4.3. Saddle-Node Bifurcation

We know that wheng
*, δ
*, α
*, α
1
*, and β
* satisfy (
C
1
) and the following equation
(42)
δ
-
(
g
-
σ
)
2
=
0
,
where σ
=
1
+
α
-
α
β
/
α
1, then there is only one interior equilibrium E
* of system (7), whose coordinates are given by
(43)
x
*
=
1
2
(
g
*
-
δ
*
-
σ
*
)
=
δ
*
σ
*
,
y
*
=
α
1
*
-
β
*
2
β
*
(
g
*
-
δ
*
-
σ
*
)
=
α
1
*
-
β
*
β
*
δ
*
σ
*
.
And the Jacobian matrix evaluated at E
* is
(44)
J
*
=
(
α
*
x
*
y
*
(
x
*
+
y
*
)
2
-
α
*
x
*
2
(
x
*
+
y
*
)
2
α
1
*
y
*
2
(
x
*
+
y
*
)
2
-
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
)
.
Next we choose the relative strength of mate-finding Allee effect δ as the bifurcation parameter; then, we have the following theorem.Theorem 10.
Assume (C
1) holds. For δ
=
δ
* and α
1
*
≠
α
*, a saddle-node bifurcation occurs at the unique positive equilibrium E
* of system (7).Proof.
First, we havedet
⁡
(
J
*
)
=
0; therefore, J
* has an eigenvalue λ
*
=
0, and if α
1
*
≠
α
*, then λ
*
=
0 is simple. Let W
1 and W
2 be the eigenvectors corresponding to the λ
*
=
0 for J
* and J
*
T, respectively. Then we obtain
(45)
W
1
=
(
1
α
1
*
-
β
*
β
*
)
:
=
(
w
11
w
12
)
,
W
2
=
(
1
-
α
*
β
*
α
1
*
(
α
1
*
-
β
*
)
)
.
From the expressions for W
1 and W
2, we get
(46)
W
2
T
F
δ
(
(
x
*
,
y
*
)
;
δ
*
)
=
-
g
*
x
*
2
(
x
*
+
δ
*
)
2
<
0
(
≠
0
)
,
W
2
T
D
2
F
(
W
1
,
W
1
)
=
-
2
x
*
x
*
+
δ
*
<
0
(
≠
0
)
,
where
(47)
F
δ
=
∂
F
∂
δ
,
D
2
F
(
W
1
,
W
1
)
=
∂
2
F
∂
x
2
w
11
w
11
+
∂
2
F
∂
x
∂
y
w
11
w
12
+
∂
2
F
∂
y
∂
x
w
12
w
11
+
∂
2
F
∂
y
2
w
12
w
12
.
Thus by Sotomayor’s theorem [26], system (7) undergoes a saddle-node bifurcation at E
* as δ passes δ
=
δ
* if α
1
*
≠
α
*.A numerical example has been carried out to demonstrate a saddle-node bifurcation of system (7) (Figure 5). The green curves represent the smaller unstable equilibria. The stability of the larger equilibria changes from stable (blue curves) to unstable (black curves) as δ passes through δ
h
(
δ
h
=
0.009724
).A saddle-node bifurcation of system (7) at δ
=
0.010932 (mark as SN), for g
=
1.44, α
=
0.8, α
1
=
0.4, and β
=
0.3. A Hopf bifurcation also occurs.
(a)
Prey population scenario
(b)
Predator population scenarioBifurcation structure of system (7) in (
δ
,
α
1
) with a Bogdanov-Takens bifurcation and phase portraits for g
=
1.44, δ
=
0.022859, α
=
0.4, α
1
=
0.4, and β
=
0.3.
(a)
Bifurcation structure
(b)
A cuspE
*


### 2.4.4. Bogdanov-Takens Bifurcation

Let us focus onE
* again. After knowing that a saddle-node bifurcation can occur at E
* if δ
=
δ
* and α
1
*
≠
α
*, we go on to consider the case of δ
=
δ
2, det
⁡
(
J
*
)
=
0, and tr
⁡
(
J
*
)
=
0. If det
⁡
(
J
*
)
=
0 and tr
⁡
(
J
*
)
=
0, then the Jacobian matrix of E
* has λ
*
1
=
λ
*
2
=
0. When det
⁡
(
J
*
)
=
0, it is easily shown that if α
1
*
=
α
* then tr
⁡
(
J
*
)
=
0, and hence we choose the relative mate-finding strength δ and the predator growing ability α
1 as two bifurcation parameters. Then the following statement holds.Theorem 11.
Let (C
1) hold. The system (7) undergoes a Bogdanov-Takens bifurcation around the equilibrium point E
* when δ
=
δ
* and α
1
*
=
α
*.Proof.
Following the ideas in [27, 28], we consider the neighborhood of (
δ
*
,
α
1
*
); that is, δ
=
δ
*
+
r
1 and α
1
=
α
1
*
+
r
2, where r
i, i
=
1,2, are sufficient small, and system (7) becomes
(48)
d
x
d
t
=
g
*
x
2
δ
*
+
r
1
+
x
-
d
*
x
(
1
+
x
)
-
α
1
*
x
y
x
+
y
,
d
y
d
t
=
(
α
1
*
+
r
2
)
x
y
x
+
y
-
β
*
y
.
We make z
1
=
x
-
x
*, z
2
=
y
-
y
* for translating E
*
=
(
x
*
,
y
*
) to the origin (
0,0
), and we obtain
(49)
z
˙
1
=
[
-
g
*
x
*
2
(
x
*
+
δ
*
)
2
r
1
+
ξ
1
(
r
1
)
]
+
[
-
2
g
*
δ
*
x
*
(
x
*
+
δ
*
)
3
r
1
+
ξ
2
(
r
1
)
+
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
]
z
1
-
α
1
*
x
*
2
(
x
*
+
y
*
)
2
z
2
+
[
g
*
δ
*
(
2
x
*
-
δ
*
)
(
x
*
+
δ
*
)
4
r
1
+
ξ
3
(
r
1
)
-
α
1
*
x
*
2
(
x
*
+
y
*
)
2
z
2
ii
+
i
+
α
1
*
y
*
2
(
x
*
+
y
*
)
3
-
x
*
x
*
+
δ
*
]
z
1
2
-
2
α
1
*
x
*
y
*
(
x
*
+
y
*
)
3
z
1
z
2
+
α
1
*
x
*
2
(
x
*
+
y
*
)
3
z
2
2
+
T
1
(
z
1
,
z
2
)
:
=
Q
1
(
z
1
,
z
2
,
r
1
,
r
2
)
,
z
˙
2
=
x
*
y
*
x
*
+
y
*
r
2
+
[
α
1
*
y
*
2
(
x
*
+
y
*
)
2
+
y
*
2
r
2
(
x
*
+
y
*
)
2
]
z
1
+
[
-
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
+
x
*
2
r
2
(
x
*
+
y
*
)
2
]
z
2
-
[
α
1
*
y
*
2
(
x
*
+
y
*
)
3
+
y
*
2
(
x
*
+
y
*
)
3
r
2
]
z
1
2
+
[
2
α
1
*
x
*
y
*
(
x
*
+
y
*
)
3
+
2
x
*
y
*
(
x
*
+
y
*
)
3
r
2
]
z
1
z
2
-
[
α
1
*
x
*
2
(
x
*
+
y
*
)
3
+
x
*
2
(
x
*
+
y
*
)
3
r
2
]
z
2
2
+
T
2
(
z
1
,
z
2
)
:
=
Q
2
(
z
1
,
z
2
,
r
1
,
r
2
)
,
where ξ
i, i
=
1,2
,
3 are polynomial functions of r
1 with r
1
p satisfying p
⩾
2, and T
1
(
z
1
,
z
2
) and T
2
(
z
1
,
z
2
) are C
∞ functions z
1
j
z
2
k satisfying j
+
k
⩾
3.
Then we apply the affine transformation(50)
v
1
=
z
1
,
v
2
=
[
-
2
g
*
δ
*
x
*
(
x
*
+
δ
*
)
3
r
1
+
ξ
2
(
r
1
)
+
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
]
z
1
-
α
1
*
x
*
2
(
x
*
+
y
*
)
2
z
2
,
to obtain
(51)
v
˙
1
=
v
2
+
[
-
x
*
x
*
+
δ
*
+
ξ
~
1
(
r
1
)
]
v
1
2
+
ξ
~
2
(
r
1
)
v
1
v
2
+
x
*
+
y
*
α
1
*
x
*
2
v
2
2
+
T
~
1
(
v
1
,
v
2
,
r
1
,
r
2
)
,
v
˙
2
=
η
~
1
(
r
1
,
r
2
)
v
1
+
η
~
2
(
r
1
,
r
2
)
v
2
+
[
-
α
1
*
x
*
2
y
*
(
x
*
+
δ
*
)
(
x
*
+
y
*
)
2
+
η
~
3
(
r
1
,
r
2
)
]
v
1
2
+
η
~
4
(
r
1
,
r
2
)
v
1
v
2
+
[
1
x
*
+
η
~
5
(
r
1
,
r
2
)
]
v
2
2
+
T
~
2
(
v
1
,
v
2
,
r
1
,
r
2
)
,
where ξ
~
i
(
r
1
)
,
i
=
1,2
,
3 and η
~
i
(
r
1
,
r
2
)
,
i
=
1,2
,
3,4
,
5 are polynomial functions with r
1
p satisfying p
⩾
1 and r
1
p
r
2
q satisfying p
+
q
⩾
1, respectively. T
~
i, i
=
1,2, are C
∞ functions in their variables with v
1
j
v
2
k satisfying j
+
k
⩾
3.
Furthermore, ifr
1
=
r
2
=
0, we can see that(BT0)
(52)
J
(
*
,
0,0
)
=
(
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
-
α
1
*
x
*
2
(
x
*
+
y
*
)
2
α
1
*
y
*
2
(
x
*
+
y
*
)
2
-
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
)
≠
θ
2
×
2
,
(BT1)
[
-
(
x
*
/
(
x
*
+
δ
*
)
)
+
ξ
~
1
(
r
1
)
+
η
~
4
(
r
1
,
r
2
)
]
|
r
1
=
0
,
r
2
=
0
=
-
x
*
/
(
x
*
+
δ
*
)
<
0,
(BT2)
[
-
(
α
1
*
x
*
2
y
*
/
(
(
x
*
+
δ
*
)
(
x
*
+
y
*
)
2
)
)
+
η
~
3
(
r
1
,
r
2
)
]
|
r
1
=
0
,
r
2
=
0
=
-
α
1
*
x
*
2
y
*
/
(
(
x
*
+
δ
*
)
(
x
*
+
y
*
)
2
)
<
0,
(BT3)
the map(53)
(
(
z
1
z
2
)
,
(
r
1
r
2
)
)
⟶
(
(
Q
1
Q
2
)
,
tr
⁡
(
∂
Q
i
∂
z
i
)
,
det
⁡
⁡
(
∂
Q
i
∂
z
i
)
)
,
i
=
1,2
is regular at(54)
(
z
1
z
2
)
=
(
0
0
)
,
(
r
1
r
2
)
=
(
0
0
)
.According to [27, 28], the system (7) undergoes a Bogdanov-Takens bifurcation around the equilibrium point E
* when δ
=
δ
* and α
1
*
=
α
*. This completes the proof.Remark 12.
Furthermore, the quantityω [27, 28] which determines the structure of the bifurcation is given by
(55)
ω
=
sign
⁡
[
α
1
*
x
*
3
y
*
(
x
*
+
δ
*
)
2
(
x
*
+
y
*
)
2
]
=
1
.The bifurcation diagram of system (7) in (
δ
,
α
1
) is presented in Figure 6(a). There are no positive equilibria of system (7) in Domain 1. In the domain bounded by Limit point cycle (green) and Hopf curve (blue), there exist a saddle E
1
* and an unstable E
2
*. In the domain bounded by Hopf curve (blue) and Homoclinic curve (dot and magenta), there is a saddle E
1
* and a stable E
2
* surrounded by an unstable limit cycle. In Domain 2, system (7) has a saddle E
1
* and a stable E
2
*. At BT point, E
* is the unique positive equilibrium of system (7) in the first quadrant, and the local dynamics near E
* is characterized by two trajectories (green curves, Figure 6(b)), which meet at E
*, but with opposite time orientation, and also behave like a cusp-like configuration.


## 2.4.1. Existence and Stability of Positive Equilibria

The interior (positive) equilibria can be evaluated by the intersections of the zero isoclines(29)
g
x
x
+
δ
-
(
1
+
x
)
=
α
y
x
+
y
,
α
1
x
x
+
y
=
β
,
in the first quadrant. The predator zero isocline, with the slope (
α
1
-
β
)
/
β, is a straight line passing through the origin, and it lies in the first quadrant if the following restriction holds(
C
1
)
α
1
>
β.Besides, we can define(30)
σ
:
=
1
+
α
-
α
β
α
1
.
Note that σ
>
1 under the restriction (
C
1
). And for x
2
+
(
δ
+
σ
-
g
)
x
+
δ
σ
=
0, we also set
(31)
Δ
2
=
(
g
-
δ
-
σ
)
2
-
4
δ
σ
.
If 0
<
δ
<
g
-
σ and Δ
2
>
0, that is,(
B
2
)
0
<
δ
<
δ
2,whereδ
2
:
=
(
g
-
σ
)
2, there are two positive equilibria E
1
*
=
(
x
1
*
,
y
1
*
) and E
2
*
=
(
x
2
*
,
y
2
*
), which also appear simultaneously, with the coordinates
(32)
x
1
*
=
1
2
(
g
-
δ
-
σ
-
(
g
-
δ
-
σ
)
2
-
4
δ
σ
)
,
x
2
*
=
1
2
(
g
-
δ
-
σ
+
(
g
-
δ
-
σ
)
2
-
4
δ
σ
)
,
(33)
y
i
*
=
α
1
-
β
β
x
i
*
,
i
=
1,2
.
Note that 0
<
δ
2
<
δ
1, it suggests that if (
B
2
) holds, then (
B
1
) holds; that is, if there exists any interior equilibrium in the first quadrant, then two predator-free equilibria will be detected. The Jacobian matrix J
i
*, i
=
1,2 evaluated at an interior equilibrium is given by
(34)
J
i
*
=
(
g
δ
x
i
*
(
x
i
*
+
δ
)
2
-
x
i
*
+
α
x
i
*
y
i
*
(
x
i
*
+
y
i
*
)
2
-
α
x
i
*
2
(
x
i
*
+
y
i
*
)
2
α
1
y
i
*
2
(
x
i
*
+
y
i
*
)
2
-
α
1
x
i
*
y
i
*
(
x
i
*
+
y
i
*
)
2
)
,
det
⁡
(
J
i
*
)
=
α
1
x
i
*
2
y
i
*
[
(
x
i
*
+
δ
)
2
-
g
δ
]
(
x
i
*
+
y
i
*
)
2
(
x
i
*
+
δ
)
2
.
Similarly, we have det
⁡
(
J
1
*
)
<
0 and det
⁡
(
J
2
*
)
>
0. Moreover,
(35)
tr
⁡
(
J
2
*
)
=
g
δ
x
2
*
(
x
2
*
+
δ
)
2
-
x
2
*
+
β
(
α
-
α
1
)
(
α
1
-
β
)
α
1
2
,
and it is interesting to find that tr
⁡
(
J
2
*
) is an increasing function of δ. In fact, we can have
(36)
d
x
2
d
δ
=
-
1
2
(
1
+
g
-
δ
+
σ
(
g
-
δ
+
σ
)
2
-
4
δ
σ
)
<
0
,
and we can define x
2
′
:
=
d
x
2
/
d
δ. Then, it follows that
(37)
dtr
(
J
2
*
)
d
δ
=
x
2
′
[
g
δ
(
x
2
+
δ
)
2
-
1
]
+
x
2
g
(
x
2
+
δ
)
3
[
x
2
-
δ
-
2
δ
x
2
′
]
>
0
,
since
(38)
-
δ
-
2
δ
x
2
′
=
-
δ
+
δ
(
1
+
g
-
δ
+
σ
(
g
-
δ
+
σ
)
2
-
4
δ
σ
)
=
δ
(
g
-
δ
+
σ
)
(
g
-
δ
+
σ
)
2
-
4
δ
σ
>
0
.
Moreover, we can have the following conclusion.Theorem 9.
Assume that(
B
2
) and (
C
1
) hold. Then one has the following.(1)
E
1
* is a saddle.
(2)
Ifα
⩽
α
1, E
2
* is stable.
(3)
Ifα
>
α
1, E
2
* can be stable or unstable depending on δ (Figure 2).


## 2.4.2. Hopf Bifurcation

Here assume that(
B
2
) and (
C
1
) hold, and we consider δ as the bifurcation parameter to discuss the Hopf bifurcation. It is easy to know that the Hopf bifurcation can occur only at E
2
* as E
1
* is always a saddle point. Then we try to show the existence of δ
=
δ
h such that tr
⁡
(
J
2
*
)
=
0 and det
⁡
(
J
2
*
)
>
0. From the expression of tr
⁡
(
J
2
*
) and some numerical attempts, we find that other parameters g, α, α
1, and β still have an important impact on the qualitative properties of the Hopf bifurcation though only δ is considered as the bifurcation parameter. The following assumption is used to guarantee the exsitence of the Hopf bifurcation with bifurcation parameter δ:(
C
2
)
α
>
α
1, (
1
-
(
β
/
α
1
)
)
(
α
-
β
+
(
α
β
/
α
1
)
)
<
g
-
1.Assumption (
C
2
) implies that tr
⁡
(
J
2
*
)
|
δ
=
0
<
0 and tr
⁡
(
J
2
*
)
|
δ
=
δ
2
>
0, and we have known that tr
⁡
(
J
2
*
) is a continuous increasing function of δ in [
0
,
+
∞
); thus, there exists δ
=
δ
h such that tr
⁡
(
J
2
*
)
=
0. Furthermore, the transversality condition (
d
/
d
δ
)
tr
⁡
(
J
2
*
)
|
δ
=
δ
h
>
0 is satisfied. Using the facts above, system (7) can undergo a Hopf bifurcation at E
2
* for δ
=
δ
h if (
C
2
) holds.Now we try to discuss the stability of the limit cycle of system (7) as a Hopf bifurcation occurs by computing the first Lyapunov coefficient l [26] at the E
2
*. Transformations u
1
=
x
-
x
2
* and u
2
=
y
-
y
2
* are used to translate the equilibrium E
2
*
(
x
2
*
,
y
2
*
) of the system (7) to the origin. Then we get
(39)
u
˙
1
=
a
10
u
1
+
a
01
u
2
+
a
20
u
1
2
+
a
11
u
1
u
2
+
a
02
u
2
2
+
a
03
u
1
3
+
a
21
u
1
2
u
2
+
a
12
u
1
u
2
2
+
a
03
u
2
3
+
R
1
(
u
1
,
u
2
)
,
u
˙
2
=
b
10
u
1
+
b
01
u
2
+
b
20
u
1
2
+
b
11
u
1
u
2
+
b
02
u
2
2
+
b
03
u
1
3
+
b
21
u
1
2
u
2
+
b
12
u
1
u
2
2
+
b
03
u
2
3
+
R
2
(
u
1
,
u
2
)
,
where
(40)
a
10
=
x
2
*
g
δ
(
x
2
*
+
δ
)
2
+
α
x
2
*
y
2
*
(
x
2
*
+
y
2
*
)
2
-
x
2
*
,
a
01
=
-
α
x
2
*
2
(
x
2
*
+
y
2
*
)
2
,
a
20
=
g
δ
2
(
x
2
*
+
δ
)
3
+
α
y
2
*
2
(
x
2
*
+
y
2
*
)
3
-
1
,
a
11
=
-
2
α
x
2
*
y
2
*
(
x
2
*
+
y
2
*
)
3
,
a
02
=
-
α
x
2
*
y
2
*
(
x
2
*
+
y
2
*
)
3
,
a
30
=
-
g
δ
2
(
x
2
*
+
δ
)
4
-
α
y
2
*
2
(
x
2
*
+
y
2
*
)
4
,
a
21
=
2
α
x
2
*
y
2
*
-
α
y
2
*
2
(
x
2
*
+
y
2
*
)
4
,
a
12
=
2
α
x
2
*
y
2
*
-
α
x
2
*
2
(
x
2
*
+
y
2
*
)
4
,
a
03
=
-
α
x
2
*
2
(
x
2
*
+
y
2
*
)
4
,
b
10
=
α
1
y
2
*
2
(
x
2
*
+
y
2
*
)
2
,
b
01
=
-
α
1
x
2
*
y
2
*
(
x
2
*
+
y
2
*
)
2
,
b
20
=
-
α
1
y
2
*
2
(
x
2
*
+
y
2
*
)
3
,
b
11
=
2
α
1
x
2
*
y
2
*
(
x
2
*
+
y
2
*
)
3
,
b
02
=
-
α
1
x
2
*
2
(
x
2
*
+
y
2
*
)
3
,
b
30
=
α
1
y
2
*
2
(
x
2
*
+
y
2
*
)
4
,
b
21
=
-
2
α
1
x
2
*
y
2
*
+
α
1
y
2
*
2
(
x
2
*
+
y
2
*
)
4
,
b
12
=
-
2
α
1
x
2
*
y
2
*
+
α
1
x
2
*
2
(
x
2
*
+
y
2
*
)
4
,
b
03
=
α
1
x
2
*
2
(
x
2
*
+
y
2
*
)
4
,
and R
i
(
x
,
y
), i
=
1,2 are power series in the powers of u
1
j
u
2
k satisfying j
+
k
⩾
4. The Lyapunov number l (as defined in [26]) is given by
(41)
l
=
-
3
π
2
a
01
Δ
3
/
2
×
(
[
a
10
b
10
(
a
11
2
+
a
11
b
02
+
a
02
b
11
)
ffffiffii
+
a
10
a
01
(
b
11
2
+
a
20
b
11
+
a
11
b
02
)
ffffiffii
+
b
10
2
(
a
11
a
02
+
2
a
02
b
02
)
-
2
a
10
b
10
(
b
02
2
-
a
20
a
02
)
fffifffii
-
2
a
10
a
01
(
a
20
2
-
b
20
b
02
)
-
a
01
2
(
2
a
20
b
20
+
b
11
b
20
)
fffffifii
+
(
a
01
b
10
-
2
a
10
2
)
(
b
11
b
02
-
a
11
a
20
)
]
ffiffifif
-
(
a
10
2
+
a
01
b
10
)
[

3
(
b
10
b
03
-
a
01
a
30
)
fffffffiffffffffffffffffffffiffff
+
2
a
10
(
a
21
+
b
12
)
ffffffffifffffffffffffffffffiffff
+
(
b
10
a
12
-
a
01
b
21
)
]
a
11
2
)
,
where Δ
=
a
10
b
01
-
a
01
b
10. Since the expression Lyapunov number for l is very complex, we fail to discuss the sign of l precisely though a subcritical Hopf bifurcation has been found numerically with the Lyapunov number l
=
5.545875
>
0 (Figure 3). Obviously, the assumption (
C
2
) is satisfied in Figure 3. We need to require δ
<
δ
1 (green line) and δ
<
δ
2 (blue line) to guarantee the existence of E
2
*; that is, (
B
1
) and (
B
2
) hold (Figure 3(a)). tr
⁡
(
J
2
*
) is a continuous increasing function of δ. If 0
<
δ
<
δ
h, system (7) has a stable E
2
* (Figure 2(a)). As δ increases, a homoclinic loop is created by joining the stable and unstable manifolds of the saddle E
1
* at some δ
~ (δ
~
<
δ
h) with a stable E
2
* inside (Figure 2(b)). Then as δ gets larger, an unstable limit cycle appears and coexists with a stable E
2
* (Figure 4). However, the limit cycle begins to shrink with the increasing δ. The shrinking unstable limit cycle (Figures 3(b) and 3(c)) disappears as δ passes through δ
h
(
δ
h is also the root of equation tr
⁡
(
J
2
*
)
=
0, δ
h
=
0.009724
). Then δ
h
<
δ
<
δ
2, E
2
* becomes unstable (Figure 2(c)).Bifurcation structure ofE
2
* for Hopf bifurcation and phase portraits of system (7) for g
=
1.44, α
=
0.8, α
1
=
0.4, and β
=
0.3.
(a)
Bifurcation structure
(b)
δ
=
0.009724
(
δ
=
δ
h
), the unstable limit cycle when a subcritical Hopf bifurcation occurs
(c)
δ
=
0.009724
(
δ
=
δ
h
), the trajectory around the unstable limit cycle when a subcritical Hopf bifurcation occursAn unstable limit cycle coexists with a stableE
2
*, for g
=
1.44, δ
=
0.0095, α
=
0.84, α
1
=
0.4, and β
=
0.3.
(a)
Inside the limit cycle
(b)
Outside the limit cycle


## 2.4.3. Saddle-Node Bifurcation

We know that wheng
*, δ
*, α
*, α
1
*, and β
* satisfy (
C
1
) and the following equation
(42)
δ
-
(
g
-
σ
)
2
=
0
,
where σ
=
1
+
α
-
α
β
/
α
1, then there is only one interior equilibrium E
* of system (7), whose coordinates are given by
(43)
x
*
=
1
2
(
g
*
-
δ
*
-
σ
*
)
=
δ
*
σ
*
,
y
*
=
α
1
*
-
β
*
2
β
*
(
g
*
-
δ
*
-
σ
*
)
=
α
1
*
-
β
*
β
*
δ
*
σ
*
.
And the Jacobian matrix evaluated at E
* is
(44)
J
*
=
(
α
*
x
*
y
*
(
x
*
+
y
*
)
2
-
α
*
x
*
2
(
x
*
+
y
*
)
2
α
1
*
y
*
2
(
x
*
+
y
*
)
2
-
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
)
.
Next we choose the relative strength of mate-finding Allee effect δ as the bifurcation parameter; then, we have the following theorem.Theorem 10.
Assume (C
1) holds. For δ
=
δ
* and α
1
*
≠
α
*, a saddle-node bifurcation occurs at the unique positive equilibrium E
* of system (7).Proof.
First, we havedet
⁡
(
J
*
)
=
0; therefore, J
* has an eigenvalue λ
*
=
0, and if α
1
*
≠
α
*, then λ
*
=
0 is simple. Let W
1 and W
2 be the eigenvectors corresponding to the λ
*
=
0 for J
* and J
*
T, respectively. Then we obtain
(45)
W
1
=
(
1
α
1
*
-
β
*
β
*
)
:
=
(
w
11
w
12
)
,
W
2
=
(
1
-
α
*
β
*
α
1
*
(
α
1
*
-
β
*
)
)
.
From the expressions for W
1 and W
2, we get
(46)
W
2
T
F
δ
(
(
x
*
,
y
*
)
;
δ
*
)
=
-
g
*
x
*
2
(
x
*
+
δ
*
)
2
<
0
(
≠
0
)
,
W
2
T
D
2
F
(
W
1
,
W
1
)
=
-
2
x
*
x
*
+
δ
*
<
0
(
≠
0
)
,
where
(47)
F
δ
=
∂
F
∂
δ
,
D
2
F
(
W
1
,
W
1
)
=
∂
2
F
∂
x
2
w
11
w
11
+
∂
2
F
∂
x
∂
y
w
11
w
12
+
∂
2
F
∂
y
∂
x
w
12
w
11
+
∂
2
F
∂
y
2
w
12
w
12
.
Thus by Sotomayor’s theorem [26], system (7) undergoes a saddle-node bifurcation at E
* as δ passes δ
=
δ
* if α
1
*
≠
α
*.A numerical example has been carried out to demonstrate a saddle-node bifurcation of system (7) (Figure 5). The green curves represent the smaller unstable equilibria. The stability of the larger equilibria changes from stable (blue curves) to unstable (black curves) as δ passes through δ
h
(
δ
h
=
0.009724
).A saddle-node bifurcation of system (7) at δ
=
0.010932 (mark as SN), for g
=
1.44, α
=
0.8, α
1
=
0.4, and β
=
0.3. A Hopf bifurcation also occurs.
(a)
Prey population scenario
(b)
Predator population scenarioBifurcation structure of system (7) in (
δ
,
α
1
) with a Bogdanov-Takens bifurcation and phase portraits for g
=
1.44, δ
=
0.022859, α
=
0.4, α
1
=
0.4, and β
=
0.3.
(a)
Bifurcation structure
(b)
A cuspE
*


## 2.4.4. Bogdanov-Takens Bifurcation

Let us focus onE
* again. After knowing that a saddle-node bifurcation can occur at E
* if δ
=
δ
* and α
1
*
≠
α
*, we go on to consider the case of δ
=
δ
2, det
⁡
(
J
*
)
=
0, and tr
⁡
(
J
*
)
=
0. If det
⁡
(
J
*
)
=
0 and tr
⁡
(
J
*
)
=
0, then the Jacobian matrix of E
* has λ
*
1
=
λ
*
2
=
0. When det
⁡
(
J
*
)
=
0, it is easily shown that if α
1
*
=
α
* then tr
⁡
(
J
*
)
=
0, and hence we choose the relative mate-finding strength δ and the predator growing ability α
1 as two bifurcation parameters. Then the following statement holds.Theorem 11.
Let (C
1) hold. The system (7) undergoes a Bogdanov-Takens bifurcation around the equilibrium point E
* when δ
=
δ
* and α
1
*
=
α
*.Proof.
Following the ideas in [27, 28], we consider the neighborhood of (
δ
*
,
α
1
*
); that is, δ
=
δ
*
+
r
1 and α
1
=
α
1
*
+
r
2, where r
i, i
=
1,2, are sufficient small, and system (7) becomes
(48)
d
x
d
t
=
g
*
x
2
δ
*
+
r
1
+
x
-
d
*
x
(
1
+
x
)
-
α
1
*
x
y
x
+
y
,
d
y
d
t
=
(
α
1
*
+
r
2
)
x
y
x
+
y
-
β
*
y
.
We make z
1
=
x
-
x
*, z
2
=
y
-
y
* for translating E
*
=
(
x
*
,
y
*
) to the origin (
0,0
), and we obtain
(49)
z
˙
1
=
[
-
g
*
x
*
2
(
x
*
+
δ
*
)
2
r
1
+
ξ
1
(
r
1
)
]
+
[
-
2
g
*
δ
*
x
*
(
x
*
+
δ
*
)
3
r
1
+
ξ
2
(
r
1
)
+
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
]
z
1
-
α
1
*
x
*
2
(
x
*
+
y
*
)
2
z
2
+
[
g
*
δ
*
(
2
x
*
-
δ
*
)
(
x
*
+
δ
*
)
4
r
1
+
ξ
3
(
r
1
)
-
α
1
*
x
*
2
(
x
*
+
y
*
)
2
z
2
ii
+
i
+
α
1
*
y
*
2
(
x
*
+
y
*
)
3
-
x
*
x
*
+
δ
*
]
z
1
2
-
2
α
1
*
x
*
y
*
(
x
*
+
y
*
)
3
z
1
z
2
+
α
1
*
x
*
2
(
x
*
+
y
*
)
3
z
2
2
+
T
1
(
z
1
,
z
2
)
:
=
Q
1
(
z
1
,
z
2
,
r
1
,
r
2
)
,
z
˙
2
=
x
*
y
*
x
*
+
y
*
r
2
+
[
α
1
*
y
*
2
(
x
*
+
y
*
)
2
+
y
*
2
r
2
(
x
*
+
y
*
)
2
]
z
1
+
[
-
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
+
x
*
2
r
2
(
x
*
+
y
*
)
2
]
z
2
-
[
α
1
*
y
*
2
(
x
*
+
y
*
)
3
+
y
*
2
(
x
*
+
y
*
)
3
r
2
]
z
1
2
+
[
2
α
1
*
x
*
y
*
(
x
*
+
y
*
)
3
+
2
x
*
y
*
(
x
*
+
y
*
)
3
r
2
]
z
1
z
2
-
[
α
1
*
x
*
2
(
x
*
+
y
*
)
3
+
x
*
2
(
x
*
+
y
*
)
3
r
2
]
z
2
2
+
T
2
(
z
1
,
z
2
)
:
=
Q
2
(
z
1
,
z
2
,
r
1
,
r
2
)
,
where ξ
i, i
=
1,2
,
3 are polynomial functions of r
1 with r
1
p satisfying p
⩾
2, and T
1
(
z
1
,
z
2
) and T
2
(
z
1
,
z
2
) are C
∞ functions z
1
j
z
2
k satisfying j
+
k
⩾
3.
Then we apply the affine transformation(50)
v
1
=
z
1
,
v
2
=
[
-
2
g
*
δ
*
x
*
(
x
*
+
δ
*
)
3
r
1
+
ξ
2
(
r
1
)
+
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
]
z
1
-
α
1
*
x
*
2
(
x
*
+
y
*
)
2
z
2
,
to obtain
(51)
v
˙
1
=
v
2
+
[
-
x
*
x
*
+
δ
*
+
ξ
~
1
(
r
1
)
]
v
1
2
+
ξ
~
2
(
r
1
)
v
1
v
2
+
x
*
+
y
*
α
1
*
x
*
2
v
2
2
+
T
~
1
(
v
1
,
v
2
,
r
1
,
r
2
)
,
v
˙
2
=
η
~
1
(
r
1
,
r
2
)
v
1
+
η
~
2
(
r
1
,
r
2
)
v
2
+
[
-
α
1
*
x
*
2
y
*
(
x
*
+
δ
*
)
(
x
*
+
y
*
)
2
+
η
~
3
(
r
1
,
r
2
)
]
v
1
2
+
η
~
4
(
r
1
,
r
2
)
v
1
v
2
+
[
1
x
*
+
η
~
5
(
r
1
,
r
2
)
]
v
2
2
+
T
~
2
(
v
1
,
v
2
,
r
1
,
r
2
)
,
where ξ
~
i
(
r
1
)
,
i
=
1,2
,
3 and η
~
i
(
r
1
,
r
2
)
,
i
=
1,2
,
3,4
,
5 are polynomial functions with r
1
p satisfying p
⩾
1 and r
1
p
r
2
q satisfying p
+
q
⩾
1, respectively. T
~
i, i
=
1,2, are C
∞ functions in their variables with v
1
j
v
2
k satisfying j
+
k
⩾
3.
Furthermore, ifr
1
=
r
2
=
0, we can see that(BT0)
(52)
J
(
*
,
0,0
)
=
(
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
-
α
1
*
x
*
2
(
x
*
+
y
*
)
2
α
1
*
y
*
2
(
x
*
+
y
*
)
2
-
α
1
*
x
*
y
*
(
x
*
+
y
*
)
2
)
≠
θ
2
×
2
,
(BT1)
[
-
(
x
*
/
(
x
*
+
δ
*
)
)
+
ξ
~
1
(
r
1
)
+
η
~
4
(
r
1
,
r
2
)
]
|
r
1
=
0
,
r
2
=
0
=
-
x
*
/
(
x
*
+
δ
*
)
<
0,
(BT2)
[
-
(
α
1
*
x
*
2
y
*
/
(
(
x
*
+
δ
*
)
(
x
*
+
y
*
)
2
)
)
+
η
~
3
(
r
1
,
r
2
)
]
|
r
1
=
0
,
r
2
=
0
=
-
α
1
*
x
*
2
y
*
/
(
(
x
*
+
δ
*
)
(
x
*
+
y
*
)
2
)
<
0,
(BT3)
the map(53)
(
(
z
1
z
2
)
,
(
r
1
r
2
)
)
⟶
(
(
Q
1
Q
2
)
,
tr
⁡
(
∂
Q
i
∂
z
i
)
,
det
⁡
⁡
(
∂
Q
i
∂
z
i
)
)
,
i
=
1,2
is regular at(54)
(
z
1
z
2
)
=
(
0
0
)
,
(
r
1
r
2
)
=
(
0
0
)
.According to [27, 28], the system (7) undergoes a Bogdanov-Takens bifurcation around the equilibrium point E
* when δ
=
δ
* and α
1
*
=
α
*. This completes the proof.Remark 12.
Furthermore, the quantityω [27, 28] which determines the structure of the bifurcation is given by
(55)
ω
=
sign
⁡
[
α
1
*
x
*
3
y
*
(
x
*
+
δ
*
)
2
(
x
*
+
y
*
)
2
]
=
1
.The bifurcation diagram of system (7) in (
δ
,
α
1
) is presented in Figure 6(a). There are no positive equilibria of system (7) in Domain 1. In the domain bounded by Limit point cycle (green) and Hopf curve (blue), there exist a saddle E
1
* and an unstable E
2
*. In the domain bounded by Hopf curve (blue) and Homoclinic curve (dot and magenta), there is a saddle E
1
* and a stable E
2
* surrounded by an unstable limit cycle. In Domain 2, system (7) has a saddle E
1
* and a stable E
2
*. At BT point, E
* is the unique positive equilibrium of system (7) in the first quadrant, and the local dynamics near E
* is characterized by two trajectories (green curves, Figure 6(b)), which meet at E
*, but with opposite time orientation, and also behave like a cusp-like configuration.


## 3. Discussion

In this paper, we have discussed a ratio-dependent predator-prey model with a (component) mate-finding Allee effect on prey. Allee effects and ratio-dependent functional response are very popular among theoretical ecologists as they can display realistic and complicated dynamical behaviors, which provide better explanations for the ecological observations. And the hyperbolic function that we concerned for mate-finding Allee effects is considered as the most appropriate analytical form from ecological point of view.


### 3.1. The Less Complicated EquilibriumE
0
(
0,0
)

We know that ratio-dependent predator-prey systems with the logistic growth of prey exhibit very complicated dynamic behaviors around the origin [5, 6, 8, 9]; for instance, the origin behaves like a stable node and an unstable saddle at the same time [5]. However, such dynamic properties can not be observed in our system, because E
0
(
0,0
) is always asymptotically stable, which implies that any trajectory, starting from a certain neighborhood of E
0
(
0,0
), goes towards the origin for all system parameters. In other words, the complicated equilibrium E
0
(
0,0
) becomes less complicated because of the mate-finding Allee effect on prey, which provides insights into pest control by introducing a mate-finding Allee effect to reduce the complexity around the origin.We have constructed diagrams with and without a mate-finding Allee effect in the same ratio-dependent predator-prey system for better understanding and comparison. It is easy to see that the origin can be a unstable node and a saddle at the same time in a ratio-dependent predator-prey system where the mate-finding Allee effect is absent in Figure7(b). It is interesting to note that, in both cases, we set α
<
α
1. Generally, the predator growing ability α
1 is assumed to be less than the consumption ability α in predator-prey systems. Hence α
<
α
1 describes the case that predators can reply on other resources, but the prey is still the limiting factor [5, 8]. So we can see obviously how the mate-finding Allee effect in the prey population influences the dynamics of the ratio-dependent predator-prey systems: it increases the extinction risk of both prey and predators, even in the case that prey population is not the only resources of the predator population.The same ratio-dependent predator-prey system with and without a mate-finding Allee effect forg
=
1.3, α
=
0.2, α
1
=
0.8, and β
=
0.3.
(a)
δ
=
0.0066
(b)
δ
=
0


### 3.2. The Impacts ofδ

In order to understand the influence of the relative strength of the mate-finding Allee effectδ on the ratio-dependent prey-predator system, we consider δ as the bifurcation parameter and find that system (7) can exhibit a saddle-node, a subcritical Hopf, and a co-2 Bogdanov-Takens bifurcations with the other bifurcation parameter α
1.It is easy to see that ifδ
>
δ
1, system (7) has no interior equilibria or other boundary equilibria except E
0
(
0,0
); then, E
0
(
0,0
) is globally asymptotically stable. In other words, if the mate-finding Allee effect on the prey population is very strong, any trajectory converges to E
0 as t tends to infinity; that is, both the prey and predators become extinct regardless of their initial density, because the prey population has great difficulty in finding mates at low density and suffers a heavy loss from predation at the same time; thus, such situation leads to a continued decline in predator population which is illustrated in Figure 1(a). Furthermore, different from general predator-prey systems with logistic growth on prey population [5, 6] or those subject to the demographic Allee effects on prey [29], the existence of predator-free boundary equilibria becomes conditional in prey-predator systems where the prey population suffers from a mate-finding Allee effect, and if such equilibria exist, they appear simultaneously with an unstable smaller equilibrium E
10.If0
<
δ
<
δ
2, α
1
>
β, there are two predator-free equilibria (an unstable node E
10 and a saddle E
20 whose stable manifolds lie on x-axis) and two positive equilibria: E
1
* (always a saddle) and E
2
* of system (7) simultaneously. In some cases, E
2
* changes from stable to unstable as the mate-finding Allee effect strength gets stronger. It is interesting that, following from Theorem 5, there do not exist any global asymptotically stable positive equilibria for any g, δ, α, α
1, and β, which means that a mate-finding Allee effect on prey increases the risk of extinction for system (7). In particular, there can exist a stable E
2
* surrounded by an unstable limit cycle for some δ; then, the other trajectories outside the unstable limit cycle are attracted to E
0
(
0,0
) except the stable manifolds of E
1
*, which indicates that the mate-finding Allee effect on prey population may destabilize the system.Though system (7) may experience unstable oscillation for some certain δ, it still needs extra assumption (
C
2
) for other system parameters apart from those basic assumptions (
A
0
) and (
C
1
), which indicates actually that the impacts δ are limit and system (7) itself is highly sensitive to its initial condition and present states. Harder to encounter a receptive mate at low density leads to fewer reproductive output directly, however, we take measures to mitigate the negative effects of mate-finding Allee effects by changing other system parameters.We also leave the co-3 bifurcation problem for future discussion. Another important and interesting issue that requires further exploration is a more sophisticated proof of the existence and stability of the limit cycles. It has been shown that, with different mathematical forms of the Allee effect in a predator-prey system, the number of limit cycles changes [20]. Recently, uniqueness of limit cycle in a Gause-type predator-prey system with an Allee effect has been proved by Olivares et al. [21]. We believe their work on the form of demographic Allee effects will be helpful for future studies on those with component Allee effects.


## 3.1. The Less Complicated EquilibriumE
0
(
0,0
)

We know that ratio-dependent predator-prey systems with the logistic growth of prey exhibit very complicated dynamic behaviors around the origin [5, 6, 8, 9]; for instance, the origin behaves like a stable node and an unstable saddle at the same time [5]. However, such dynamic properties can not be observed in our system, because E
0
(
0,0
) is always asymptotically stable, which implies that any trajectory, starting from a certain neighborhood of E
0
(
0,0
), goes towards the origin for all system parameters. In other words, the complicated equilibrium E
0
(
0,0
) becomes less complicated because of the mate-finding Allee effect on prey, which provides insights into pest control by introducing a mate-finding Allee effect to reduce the complexity around the origin.We have constructed diagrams with and without a mate-finding Allee effect in the same ratio-dependent predator-prey system for better understanding and comparison. It is easy to see that the origin can be a unstable node and a saddle at the same time in a ratio-dependent predator-prey system where the mate-finding Allee effect is absent in Figure7(b). It is interesting to note that, in both cases, we set α
<
α
1. Generally, the predator growing ability α
1 is assumed to be less than the consumption ability α in predator-prey systems. Hence α
<
α
1 describes the case that predators can reply on other resources, but the prey is still the limiting factor [5, 8]. So we can see obviously how the mate-finding Allee effect in the prey population influences the dynamics of the ratio-dependent predator-prey systems: it increases the extinction risk of both prey and predators, even in the case that prey population is not the only resources of the predator population.The same ratio-dependent predator-prey system with and without a mate-finding Allee effect forg
=
1.3, α
=
0.2, α
1
=
0.8, and β
=
0.3.
(a)
δ
=
0.0066
(b)
δ
=
0


## 3.2. The Impacts ofδ

In order to understand the influence of the relative strength of the mate-finding Allee effectδ on the ratio-dependent prey-predator system, we consider δ as the bifurcation parameter and find that system (7) can exhibit a saddle-node, a subcritical Hopf, and a co-2 Bogdanov-Takens bifurcations with the other bifurcation parameter α
1.It is easy to see that ifδ
>
δ
1, system (7) has no interior equilibria or other boundary equilibria except E
0
(
0,0
); then, E
0
(
0,0
) is globally asymptotically stable. In other words, if the mate-finding Allee effect on the prey population is very strong, any trajectory converges to E
0 as t tends to infinity; that is, both the prey and predators become extinct regardless of their initial density, because the prey population has great difficulty in finding mates at low density and suffers a heavy loss from predation at the same time; thus, such situation leads to a continued decline in predator population which is illustrated in Figure 1(a). Furthermore, different from general predator-prey systems with logistic growth on prey population [5, 6] or those subject to the demographic Allee effects on prey [29], the existence of predator-free boundary equilibria becomes conditional in prey-predator systems where the prey population suffers from a mate-finding Allee effect, and if such equilibria exist, they appear simultaneously with an unstable smaller equilibrium E
10.If0
<
δ
<
δ
2, α
1
>
β, there are two predator-free equilibria (an unstable node E
10 and a saddle E
20 whose stable manifolds lie on x-axis) and two positive equilibria: E
1
* (always a saddle) and E
2
* of system (7) simultaneously. In some cases, E
2
* changes from stable to unstable as the mate-finding Allee effect strength gets stronger. It is interesting that, following from Theorem 5, there do not exist any global asymptotically stable positive equilibria for any g, δ, α, α
1, and β, which means that a mate-finding Allee effect on prey increases the risk of extinction for system (7). In particular, there can exist a stable E
2
* surrounded by an unstable limit cycle for some δ; then, the other trajectories outside the unstable limit cycle are attracted to E
0
(
0,0
) except the stable manifolds of E
1
*, which indicates that the mate-finding Allee effect on prey population may destabilize the system.Though system (7) may experience unstable oscillation for some certain δ, it still needs extra assumption (
C
2
) for other system parameters apart from those basic assumptions (
A
0
) and (
C
1
), which indicates actually that the impacts δ are limit and system (7) itself is highly sensitive to its initial condition and present states. Harder to encounter a receptive mate at low density leads to fewer reproductive output directly, however, we take measures to mitigate the negative effects of mate-finding Allee effects by changing other system parameters.We also leave the co-3 bifurcation problem for future discussion. Another important and interesting issue that requires further exploration is a more sophisticated proof of the existence and stability of the limit cycles. It has been shown that, with different mathematical forms of the Allee effect in a predator-prey system, the number of limit cycles changes [20]. Recently, uniqueness of limit cycle in a Gause-type predator-prey system with an Allee effect has been proved by Olivares et al. [21]. We believe their work on the form of demographic Allee effects will be helpful for future studies on those with component Allee effects.



---

*Source: 673424-2014-10-20.xml*