DYNAMICS IN A DIFFUSIVE PREDATOR-PREY SYSTEM WITH STAGE STRUCTURE AND STRONG ALLEE EFFECT

. In this paper, we consider the dynamics of a diffusive predator- prey system with stage structure and strong Allee effect. The upper-lower solution method and the comparison principle are used in proving the nonneg- ativity of the solutions. Then the stability and the attractivity basin of the boundary equilibria are obtained, by which we investigated the bistable phe- nomena. The existence and local stability of the positive constant steady-state are investigated, and the existence of Hopf bifurcation is studied by analyzing the distribution of eigenvalues. On the center manifold, we studied the crit-icality of the Hopf bifurcation by the normal form theory. Some numerical simulations are carried out for illustrating the theoretical results.

1. Introduction. Stage structure in a predator-prey model may improve the assumption that each individual predator has the same ability to capture prey. In fact, most immature predators do not have the ability to attack prey. Thus, immature individuals and mature individuals of the predator should be divided by a fixed age since the ability of attacking prey are obviously different. There has been a fair amount of work on modeling stage structured model consisting of immature and mature individuals for species (see [5,6,22,18] and the references therein). Wang and Chen [19] found that an orbitally asymptotically stable periodic orbit exists in that model. They also established the condition for the permanence of the populations.
Saito and Takeuchi [14] proposed a stage structured predator-prey model as follows: Here, x(t), y(t), Y (t) denote the density of the prey, the immature and mature predator at time t, respectively. The immature individuals and mature individuals are separated by age τ . The term a 31 e −r2τ x(t − τ )Y (t − τ ) in the second equation represents the immature individuals (with the mature birth rate a 31 ) born at time t − τ and alive time t (with the immature death rate −r 2 ). They obtained the conditions for the local asymptotic stability and global attractivity of the interior equilibrium of the model. Apart from this, the conditions under which the Hopf bifurcation may occur in system (1) are derived in [12].
Since the prey and predator always distribute inhomogeneously in different locations, the diffusion has been taken into considered in many ecological models (see [3,27,8,24]). Xu and Wei [23] considered a diffusive budworm model with a structured population model subjected to the Neumann boundary conditions as follows: where u(x, t) is the mature budworm density at location x and time t, τ is the maturation delay. They investigated the stability of the steady state and the existence of Hopf bifurcation. The global existence of the periodic solutions is also established.
In recent years, the Allee effect growth of a population has been studied extensively [17,16]. A strong Allee effect refers to the phenomenon that one population has a negative growth when the size of the population is below the threshold value [15]. This is a universal phenomenon in nature since many preys are social animals, they can only live in groups, like ants, birds, rabbits and so on [1].
Therefore, in order to describe the predator-prey model more accurately, we combine the stage structure in predator and the strong Allee effect growth in prey, and impose a no-flux boundary condition, so it is a closed ecosystem. Such a diffusive system with strong Allee effect and stage structure takes the following form: where u(t, x) and v(t, x) are the population densities of the prey and the mature predator, respectively. d i (i = 1, 2) is the diffusion coefficient. k represents the predation rate. r is the average mortality rate of the mature predator. θ is the conversion rate. √ b denotes the half saturation coefficient. τ is the maturation time delay of the immature predator. r 0 is the average mortality rate of the immature predator. The Laplacian operator ∆ shows the diffusion effect. The habitat Ω ∈ R n is a bounded domain with smooth boundary for n ≥ 1, and ν is the outer normal direction for x ∈ ∂Ω. d 1 , d 2 , α, k, r, θ, b and r 0 are all positive constants, M represents the Allee threshold of prey with 0 < M < N , and N denotes the carrying capacity of the prey. The term represents the number of immature predators that were born at time t − τ which still survive at time t and are transferred from the immature stage to the mature stage at time t.
With a nondimensionalized change of variables: Dropping the tildes of u and v, we obtain the simplified dimensionless system of equations: In the remaining part of this article, we focus on system (4). For the new parameters, K is a rescaled carrying capacity of the prey. r is the mortality rate of the predator, and m is the strength of the interaction. To maintain the meaningful of the system, we assume that m > r.
In this paper, we have proved the global existence and non-negativity of the solutions to (4). Several groups of conditions are given under which global asymptotical behaviour of the solutions can be determined. The attractivity basin of the extinction equilibrium and the local stability of boundary equilibria is investigated. Particularly, we show that large amount of predator will always drives both predator and prey into extinction, which is called over-predation. This is a character of predator-prey system with strong Allee effect. Over-predation exists widely in nature. For example, people have driven Emberiza aureola, Manis pentadactyla and Chinese sturgeon to the edge of extinction by capturing and eating them. We also proved that the bistable phenonmenon may occur which means the system dynamics has strong selectivity with respect to the initial conditions. Besides, we investigated Hopf bifurcations at the positive steady state: when the delay τ decreases passing through a critical value, the positive steady state loses its stability and the stable or unstable time-periodic oscillations appear. This explains why the number of the predator and prey vary periodically, as reported in [19,12]. Finally, by using a group of parameters, we illustrate the rich dynamical behaviours of our model.
Methods for analyzing diffusive predator-prey systems have been developed since past decades. In this paper we apply some classical techniques like upper-lower solution method, comparison methods, geometry criterion and bifurcation theory to proceed the analysis. We emphasize several points in our derivations: for one thing, the proof of over-predation requires subtle analysis, since stage structure makes (4) a diffusive functional differential equation; for another, stability analysis of the equilibria is not trivial since the characteristic equation is a series of equations with coefficients dependant with delay. Besides, when seeking for the conditions under which the positive equilibrium E 3 is stable, we searched the stable region and then the critical bifurcation point with the decrease of τ and draw the conclusion.
The rest of the paper are structured in the following way. In section 2, we prove the nonnegativity of the solutions in system (4). In section 3, by analyzing the distribution of the roots of the associated characteristic equation, we study the stability of the boundary equilibria and the positive equilibrium. Taking τ as a parameter, using the approach of Beretta and Kuang [2], we show that the positive steady-state can be destabilized through a Hopf bifurcation. In section 4, we illustrate our results by numerical simulations and find that the system may exhibit stable spatially homogeneous periodic solutions or transient periodic solutions near different Hopf bifurcation points. In section 5, a conclusion is given to complete this paper.
2. Nonnegativity of solutions. In this section, the nonnegativity of the solutions is proved. Let Ω ∈ R N be a bounded domain with smooth boundary. Define with the standard inner product ·, · , and Now we give the conclusion as follows.
On the other hand, v(x, t) in system (4) can be rewritten as: Let Since v 0 (x, t) ≥ 0 and v 0 (x, 0) ≡ 0, using the comparison principle discussed above, we can get that the solution W (x, t) of Eq.
From the proof above, we may get the following remark.

