Bogdanov–Takens bifurcation in a predator–prey model

In this paper, we investigate a class of predator–prey model with age structure and discuss whether the model can undergo Bogdanov–Takens bifurcation. The analysis is based on the normal form theory and the center manifold theory for semilinear equations with non-dense domain combined with integrated semigroup theory. Qualitative analysis indicates that there exist some parameter values such that this predator–prey model has an unique positive equilibrium which is Bogdanov–Takens singularity. Moreover, it is shown that under suitable small perturbation, the system undergoes the Bogdanov–Takens bifurcation in a small neighborhood of this positive equilibrium.


Introduction
In this article, we will analyze the following predator-prey system with an age structure ∂u(t, a) ∂t + ∂u(t, a) ∂a = −μu(t, a), for a ≥ 0, , u(0, ·) = u 0 ∈ L 1 + ((0, +∞) , R) and V (0) = V 0 ≥ 0, (1.1) where the number μ denotes the death rate of the predator, t is time, and a is the chronological age (i.e., the age of individuals since they were born). η denotes the maximal growth rate of predator. The function u(t, a) is the density of population of predators with respect to the age at time t. This means that the number of predator with age between a 1 and a 2 is a2 a1 u(t, a)da.
The number V (t) is the total number of prey at time t. The parameter r = b − d is the intrinsic growth rate of the prey (where b and d denote the birth rate and death rate, respectively) and K > 0 is the carrying capacity of the prey. The function β ∈ L ∞ + ((0, +∞) , R) represents a fraction of predator which can produce new individuals by eating the prey. In particular, when β(a) ≡ 1 system (1.1) becomes the ordinary differential equation Note that the bifurcations of system (1.2) have been analyzed in Xiao and Ruan [52] and it has been shown that system (1.2) can undergo Bogdanov-Takens bifurcation. Namely, they prove the existence of a series of bifurcations including saddle-node bifurcation, Hopf bifurcation, and homoclinic bifurcation. The functional response V (t)U (t) h+V (t) 2 is a so-called non-monotonic functional response. We refer to Ruan and Xiao [38] for a detailed presentation of model (1.2) and an overview on this topic. We also refer to Cushing and Saleem [13] for a presentation of prey-predator systems with an age structure.
Here it makes sense to assume that the function β(a) is not identically equal to 1 because the young predator cannot reproduce even if they eat some preys. Therefore, in order to take care of this fact we will assume that β is a step function of the form Moreover, it will be convenient for the computation to assume that +∞ 0 β(a)e −μa da = 1.
One may observe that without loss of generality we can make this assumption. Because by replacing β(a) by β(a) +∞ 0 β(s)e −μs ds and replacing η by η +∞ 0 β(s)e −μs ds the system is unchanged.
Age-structured models, described by hyperbolic partial differential equations, have been studied by many researchers (see the monographs of Cushing [12], Iannelli [24], Webb [49], and the references cited therein). Various approaches have been developed to study these models. In order to analyze the qualitative properties of such a system, it is usually convenient to combine (a) integration along the characteristics combined with Volterra integral equations (Webb [49] , Iannelli [24]) and (b) integrated semigroups method (Thieme [44,46,47], Magal [31], Magal and Ruan [32,34]). The classical study on age-structured models focused on the existence, boundedness and stability of solutions. It has been shown that some age-structured models exhibit non-trivial periodic solutions induced by Hopf bifurcation (see Prüss [37], Cushing [12], Swart [40], Kostava and Li [26], and Bertoni [2]). Recently, there has been great interest in bifurcation analysis on degenerate equilibria of age-structured models. Since age-structured models can be rewritten as abstract semilinear equations with non-dense domain (Thieme [44,46,47], Magal and Ruan [32,34]), Magal and Ruan [33] developed the center manifold theory for abstract semilinear Cauchy problems with non-dense domain. Based on the center manifold theorem proved in Magal and Ruan [33], a Hopf bifurcation theorem has been presented for abstract non-densely defined Cauchy problem in Liu et al. [28]. These theorems have been successfully applied to study the existence of Hopf bifurcation for some age-/size-structured models (see [3,[9][10][11]34,35,39,48]). Recently, much attention has been focused on high-codimensional bifurcations, since they may reveal some complex dynamical behaviors.
Bogdanov-Takens bifurcation is one of the most important high-codimensional bifurcations in nonlinear dynamics. It is a well-studied example of a bifurcation with codimension two and named after Bogdanov [4,5] and Takens [41,42], who independently and simultaneously described this bifurcation. For more results about Bogdanov-Takens bifurcation, see, for example, Chow and Hale [6], Dumortier et al. [16], Arnold [1], Chow et al. [7], Guckenheimer and Holmes [22], Kuznetsov [27], Xiao and Ruan ([38,51]), and others. Bogdanov-Takens bifurcation occurs also in infinite dimensional differential equations associated with functional differential equations and a few other partial differential equations. In the context of functional differential equations, we refer to the book by Hale et al. [23] and papers [18,25,50,52,53]. In the context of partial differential equations with delays, we refer to [19][20][21]30]. But to the best of our knowledge, there is no result on Bogdanov-Takens bifurcation for an age-structured model which is an hyperbolic partial differential equations with nonlinear and non-local boundary conditions.
In this article, our aim is to study whether system (1.1) also undergoes the Bogdanov-Takens bifurcation near a positive equilibrium. Hence, we first look for some conditions of parameters which guarantee system (1.1) has a positive equilibrium who is Bogdanov-Takens singularity by applying the normal form theory developed by Liu, Magal, and Ruan [29] for the non-densely defined abstract Cauchy problems, then we choose suitable small perturbation of parameters such that system (1.1) can undergo the Bogdanov-Takens bifurcation. More precisely, we obtain that there exists a unique positive equilibrium which is Bogdanov-Takens singularity of system (1.1) if η = 2 √ h and K = 2 √ h, and prove the following small perturbation system is a versal unfolding of this Bogdanov-Takens singularity, where ( α 1 , α 2 ) is a very small parameter vector, other parameters h, r, and μ are any fixed, and η = 2 √ h and K = 2 √ h. Here we take the parameters α 1 , α 2 as the bifurcation parameters. Of course, the original system (1.1) corresponds to the case where α 1 = 0 and α 2 = 0. Here we will study the bifurcations of system (1.3) near the unique positive equilibrium as parameters α 1 and α 2 change in the small neighborhood of the origin. The aim of the present paper is to show how to apply the normal form theory in the context of age-structured models. In this paper, we will prove that Bogdanov-Takens bifurcation occurs for system (1.3). Some more biologically meaningful examples are left for further investigation.
The paper is organized as follows. In Sect. 2, we formulate an age-structured model based on system (1.1) as a non-densely defined Cauchy problem. Then we consider the existence of the positive equilibrium, linearize the system around the positive equilibrium, and study the spectral properties of the linearized equation in Sect. 3. In Sect. 4, the eigenvalue problem for the linearized system of (1.1) around the unique positive equilibrium is investigated and the normal form and center manifold theory for semilinear equations with non-dense domain is used to carry out the analysis of the Bogdanov-Takens bifurcation. A brief conclusion is provided in the last section.

