BOGDANOV-TAKENS BIFURCATION IN PREDATOR-PREY SYSTEMS

In this paper, we consider Bogdanov-Takens bifurcation in two predator-prey systems. It is shown that in the full parameter space, BogdanovTalens bifurcation can be codimension 2, 3 or 4. First, the simplest normal form theory is applied to determine the codimension of the systems as well as the unfolding terms. Then, bifurcation analysis is carried out to describe the dynamical behaviour and bifurcation property of the systems around the critical point.

1. Introduction. Since Lotka and Volterra created the so-called Lotka-Volterra predator-prey model [18,22], more sophisticated predator-prey models have been developed, given in the following general form of ordinary differential equations [6], where g(x, K) is a continuous and differentiable function, describing the specific growth rate of the prey in the absence of predators. A specific form of g, the logistic growth g(x, K) = r(1 − x K ), is usually considered as a prototype, satisfying g(0, K) = r > 0, g(K, K) = 0, g x (K, K) < 0, g x (x, K) ≤ 0, and g K (x, K) > 0 for any x > 0.
The functional response p(x) of predators to the prey is a continuous and differentiable function, satisfying p(0) = 0, which describes the change in the density of the prey attacked per unit time per predator as the prey density changes. The functional response functions which have been extensively used in modeling population dynamics have the following three forms. (a) Lotka-Volterra type: p(x) = mx, where m > 0, which is a linear function with lim p(x) = mx 2 ax 2 +bx+1 , where m > 0, a > 0 and b is a constant. When b = 0, it is called the Holling type III response function. When b > −2 √ a (so that ax 2 +bx+1 > 0 and hence p(x) > 0), it is called the generalized Holling type III or sigmoidal functional response [1].
The function q(x) in system (1), on the other hand, describes how predator convert the consumed prey into the growth of predators, and the parameter c represents the efficiency of predators in converting consumed prey into their growth, while d denotes the predator mortality rate. In the literature, the following three typical forms of q(x) have been used. (i) q(x) = p(x), which is used in most classical predator-prey models. (ii) q( x y ) = p( x y ), which depends upon the ratio of the prey to predators, not the prey density. (iii) q( y x ), which is based on the ratio of the predators to their prey, and the typical second equation of (1) is given bẏ y = yq( y x ) = sy(1 − y hx ) while the first equation of (1) still takes p(x). Now, we combine the three types of function q(x) with the three types of function p(x) to obtain nine systems: A ii : ẋ = r x 1 − x K − m x y a+x , y = y m c x a+x − d ; (3) B iii : ẋ = r x 1 − x K − m x 2 y ax 2 +bxy+y 2 , y = y mcx 2 ax 2 +bxy+y 2 − d ; (7) In order to simplify the analysis on the above models, we first use state scaling and time scaling τ = r t to transform (2)-(10) into the following dimensionless systems (through x = KX with D = d r and S = s r ): System A ii : System A iii : System B i : System B ii : System C ii : We still use the dot to denote the differentiation with respect to the new time τ . Note that all the parameters should take positive values, except for B (or b) which may also take zero or negative values, provided B > −2 √ A (or b > −2 √ a). If we add a negative, constant term to the second equation of the above systems, which measures the rate of harvesting or removal [23] for the systems, we can study the general effect of harvesting on these models. This investigation is left for future work.
The dynamics and bifurcations of systems (2)-(6) (or (11)-(15)) have been studied completely in [14], except for the Bogdanov-Takens (BT) bifurcation for system (4) (or (13)). Similarly, in [15] the dynamics and bifurcations of system (7) (or (16)) have been completely analyzed, but the BT bifurcation. Thus, in this paper we will particularly investigate the BT bifurcations in systems (4) and (7) (or (13) and (16)). General theory of BT bifurcation for codimension two was established by Bogdanov [2] and Takens [21], which has been presented and discussed in many books (for example, see [4,10,16]). Recently, Han et al. revisited this problem and gave a fairly complete study [11]. For codimension three and higher were studied by Dumortier et al. [5] and Mardesic [19]. The cusp case was studied by Han [12] to obtain the normal form with unfolding up to any higher order. Suppose the normal form of BT bifurcation without unfoliding is given bẏ where k = 3(n−1) 2 (n ≥ 2). It has been shown in [12] that the general unfolding of (20) can be put in the form, However, in reality not all practical problems can have a standard ("perfect") unfolding due to limitation on the physical parameters. In the following sections, the BT bifurcation analysis is given to the two systems (13) and (16) by using the simplest normal form (SNF) theory (e.g., see [7,8,24]) and the parametric simplest normal form (PSNF) theory (e.g., see [9,26]). We will show that both the two systems can exhibit the cusp case of BT bifurcation, but the codimension for system (13) is not over three, while that for system (16) is not over four which does not have a standard unfolding.
2. BT bifurcation for system (13). In this section, we present BT bifurcation analysis for system (13). It is noted that global stability for predator-prey systems with p(x) = q(x) = mx n a+x n (m > 0) has been investigated by Cheng et al. [3], with p(x) = q(x) = x a+x 2 has been considered by Ruan and Xiao [20]. Also, global stability for a class of predator-prey systems, which consist of the first equation in (13) and the second equation in (15), has been studied in [13]. However, to our best knowledge, global stability and dynamics for system (13) have not been completely studied.
2.1. Equilibria and BT critical point. System (13) has three equilibria: where X 2 is determined from the following quadratic polynomial equation: The condition 0 < X 2 < 1 guarantees that the equilibrium E 2 is an interior (positive) equilibrium.
It is easy to show that E 0 is always a saddle. Define Then, when B > max{−2 √ A, −2}, E 1 is globally asymptotically stable (GAS) for C ∈ (0, C t ] and unstable for C > C t ; when −2 √ A < B ≤ −2, E 1 is GAS for C ∈ (0, C * ] and locally asymptotically stable (LAS) for C ∈ (C * , C t ), and unstable for C > C t . No Hopf bifurcation can occur from E 1 , but transcritical bifurcation can happen between E 1 and E 2 , which may be a forward bifurcation if B > max{−2 It has been shown in [14] that either none or two Hopf bifurcations can occur from E 2 when B > max{−2 √ A, −2}; while when −2 √ A < B ≤ −2, saddle-node bifurcation, none, one or two Hopf bifurcations may occur from E 2 . Thus, bistable or even tristable phenomenon (involving equilibria and stable limit cycles) may occur. In addition, for −2 √ A < B ≤ −2, when the saddle-node bifurcation and a Hopf bifurcation coincide, a BT bifurcation occurs. The condition for a BT bifurcation to occur can be found from the characteristic polynomial: where J 2 is the Jacobian of system (13) evaluated at the equilibrium E 2 , yielding The necessary and sufficient conditions for system (13) to have a BT bifurcation from E 2 are Tr(J 2 ) = det(J 2 ) = 0, which gives (note that X 2 = 1 is not a solution) In the following, we will apply the simplest normal form (SNF) theory (e.g., see [7,8,24]) to determine the codimension for BT bifurcation of system (13) and the parametric simplest normal form (PSNF) theory (e.g., see [9,26]) to obtain the normal form of system (13) with unfolding terms.

2.2.
The SNF and codimension of the BT bifurcation. In order to obtain the normal form for the BT bifurcation, we first need to determine the codimention of system (13). Define the BT critical point as (µ 1 , µ 2 ) = (0, 0), where . Now, we assume (µ 1 , µ 2 ) = (0, 0), and apply the simplest normal form theory [7,8,24] to determine the codimension of system (13). To achieve this, introducing the following transformation: into (13), we obtain the following system: Then, we expand the above system around (u, v) = (0, 0) and apply the SNF theory [8,24] to the resulting system to obtain the following results. First, consider C = 4D. We introduce the following transformation, (13) to obtain the SNF up to 2nd-order terms: Further, introducing the transformation, into the above equations, we obtain Next, consider C = 4D, for which the coefficient of the term x 1 x 2 in (32) becomes zero. Thus, we need to find higher-order terms in the SNF. To achieve this, applying the change of state variables: where a ij and b ij are functions in D, and the time rescaling, into (13) yields the SNF up to 7th-order terms: Finally, applying the transformation into the above system we obtain The following theorem directly follows from (32) and (31).
The BT bifurcation is codimension 2 if C = 4D, and codimesion 3 if C = 4D. No codimension higher than 3 can happen.
In the following two sections, we will use the PSNF theory to obtain the normal forms with unfolding terms up to 2nd-order terms for the codimension-2 BT bifurcation and up to 4th-order terms for the codimension-3 BT bifurcation, and give a summary on the bifurcation analysis. Note that in the literature usually a number of nonlinear transformations are used to obtain the normal form with unfolding terms (e.g., see [5,17]). However, here we apply the PSNF theory to obtain the normal form with unfolding terms via only one nonlinear transformation.
2.3. The PSNF of the codimension-2 BT bifurcation and bifurcation analysis. Now, suppose C = 4D. To obtain the normal form with unfolding, we introduce the parametric transformation, together with the change of state variables (28), into (13) to obtain Then, we expand the above system around the point (u, v, µ 1 , µ 2 ) = (0, 0, 0, 0) and apply the PSNF theory, with the change of state variables: and the parametric transformation: to obtain the standard normal form with unfolding: Finally, introducing the transformation into the above system yields the normal form with unfolding up to 2nd-order terms: Now, we use the normal form (39) to analyze the codimension-2 BT bifurcation. Note that the normal form (39) is in the standard form given in [10]. Thus, we follow the approach described in [10] to obtain the following theorem. (1) Saddle-node bifurcation occurs from the bifurcation curve: (2) Hopf bifurcations occur from the bifurcation curve: (3) Homoclinic orbits occur from the bifurcation curve: The above formulas for bifurcation curves can be expressed in terms of the original perturbation parameters µ 1 and µ 2 by using (38).

2.4.
The PSNF of the codimension-3 BT bifurcation and bifurcation analysis. To analyze the codimension-3 BT bifurcation under the condition C = 4D, at which E 2 becomes (X, Y ) = ( 1 2 , 1), we need to find the corresponding normal form with unfolding. To achieve this, let We denote µ = (µ 1 , µ 2 , µ 3 ) and then apply (40) together with into (13) to obtain the following system: Then, we expand the vector field of (42) around the point (u, v, µ) = (0, 0, 0) and employ the PSNF theory [9,26] to obtain the parametric normal form. More precisely, applying the following change of state variables: where U 1ijklm and U 2ijklm are functions in D, with the time rescaling: and then the parametrization: into (13) yields where β = (β 1 , β 2 , β 3 ). Finally, introducing the transformation into the above equations yields the normal form with unfolding up to 4th-order terms: It is easy to use (45) to verify that det ∂(µ1,µ2,µ3) ∂(β1,β2,β3) = 128 + O(β) = 0 for small β with D > 0, which shows that near the critical point (A, B, C) = (8, − 4, 4D), system (13) has the same bifurcation set with respect to µ as system (47) has with respect to β, up to a homeomorphism in the parameter space. Now based on the results of [5], we apply the method of normal forms and Abelian integral (or Melnikov function) to obtain the following result (proof is omitted due to page limit). (1) Saddle-node bifurcation occurs from the critical surface: (2) Hopf bifurcation occurs from the critical surface: (3) Homoclinic bifurcation occurs from the critical surface: (4) Generalized Hopf bifurcation occurs from the critical curve: Degenerate homoclinic bifurcation occurs from the critical curve: 3. BT bifurcation for system (16). In this section, we present BT bifurcation analysis for system (16).
3.1. Equilibria and BT critical point. System (16) has also the exact three equilibria given in (22), but now X 2 is determined from the following quadratic polynomial equation: Again, the condition 0 < X 2 < 1 guarantees that the equilibrium E 2 is an interior (positive) equilibrium. Define where It is straightforward to prove that E 0 is always a saddle. Moreover, it is easy to see that E 1 is a stable node for C < AD, and a saddle for C > AD. Furthermore, we can show that E 1 is GAS for C ≤ AD. There is no Hopf bifurcation from E 1 . For the stability of E 2 , the negative part (with X + 2 < 0 or X − 2 < 0) is a saddle (only mathematically meaningful), and the positive equilibrium E + 2 (with X + 2 > 0) is also a saddle; while the positive equilibrium E − 2 (with X − 2 > 0) may lose stability at Hopf critical points, whose number must be two if they exist. So the only possibility for a BT bifurcation to occur is from E 2 . The condition for a BT bifurcation to occur can be found from the characteristic polynomial (25) for which J 2 represents the Jacobian of system (16) evaluated at the equilibrium E 2 , with The necessary and sufficient conditions for system (16) to have a BT bifurcation from E 2 are Tr(J 2 ) = det(J 2 ) = 0. Solving the three equations: F 2 = Tr(J 2 ) = det(J 2 ) = 0, yields (noting that X 2 = 1 is not a solution) a family of biologically meaningful solutions for system (16) to have a BT bifurcation from E 2 , given by for which the equilibrium E 2 becomes the same one given in (27): E 2 : 1 2 , C 4D . Similarly, we will first apply the SNF theory to determine the codimension for the BT bifurcation of system (16) and then employ the PSNF theory to obtain the normal form of system (16) with unfolding terms.

The SNF and codimension of the BT bifurcation. Define
Then, the BT critical point is defined by (µ 1 , µ 2 ) = (0, 0). Now, we assume (µ 1 , µ 2 ) = (0, 0), and apply the SNF theory [7,8,24] to determine the codimension of system (16). To achieve this, introducing the same transformation (29) into (16) yields the following system: Then, we expand the above system around (u, v) = (0, 0) and apply the SNF theory [8,24] to the resulting system to obtain the following results. For C = 4D 1+2D , we apply the following transformation: into (16) to obtain the following SNF up to 2nd-order terms: Further, applying the transformation (31) into the above system we obtain If C = 4D 1+2D , we apply the change of state variables, where a ij and b ij are functions in C and D, and the time rescaling, , to obtain the SNF up to 7th-order terms as follows: Further, introducing the transformation (33) into the above equations we obtain Based on (55) and (56), we have the following theorem. Similarly, in the following two sections, we will use the PSNF theory to obtain the normal forms with unfolding terms for the BT bifurcations, and present a summary on the bifurcation analysis. More precisely, we will show that due to the limitation on the parameters, the codimension of Hopf bifurcation is 2, which implies that the possible condimension-4 BT bifurcation is degenerate. Again, we will obtain the normal forms with unfolding terms through only one nonlinear transformation.

The PSNF of the codimension-2 BT bifurcation and bifurcation anal-
1+2D . To obtain the normal form with unfolding up to 2nd-order terms, we first apply the parametric transformation (53) together with the change of state variables (29) into (16) and then expand the resulting equations around the point (u, v, µ 1 , µ 2 ) = (0, 0, 0, 0) to obtain the following equations: 3264 BING ZENG, SHENGFU DENG AND PEI YU Next, we apply the PSNF theory, with the change of state variables: and the parameter transformations: to obtain the standard normal form with unfolding: Finally, introducing the transformation: yields the normal form with unfolding up to 2nd-order terms: Again, the normal form (60) is in the standard form given in [10], which yields the following theorem. (1) Saddle-node bifurcation occurs from the bifurcation curve: (2) Hopf bifurcations occur from the bifurcation curve: (3) Homoclinic orbits occur from the bifurcation curve: 1+2D ), stable.
The above formulas for bifurcation curves can be expressed in terms of the original perturbation parameters µ 1 and µ 2 by using (59).
It is seen from (67) that the system has three unfolding parameters and the coefficient of the term x 3 1 x 2 contains the parameter β 3 . Since the coefficient of the term x 4 1 x 2 is equal to − 32 D 4 (1 + D)(1 + 2D) 2 = 0, it seems that the system may have codimension-4 BT bifurcation. However, we will show that due to limitation on the parameters, codimension-3 of Hopf bifurcation is not possible. To achieve this, first it is easy to use (67) to obtain the trace of the Jacobian of system (67), J 2 , given by Tr(J 2 ) = β 2 − (−β 1 ) under which the determinant of J 2 is det(J 2 ) = 2 −β 1 > 0, (β 1 < 0),