Stability and bifurcation analysis of the constant steady states.
In the following of this paper, we consider system (4) on Ω = (0, lπ). Then Clearly, the constant equilibria of (4) are determined by the following equations Solving Eq.(11), system (4) has three non-negative boundary constant steady states E 0 , E 1 and E 2 and one possible positive constant steady state E 3 . The steady states are given by E 0 = (0, 0), Notice that K = N M > 1 and m > r, we have the following result. Lemma 3.1. The positive equilibrium E 3

Characteristic equation.
In this section, we aim at the stability analysis of the equilibria E 0 , E 1 , E 2 and E 3 . Denote , X), system (4) can be written as an abstract differential equation as follows where and F : C → X is given by Linearizing system (12) at one of the constant steady states, for simplicity, denote the constant steady state as (u 0 , v 0 ), thus we obtain where L : C → X is given by

YUYING LIU, YUXIAO GUO AND JUNJIE WEI
and From Wu [21], we obtain that the characteristic equation for the linear system (13) satisfying It is well known that the operator u −→ ∆u with ∂ ν u = 0 at x = 0 and lπ has eigenvalues − n 2 l 2 (n ∈ N 0 ) with the corresponding eigenfunctions cos be an eigenfunction of the operator D∆ + L with eigenvalue λ, see [25,26], and substituting this eigenfunction into (14), we get Therefore, the characteristic equation (14) at the equilibrium (u 0 , v 0 ) is equivalent to where The stability of the steady state (u 0 , v 0 ) can be determined by the distribution of the roots of Eq. (15): if all the roots of Eq.(15) have negative real parts, then (u 0 , v 0 ) is locally asymptotically stable; if at least one root of Eq.(15) has positive real parts, then (u 0 , v 0 ) is unstable; if any root has zero real part and other roots all have negative real parts, then the stability of (u 0 , v 0 ) cannot be determined by the linearized system directly.
For simplicity, denote Then, the characteristic Eq.(15) at (u 0 , v 0 ) turns out to be