Preliminary
We first reformulate the model (1.1) as a semilinear Cauchy problem with non-dense domain. By setting we can rewrite the second equation in system (1.1) as follows Thus, by setting w(t, a) = u(t, a) v(t, a) , we obtain the equivalent system of system (1.1) Following the results developed in Thieme [44] and Magal [31], we consider the Banach space We observe that L is non-densely defined since  Then setting , we can rewrite system (1.1) as the following non-densely defined abstract Cauchy problem . From the results of Magal [31] and Magal and Ruan [34], we obtain the global existence, uniqueness, and positive of solutions for system (2.1).

Equilibria and linearized equation
The equilibrium solutions of Eq. (2.1) are obtained by solving the equation and we obtain the following lemma.
Lemma 3.1. The system (2.1) has always the equilibria Furthermore, there exists a unique positive equilibrium of system (2.1) In the following, we assume that η = 2 √ h and K = √ h. Now we make the following change of variable Therefore, the linearized equation of (3.1) around the equilibrium 0 is given by Then (3.1) can be written as is a linear operator and satisfying H(0) = 0 and DH(0) = 0. Next we will study the spectral properties of the linearized equation of (3.1) in order to investigate the dynamical behavior for system (3.1).
By applying the results of Liu et al. [28], we obtain the following result. Lemma 3.2. If λ ∈ Ω, then λ ∈ ρ(L) and From the above we have the following result. otherwise.

Now it remains to precise the spectral properties of
Then we obtain the system From the formula of DB(w) we know that Lemma 3.4. The following results hold.
Proof. Assume that λ ∈ Ω and det(Δ(λ)) = 0. From the above discussion, we have By Lemma 3.2 we obtain (3.5) and (ii) follows. Therefore, we have {λ ∈ Ω : det(Δ(λ)) = 0} ⊂ ρ(L+C)∩Ω, and Conversely, assume that λ ∈ Ω and det(Δ(λ)) = 0. We claim that we can find In fact, set We can find a nonzero solution of (3.7) if and only if we can find From the above arguments, we can obtain Since by assumption det(Δ(λ)) = 0, we can find δ = 0 such that Δ(λ)δ = 0. So we can find (3.7), and thus λ ∈ σ P (L + C). Hence, From the explicit formula for the resolvent, we deduce the following lemma. Since L is a Hille-Yosida operator, and DF(x) is bounded, A is also a Hille-Yosida operator. Consequently, A 0 generates a strongly continuous semigroup {T A0 (t)} on X 0 . In order to apply the center manifold theorem and normal form theory, we need to study the essential growth bound of A. The essential growth bound ω 0,ess (A) ∈ [−∞, +∞) of A is defined by where T A (t) ess is the essential norm of T A (t) defined by is the Kuratovsky measure of non-compactness. By using a perturbation result we obtain the following estimation.
Proposition 3.6. The essential growth rate of the strongly continuous semigroup generated by A 0 is strictly negative, that is, Since DF(x) is a compact bounded linear operator, we can apply the perturbation results in Thieme [43] or Ducrot et al. [14] to deduce that
(λI+Q)dl ϕ(s)ds . We need to differentiate with respect to λ the above formula. We first have Thus, we obtain where ψ is defined in (4.4).
From the above results, we obtain a state space decomposition with respect to the spectral properties of the linear operator A. More precisely, the projector on the linear center manifold is defined by (4.6) we obtain the basis {χ 1 , χ 2 } of X c defined by