3.2.
Local stability of E 0 , E 1 and E 2 . First we investigate the local stability of E 0 (0, 0) and E 1 (1, 0). By direct calculation on the roots of Eq.(15) at E 0 , E 1 , we can get the following theorem.
As to the stability of E 2 , first we introduce a lemma.
Applying Lemma 3.3, we have the conclusion as follows. where Since λ 1(n) = − a < 0, to get the stability of E 2 (K, 0), we only need to check the sign of real part of the roots in λ + d − f e −λτ = 0. When The roots λ 2(n) , as a function of n 2 , reaches the maximum when n = 0 and λ 2(0) = mK 2 b+K 2 − r.
Clearly a, b are real numbers. From Lemma 3.3, we get that all the roots of (z + a)e z + b = 0 have negative real parts if and only if Eq.(16) holds, that means, where ξ = − dτ tan ξ, 0 < ξ < π.
3.3. The attractivity of E 0 . In this section we investigate the attractivity of To prove v(x, t) → 0 when supΩ u 0 (x, 0) < 1, we first introduce a lemma.
with homogeneous Neumann boundary conditions and non-negative initial conditions converges to zero (uniformly in x) as t → ∞ if and only if f (y) < δy for all y > 0.
From Lemma 3.5, every solution Q(x, t) of Eq. (20) converges to Q = 0 (uniformly in x) as t → ∞, so does the solution v(x, t) of equation (19). Combining the local stability of E 0 , we have the following theorem.  (u(x, t), v(x, t)) tends to (0, 0) as t → ∞ for all x ∈ Ω. That is, E 0 (0, 0) has a basin of attraction The following theorem states that for a given initial value t) is big enough, then the corresponding solution (u(x, t), v(x, t)) of system (4) originating from (u 0 (x, t), v 0 (x, t)) will also be attracted by (0, 0). In fact, this describes the over-predation phenomenon. Over-predation devotes the extinction of the prey.
In the following, we will prove when sup Since sup x∈Ω u(x, T 0 ) ≥ 1, it is easy to get u(x, T 0 ) ≡ 0, then from Theorem 2.1, we have u(x, t) > 0 for t > T 0 and x ∈Ω. For convenience, take Then v(x, t) ≥ v 1 (x, t) from the comparison principle of parabolic equation for any where c n satisfies g(x, t).
In the following, we will prove if v * 0 and T 2 are large enough, then u(x, T 1 + T 2 ) < 1.

Direct calculation implies that if we choose
and for any x ∈Ω, This contradicts with the assumption that u(x, t) ≥ 1 for all t > T 1 , which means for a given u 0 (x, t), if v * 0 is large enough and v 0 (x, t) ≥ v * 0 , the solution u(x, t) will be smaller than 1 eventually. Therefore, (u(x, t), v(x, t)) tends to (0, 0) as t → ∞ from Theorem 3.6. Since ϵ is chosen arbitrarily, it is clear that v * 0 depends only on the fixed parameters and T 0 which depends on u 0 (x, t).
The above theorem implies that (0, 0) is always a locally stable steady state with basin of attraction including all large v 0 (x, t) for a given u 0 (x, t). Thus system (4) is bistable if there is another locally stable steady state or periodic solutions. We have the following conclusion. (4) is bistable when re dτ > mK 2 b+K 2 in the following sense: if re dτ > mK 2 b+K 2 , the die-out equilibrium E 0 (0, 0) of system (4) is locally asymptotically stable; meanwhile, the boundary equilibrium E 2 (K, 0) is locally asymptotically stable.

Local stability and Hopf bifurcation of
. For convenience, we rewrite the hypotheses in the following: From Lemma 3.1 we know that (H1) guarantees the existence of the positive equilibrium E 3 .

Remark 2.
(H1) holds provided that either of the following conditions is satisfied: Since case (ii) can be dealed in the same way as case (i), in the following of this paper, we focus on the case (i).
For convenience, we further make the following hypotheses: First, we analyze the local stability of E 3 when τ is close enough to τ max . For E 3 (u * (τ ), v * (τ )), the characteristic equation becomes where Theorem 3.9. E 3 (u * , v * ) is locally asymptotically stable, provided that Proof. For simplicity, denote s 1 = . Then the characteristic equation (22) becomes (23) with nonnegative real part, then there must exists an integer k such that λ p = a + ib is the root of the (k + 1)th equation in (23). Substituting λ p = a + ib into the (k + 1)th equation of (23), we have the right side of the Eq. (24) is re −(a+ib)τ , thus we have The left side of Eq.(24) becomes since s 1 , s 2 > 0 and a ≥ 0, we have This contradicts with the assumption. Thus, all roots of Eq.(22) have negative real part, therefore, E 3 (u * , v * ) is locally asymptotically stable.

Theorem 3.10. Letτ be such that
This can be proved easily by the fact that u * (τ ) is an increasing function of τ .
Here we omit the proof.

Remark 3.
Ifτ ≤ 0, then E 3 (u * , v * ) is stable for all τ ∈ (0, τ max ). Forτ > 0 and taking τ as the bifurcation parameter, we may investigate the existence of Hopf bifurcation between 0 andτ . Without loss of generality, here we only consider the case whenτ > 0. Since is an increasing function of τ , therefore, when τ reduces fromτ , there must exists a critical value τ c such that when τ ∈ (τ c ,τ ), and when τ = τ c , Let τ = max{τ c , 0}, now we can investigate the Hopf bifurcation when τ ∈ (τ ,τ ]. Obviously When τ ∈ (τ ,τ ], we can apply the geometry criterion due to Beretta and Kuang [2] to study the stability of E 3 . Here we need to verify the following properties for n ∈ N 0 and τ ∈ (τ ,τ ]: Proof. Given the fact that which means (iii) follows. Let F n (ω, τ ) be defined as in (iv) with the following expression It is obvious that property (iv) is satisfied, and by the Implicit Function Theorem, (v) is also satisfied.
Since a 1 (τ ) > 0, we get that F n (ω, τ ) = 0 has positive root only and if only a 2 (τ ) < 0. In fact we have a 2 (τ ) = ( a d − a f − b c)( a d + a f + b c) and a d − a f − b c > 0. Thus, F n (ω, τ ) = 0 has positive root only and if only a d + a f + b c < 0. Denote the only positive root as ω n (τ ), then ω n (τ ) satisfies This means that ω n (τ ) makes sense if and only if a 2 1 (τ ) − 4a 2 (τ ) > a 1 (τ ), which requires a d + a f + b c < 0. For simplicity, denote H n (τ ) = a d + a f + b c.
One can check that iω n (τ * )(ω n (τ * ) > 0) is a purely imaginary root of Eq. (22) if and only if τ * is a root of the function S m n , defined by The following is the result introduced by Beretta and Kuang [2].
Lemma 3.11. For a fixed n 0 ∈ N 0 , assume that the function S m n0 (τ ) has a simple positive root τ * ∈ I n0 for some m ∈ N 0 , then a pair of simple purely imaginary roots ±iω n0 (τ * ) of equation (22) exists at τ = τ * and Since Therefore, this pair of simple conjugate purely imaginary roots crosses the imaginary axis from left to right if δ(τ * ) = 1 and from right to left if δ(τ * ) = −1.

Remark 4.
It can be easily observed that I n+1 ⊂ I n , S m n (τ ) > S m+1 n (τ ) for all τ ∈ I n , S m n (τ ) > S m n+1 (τ ) for all τ ∈ I n+1 and S m n (τ ) approaches to negative infinity atτ . Thus, if S m0 n0 has no zero in I n0 , then for any n ≥ n 0 , m ≥ m 0 , τ ∈ I n , we have S m n (τ ) ≤ S m0 n0 (τ ) and S m n have no zeros in I n .
Applying Corollary in [13], we can draw the conclusion: If m b+1 < r < mK 2 b+K 2 and τ ∈ (τ , τ max ), then all roots of Eq.(22) have negative real parts when τ ∈ (τ k , τ max ) and at least a pair of roots have positive real parts when τ ∈ (τ k−1 , τ k ). Furthermore, all other roots of Eq. (22), except a pair of purely imaginary roots, have negative real parts when τ = τ k . Now we can state the following theorem on the existence of a Hopf bifurcation at the positive steady state. (ii) If J = ∅, then the steady state E 3 of Eq.(4) is locally asymptotically stable for τ ∈ (τ k , τ max ) and unstable for τ ∈ (τ k−1 , τ k ) with a Hopf bifurcation occurring at E 3 when τ = τ i ∈ J. Theorem 3.12 gives some sufficient conditions to ensure that Eq.(4) undergoes a Hopf bifurcation at E 3 . Next, under the conditions of Theorem 3.12(ii), we shall use the center manifold and normal form theories presented by Wu [21] and Faria [4] to study the direction of Hopf bifurcation and the stability of the bifurcating periodic solutions from E 3 . As the details are given in the Appendix, we summarize the results in the following theorem.
Under this parameter set, one can easily see that m b+1 ≈ 0.0513, mK 2 b+K 2 ≈ 1.5220 and thus m b+1 < r < mK 2 b+K 2 . By calculation, we obtain that τ max ≈ 42.0034,τ ≈ 17.6558.    In order to gain the set I n , we picture the graphs of (s 1 + s 2 )(u) and H n (τ ) in Figure 1. From the graph of (s 1 + s 2 )(u) in Figure 1, we get that u(τ c ) ≈ 5.9108, by calculation, we can get that τ ≈ 4.2894. Note that τ ∈ I n is equivalent to H n (τ ) < 0. Hence we have I 0 = (4.29, 11.51), I 1 = (4.29, 10.13) I 2 = (4.29, 5.37), I n = ∅, n ≥ 3. Accordingly, the pictures of S m n on I n can be drawn clearly (see Figure 2). From . Therefore, the positive steady state E 3 is asymptotically stable when τ ∈ (τ 1 , τ max ), and unstable when τ ∈ (τ 0 , τ 1 ), as well as Hopf bifurcation takes place when τ = τ i ∈ J. The stability of E 3 is illustrated by Figure 3, the bifurcating periodic solutions are illustrated by Figure 4 and Figure 5.    respectively. It follows that the direction of Hopf bifurcation is forward at τ 0 and backward at τ 1 . The bifurcating periodic solutions exist for τ < τ 1 and close to τ 1 , these bifurcating periodic solutions are orbitally asymptotically stable. The bifurcating periodic solutions exist for τ > τ 0 and close to τ 0 and these bifurcating periodic solutions are unstable. These are illustrated in Figure 4 and Figure 5. When τ > τ max , the positive steady state E 3 disappears, and by Theorem 3.4, the steady state E 2 is local stable (see Figure 6). By theorem 3.7, when the initial predator value is large enough, the corresponding solution of system (4) will tend to (0, 0) uniformly for x ∈ Ω as t → ∞ (see Figure 7). Figure 7 illustrates the excessive predation phenomenon.