Note that
The matrix of A c in the basis {χ 1 , χ 2 } of X c is given by or by using a usual symbolic matrix computation formula We now consider the normal form of system We decompose V 2 (X c , X c ) into the direct sum where R c 2 := R(Θ c 2 ) is the range of Θ c 2 , and C c 2 is some complementary space of R c 2 into V 2 (X c , X c ). Let D(A)). (4.10) By applying Theorem 4.4 in Liu, Magal, and Ruan [29], we obtain that after the following change of variable locally around 0 H(y(t)), for t ≥ 0, (4.12) and the reduced equation of system (4.12) has the following form which is in the normal form at the order 2, and the remainder term R c ∈ C 3 (X c , X c ) is a rest of order 3, that is to say that D j R c (0) = 0 for each j = 1, 2.
In the following, we will compute

d(a−s) ϕ(s)ds
and We observe that for each By projecting on X c , we obtain (4.14) Set We shall compute (4.12) expressed in terms of the basis {χ 1 , χ 2 }. Consider V 2 (R 2 , R 2 ) which denotes the linear space of the homogeneous polynomials of degree 2 in two real variables, x 1 , x 2 with coefficients in R 2 . The operators Θ c 2 considered in (4.8) now act in the spaces V 2 (R 2 , R 2 ) and satisfies Θ c 2 (G 2,c ) It is known [7,22,27,42] that a normal form for Bogdanov-Takens singularity which gives the flow on the center manifold is . In order to obtain from equation (4.13) the second-order terms in the above normal form for Bogdanov-Takens singularity, we need to choose a complementary space for R(Θ c 2 ) in V 2 (R 2 , R 2 ). The canonical basis of V 2 (R 2 , R 2 ) has six elements: Their images under Θ c 2 are, respectively, We choose the following Note that in the above formula we are using a matrix symbolic computation.
Remember we assumed that K = 2 √ h; therefore, we obtain that (4.13) expressed in terms of the basis From the above analysis, we obtain the following theorem.

Bogdanov-Takens bifurcation
To observe whether system (4.1) can undergo Bogdanov-Takens bifurcation under a small perturbation, we rewrite the linear operator of system (4.1) We introduce two small parameters α = ( α 1 , α 2 ) by setting α 1 = α 1 , α 2 = β 1 − 1 to system (4.1). Then the following system is a small perturbation system of system and with α 1. Note that O(0, 0) is an equilibrium of system (5.1) and It is easy to check that when α 1 = 0 and α 2 = 0, λ = 0 is a root of det Δ * (λ) = 0 with multiplicity 2. From Theorem 4.3, we obtain that the equilibrium O(0, 0) of system (5.1) is a cusp of codimension 2. To determine that system (5.1) is the versal unfolding of system (4.1) with Bogdanov-Takens singularity, we consider the following suspension system In order to rewrite (5.4) as an abstract Cauchy problem, we set We consider the linear operator A : Then Since A is a Hille-Yosida operator, we can prove that A is also a Hille-Yosida operator. Consider F : D (A) → X the nonlinear map defined by Then we have Then So λ ∈ σ (A) .

Lemma 5.2. Let Assumption 4.1 be satisfied. Then we have
and Furthermore, Thus, for each λ > 0 large enough, It follows that T A0 (t) and S A (t) are defined, respectively, by (5.6) and (5.7). By using formula (5.6), we deduce that T A0 (t) ess = T A0 (t) ess , for all t ≥ 0, and it follows that This completes the proof.
Now we compute the projectors on the generalized eigenspace associated to eigenvalue 0 of A. Note that where S C (0, ε) + is the counterclockwise oriented circumference |λ| = ε for sufficiently small ε > 0 such that |λ| ≤ ε does not contain other point of the spectrum than 0. Since where ψ 1 and ψ 2 are defined in (4.5) and (4.6). Set Then The projector on the generalized eigenspace of A associated to 0 is given in the following lemma.
Lemma 5.3. 0 is a pole of order 2 of the resolvent of A, and the projector on the generalized eigenspace of A associated to 0 is given by From the above results, we obtain a state space decomposition with respect to the spectral properties of the linear operator A. More precisely, the projector on the linear center manifold is defined by We denote by We observe that for each By projecting on X c , we obtain precisely, on the curve SN , the model has only a unique equilibrium which is a saddle-node; when the parameters ( α 1 , α 2 ) lie between the curve SN and the curve H, the model has a saddle and a stable focus and no periodic orbit; on the HL curve, the model has an unstable focus and a stable homoclinic loop; the model has a saddle, an unstable focus, and no periodic orbits when the parameters lie between the curve HL and the curve SN ; however, when the parameters lie between the curve H and the curve HL, the model has a unique stable limit cycle (see Fig. 1).