Conclusion.
In this paper, a diffusive predator-prey system with strong Allee effect and stage structure is investigated. Theoretical analysis indicates that both strong Allee effect and stage structure would make the dynamics of the system more complicated. These contains over-predation phenomenon, bistable phenomenon caused by the strong Allee effect and time periodic orbits caused by the stage structure. The attractivity basin of the die-out equilibrium and the local stability of boundary equilibria are investigated, suggesting that under certain condition, the initial value of the prey is important as well as the Allee threshold. If the initial value of the prey is smaller than the Allee threshold, the prey species will die out, and so does the predator species, see Theorem 3.6; if the initial value of the prey is larger than the Allee threshold, and when the initial value of the predator is large enough, both the prey and the predator will die out, see Theorem 3.7. Moreover, the existence and stability of the positive steady state are studied, the Hopf bifurcation at the positive steady state is investigated. We show that when the delay τ decreases to a critical value, the positive steady state will lose its stability and Hopf bifurcation will occur. Since the coefficients of the corresponding characteristic equation depend on the delay τ , we prove that all roots of the characteristic equation have negative real parts when τ is close enough to τ max , see Theorem 3.9 and 3.10. By the theory of normal form and center manifold method, the conditions for determining the direction of the bifurcation and the stability of bifurcating periodic solution are derived in the appendix. We conclude that the time delay τ can affect the coexistent of the prey and the predator. Under certain condition, when τ k < τ < τ max , these two species will coexistent; when τ ∈ (τ k−1 , τ k ), the prey and the predator will exhibit oscillatory behaviour.
Appendix. In this appendix, by applying the center manifold theorem and normal form theory of the partial differential equations with delay [4,21], we study the direction of Hopf bifurcation and stability of bifurcating periodic solutions from the positive equilibrium E 3 near τ = τ * under the conditions of Theorem 3.12(ii). Denote while With the center manifold theorem, there exists a center manifold C 0 and we can write W in the following form on C 0 nearby (0, 0): For solution U t ∈ C 0 , denote Therefore the system restricted to the center manifold is given bẏ by calculation, we can get that g 20 = τ * M From (36) and (44), for θ ∈ [−1, 0), we have in which E 1 and E 2 are both 2-dimension vectors in X, and can be determined by setting θ = 0 in H. Set θ = 0 and by (44) and (45), we have The termsF 20 andF 11 are elements in the space C , and they can be denoted bỹ Denote substituting E 1 , E 2 into the equation (46), then for n = 0, 1, · · ·, we have (A 0 − 2iω n0 τ * I)E n 1 e 2iωn 0 τ * θ | θ=0 = − F 20 , β n b n , A 0 E n 2 | θ=0 = − F 11 , β n b n . Thus, E n 1 and E n 2 could be calculated by F 20 = 2a 1 + 2a 2 q 1 2a 5 q 1 e −2iωn 0 τ * + 2a 6 e −2iωn 0 τ * , F 11 = 2a 1 + a 2 (q 1 +q 1 ) a 5 (q 1 + q 1 ) + 2a 6 .