Hopf Bifurcation of a Differential-Algebraic Bioeconomic Model with Time Delay

Abstract

We investigate the dynamics of a differential-algebraic bioeconomic model with two time delays. Regarding time delay as a bifurcation parameter, we show that a sequence of Hopf bifurcations occur at the positive equilibrium as the delay increases. Using the theories of normal form and center manifold, we also give the explicit algorithm for determining the direction of the Hopf bifurcations and the stability of the bifurcating periodic solutions. Numerical tests are provided to verify our theoretical analysis.

1. Introduction

Predator-prey systems with delay play an important role in population dynamics since Vito Volterra and James Lotka proposed the seminal models of predator-prey in the mid of 1920s. In recent years, predator-prey systems with delay have been studied extensively due to their theoretical and practical significance.

In 1973, May [1] first proposed the delayed predator-prey system 𝑥̇=𝑥𝑟₁−𝑎₁₁𝑥(𝑡−𝜏)−𝑎₁₂𝑦,𝑦̇=𝑦−𝑟₂+𝑎₂₁𝑥−𝑎₂₂𝑦, (1.1) where 𝑥 and 𝑦 can be interpreted as the population densities of prey and predator at time 𝑡, respectively. 𝜏>0 is the feedback time delay of prey species to the growth of species itself. 𝑟₁>0 denotes intrinsic growth rate of prey, and 𝑟₂>0 denotes the death rate of predator. The parameters 𝑎_𝑖𝑗(𝑖,𝑗=1,2) are all positive constants.

Recently, Song and Wei [2] further considered the existence of local Hopf bifurcations of system (1.1). They regarded the feedback time delay 𝜏 as a bifurcation parameter and investigated the stability and the direction of periodic solutions bifurcating from Hopf bifurcations, by applying the normal form theory and the center manifold reduction developed by Hassard et al. [3].

Considering the feedback time delay of predator species to the growth of species itself and also with the delay 𝜏, Yan and Li [4] studied Hopf bifurcation and global periodic solutions of the following modified delayed predator-prey system: 𝑥̇=𝑥𝑟₁−𝑎₁₁𝑥(𝑡−𝜏)−𝑎₁₂𝑦,𝑦̇=𝑦−𝑟₂+𝑎₂₁𝑥−𝑎₂₂𝑦(𝑡−𝜏). (1.2)

They established the global existence results of periodic solutions bifurcating from Hopf bifurcations using a global Hopf bifurcation result due to Wu [5], which was different from that used in Song and Wei [2]. In [6], Yuan and Zhang also investigated the stability of the positive equilibrium and existence of Hopf bifurcation of the model (1.2).

Noting that, in real situations, the feedback time delay of the prey to the growth of the species itself and the feedback time delay of the predator to the growth of the species itself are different. If the time delay of the predator to the growth of the species itself is zero, (1.2) becomes (1.1). Meanwhile, if the time delay of the prey to the growth of the species itself is zero, (1.2) is simplified as 𝑥̇=𝑥(𝑡)𝑟₁−𝑎₁₁𝑥(𝑡)−𝑎₁₂𝑦(𝑡),𝑦̇=𝑦(𝑡)−𝑟₂+𝑎₂₁𝑥(𝑡)−𝑎₂₂𝑦(𝑡−𝜏). (1.3)

In 1954, Gordon [7] studied the effect of harvest effort on ecosystem from an economic perspective and proposed the following economic theory: NetEconomicRevenue(NER)=TotalRevenue(TR)−TotalCost(TC). (1.4)

This provides theoretical fundament for the establishment of bioeconomic systems by differential-algebraic equations (DAEs).

Let 𝐸 and 𝑦 represent the harvest effort and the density of harvested population, respectively. Then TC=𝑐𝐸 and TR=𝛼𝐸𝑦, where 𝛼 is unit price of harvested population and 𝑐 is the cost of harvest effort. Considering the economic interest 𝑚 of the harvest effort on the predator, we can establish the following algebraic equation: 𝐸𝛼𝑦−𝑐=𝑚. (1.5)

Recently, a class of differential-algebraic bioeconomic models were proposed and analyzed in [8–14]. For example, [8] discussed the problems of chaos and chaotic control for a differential-algebraic system. The existence and stability of equilibrium points, the sufficient conditions of existence for various bifurcations (transcritical bifurcation, singular induced bifurcation, and Hopf bifurcation) are invested in [9–14]. However, the stability and the direction of periodic solutions bifurcating from Hopf bifurcations have not been discussed there. Very recently, Zhang et al. [15] analyzed the stability and the Hopf bifurcation of a differential-algebraic bioeconomic system with a single time delay and Holling II type functional response and first investigated the direction of Hopf bifurcation and the stability of the bifurcating periodic solutions.

Based on the above economic theory and system (1.3), in this paper, we study the Hopf bifurcation of the following differential-algebraic bioeconomic system: 𝑥̇=𝑥𝑟₁−𝑎₁₁𝑥(𝑡)−𝑎₁₂𝑦,𝑦̇=𝑦−𝑟₂+𝑎₂₁𝑥−𝑎₂₂𝑦(𝑡−𝜏)−𝐸𝑦,0=𝐸𝛼𝑦−𝑐−𝑚. (1.6)

The rest of this paper is organized as follows. In Section 2, regarding time delay 𝜏 as bifurcation parameter, first, we study the stability of the equilibrium point of system based on the new normal form approach proposed by Chen et al. [16]. Then, following the normal form approach theory and the center manifold theory introduced by Hassard et al. [3], we compute the normal form for the Hopf bifurcation of system (1.6) and analyze the direction and stability of the bifurcating periodic orbits of the system. Numerical examples are given in Section 3 to illustrate our theoretical results.

2. Hopf Bifurcation

In this section, we analyze the stability and the Hopf bifurcation of the differential-algebraic bioeconomic system (1.6). Furthermore, we will derive the formula for determining the properties of the Hopf bifurcation.

2.1. Existence of Hopf Bifurcation

For 𝑚>0, the system (1.6) has positive equipment point 𝑃_∗=(𝑥_∗,𝑦_∗,𝐸_∗), where 𝑥_∗=(𝑟₁−𝑎₁₂𝑦_∗)/𝑎₁₁, 𝐸_∗=𝑚/(𝛼𝑦_∗−𝑐) and 𝑦_∗ is the positive root of the following equation: 𝛼𝑎₁₁𝑎₂₂+𝑎₁₂𝑎₂₁𝑦²−𝑐𝑎₁₁𝑎₂₂+𝑐𝑎₁₂𝑎₂₁+𝛼𝑟₁𝑎₂₁−𝛼𝑎₁₁𝑟₂𝑦+𝑚𝑎₁₁+𝑐𝑟₁𝑎₂₁−𝑐𝑎₁₁𝑟₂=0. (2.1)

In order to guarantee the existence of the positive equipment point 𝑃_∗, some restrictions must be satisfied also: 𝑐𝛼<𝑦_∗<𝑟₁𝑎₁₂,𝑐𝑎₁₁𝑎₂₂+𝑐𝑎₁₂𝑎₂₁+𝛼𝑟₁𝑎₂₁>𝛼𝑎₁₁𝑟₂, (2.2) and 𝑚𝑎₁₁+𝑐𝑟₁𝑎₂₁−𝑐𝑎₁₁𝑟₂<0.

In order to analyze the local stability of the positive equilibrium point, we first use the linear transformation (𝑥,𝑦,𝐸)^𝑇=𝑄⋅(𝑢,𝑣,𝐸)^𝑇, where 𝑄=1000100−𝛼𝐸_∗𝛼𝑦_∗1−𝑐. (2.3)

Then the system (1.6) yields 𝑢̇=𝑢𝑟₁−𝑎₁₁𝑢−𝑎₁₂𝑣,𝑣̇=𝑣−𝑟₂+𝑎₂₁𝑢−𝑎₂₂𝑣(𝑡−𝜏)−𝐸𝑣+𝛼𝐸_∗𝑣²𝛼𝑦_∗,−𝑐0=𝐸−𝛼𝐸_∗𝑣𝛼𝑦_∗)−𝑐(𝛼𝑣−𝑐−𝑚. (2.4)

According to the literature [16], considering the local parametrization of the third equation of the system (2.4), and we can obtain the following parametric system of it: 𝑥̇₁=𝑥₁+𝑢_∗𝑟₁−𝑎₁₁𝑥₁+𝑢_∗−𝑎₁₂𝑥₂+𝑣_∗,𝑥̇₂=𝑥₂+𝑣_∗−𝑟₂+𝑎₂₁𝑥₁+𝑢_∗−𝑎₂₂𝑥₂(𝑡−𝜏)+𝑣_∗−𝐸_∗+ℎ𝑥₂+𝑣_∗+𝛼𝐸_∗𝑥₂+𝑣_∗²𝛼𝑦_∗,−𝑐 (2.5) where 𝑢=𝑢_∗+𝑥₁,𝑣=𝑣_∗+𝑥₂,𝐸=𝐸_∗+ℎ𝑥₁,𝑥₂, (2.6)ℎ∶𝑅²→𝑅 is a smooth mapping, and (𝑢_∗,𝑣_∗,𝐸_∗)^𝑇 satisfies 𝑥_∗,𝑦_∗,𝐸_∗^𝑇=𝑄⋅𝑢_∗,𝑣_∗,𝐸_∗^𝑇. (2.7)

Thus the linearized system of parametric system (2.5) at (0,0) is given as below: 𝑥̇₁=−𝑎₁₁𝑥_∗𝑥₁−𝑎₁₂𝑥_∗𝑥₂,𝑥̇₂=𝑎₂₁𝑦_∗𝑥₁+𝛼𝐸_∗𝑦_∗𝛼𝑦_∗𝑥−𝑐²−𝑎₂₂𝑦_∗𝑥₂(𝑡−𝜏). (2.8)

The characteristic equation of the linearized system of parametric system (2.8) at (0,0) is given by 𝜆²+𝑝𝜆+𝑟+𝑠𝜆+𝑞𝑒^−𝜆𝜏=0, (2.9) where 𝑝=𝑎₁₁𝑥_∗−𝛼𝑦_∗𝐸_∗/(𝛼𝑦_∗−𝑐),𝑟=𝑥_∗𝑦_∗(𝑎₂₁𝑎₁₂−𝛼𝑎₁₁𝐸_∗/(𝛼𝑦_∗−𝑐)),𝑠=𝑎₂₂𝑦_∗, and 𝑞=𝑎₁₁𝑎₂₂𝑥_∗𝑦_∗.

The second-degree transcendental polynomial equation (2.9) has been extensively studied and applied by many researchers [2, 17–20].

Let(H1)𝑟₁>𝛼𝑦_∗𝐸_∗/(𝛼𝑦_∗−𝑐);(H2)𝑎₁₂𝑎₂₁+𝑎₁₁𝑎₂₂>𝛼𝑎₁₁𝐸_∗/(𝛼𝑦_∗−𝑐);(H3) either 2𝑥_∗𝑦_∗𝑎₂₁𝑎₁₂+𝑎²₂₂𝑦²_∗<𝑎²₁₁𝑥²_∗+(𝑦_∗𝐸_∗𝛼/(𝛼𝑦_∗−𝑐))² and (𝑎₂₁𝑎₁₂−𝑎₁₁𝐸_∗𝛼/(𝛼𝑦_∗−𝑐))²>(𝑎₁₁𝑎₂₂)² or 2𝑥_∗𝑦_∗𝑎₂₁𝑎₁₂+𝑎²₂₂𝑦²_∗−𝑎²₁₁𝑥²_∗−𝑦_∗𝐸_∗𝛼𝛼𝑦_∗−𝑐²²<4𝑥²_∗𝑦²_∗𝑎₂₁𝑎₁₂−𝑎₁₁𝐸_∗𝛼𝛼𝑦_∗−𝑐²−𝑎²₁₁𝑎²₂₂; (2.10)(H4) either (𝑎₂₁𝑎₁₂−𝑎₁₁𝐸_∗𝛼/(𝛼𝑦_∗−𝑐))²<(𝑎₁₁𝑎₂₂)² or 2𝑥_∗𝑦_∗𝑎₂₁𝑎₁₂+𝑎²₂₂𝑦²_∗>𝑎²₁₁𝑥²_∗+(𝑦_∗𝐸_∗𝛼/(𝛼𝑦_∗−𝑐))², and 2𝑥_∗𝑦_∗𝑎₂₁𝑎₁₂+𝑎²₂₂𝑦²_∗−𝑎²₁₁𝑥²_∗−𝑦_∗𝐸_∗𝛼𝛼𝑦_∗−𝑐²²=4𝑥²_∗𝑦²_∗𝑎₂₁𝑎₁₂−𝑎₁₁𝐸_∗𝛼𝛼𝑦_∗−𝑐²−𝑎²₁₁𝑎²₂₂; (2.11)(H5)(𝑎₂₁𝑎₁₂−𝑎₁₁𝐸_∗𝛼/(𝛼𝑦_∗−𝑐))²>(𝑎₁₁𝑎₂₂)², 2𝑥_∗𝑦_∗𝑎₂₁𝑎₁₂+𝑎²₂₂𝑦²_∗>𝑎²₁₁𝑥²_∗+(𝑦_∗𝐸_∗𝛼/(𝛼𝑦_∗−𝑐))², and 2𝑥_∗𝑦_∗𝑎₂₁𝑎₁₂+𝑎²₂₂𝑦²_∗−𝑎²₁₁𝑥²_∗−𝑦_∗𝐸_∗𝛼𝛼𝑦_∗−𝑐²²>4𝑥²_∗𝑦²_∗𝑎₂₁𝑎₁₂−𝑎₁₁𝐸_∗𝛼𝛼𝑦_∗−𝑐²−𝑎²₁₁𝑎²₂₂. (2.12)

From Lemma 2.1 in [2], we easily obtain the following results about the stability of the positive equilibrium and the Hopf bifurcation of (1.6).

Theorem 2.1. For system (1.6). One has the following.
(i) If (H1)–(H3) hold, then the equilibrium 𝐸_∗ of the system (1.6) is asymptotically stable for all 𝜏≥0.
(ii) If (H1), (H2), and (H4) hold, then the equilibrium 𝐸_∗ of the system (1.6) is asymptotically stable when 0≤𝜏<𝜏⁺₀ and unstable when 𝜏>𝜏⁺₀. System (1.6) undergoes a Hopf bifurcation at 𝐸_∗ when 𝜏=𝜏⁺_𝑗.
(iii) If (H1), (H2), and (H5) hold, then there is a positive integer 𝑘, such that the equilibrium 𝐸_∗ switches 𝑘 times from stability to instability to stability; that is, 𝐸_∗ is asymptotically stable when 𝜏∈0,𝜏⁺₀∪𝜏⁻₀,𝜏⁺₁∪⋯∪𝜏⁻_𝑘−1,𝜏⁺_𝑘 (2.13) and unstable when 𝜏∈𝜏⁺₀,𝜏⁻₀∪𝜏⁺₁,𝜏⁻₁∪⋯∪𝜏⁺_𝑘−1,𝜏⁻_𝑘−1∪𝜏⁺_𝑘−1,+∞. (2.14)
Here, 𝜔_±=√222𝑥_∗𝑦_∗𝑎₂₁𝑎₁₂+𝑎²₂₂𝑦²_∗−𝑎²₁₁𝑥²_∗−𝑦_∗𝐸_∗𝛼𝛼𝑦_∗−𝑐²±2𝑥_∗𝑦_∗𝑎₂₁𝑎₁₂+𝑎²₂₂𝑦²_∗−𝑎²₁₁𝑥²_∗−𝑦_∗𝐸_∗𝛼𝛼𝑦_∗−𝑐²²+4𝑥²_∗𝑦²_∗𝑎₂₁𝑎₁₂−𝑎₁₁𝐸_∗𝛼𝛼𝑦_∗−𝑐²−𝑎₁₁𝑎₂₂^21/2,𝜏^±_𝑗=1𝜔_±𝑎arccos²²𝐸_∗𝛼/𝛼𝑦_∗−𝑐𝜔²_±−𝑎₁₁𝑎₂₂𝑥²_∗𝑎₂₁𝑎₁₂−𝑎₁₁𝐸_∗𝛼/𝛼𝑦_∗−𝑐𝑎²₂₂𝜔²_±+𝑎²₁₁𝑥²_∗+2𝑗𝜋𝜔_±,𝑗=0,1,2,…. (2.15)

Denote 𝜆_𝑗(𝜏)=𝛼_𝑗(𝜏)+𝑖𝜔_𝑗(𝜏),𝑗=0,1,2,…; then the following transversality conditions hold: 𝑑𝑑𝜏Re𝜆_𝑗𝜏⁺_𝑗𝑑>0,𝑑𝜏Re𝜆_𝑗𝜏⁻_𝑗<0. (2.16)

2.2. Direction and Stability of the Hopf Bifurcation

In the previous section, we obtain the conditions for Hopf bifurcations at the positive equilibrium 𝑃_∗ when 𝜏=𝜏^±_𝑛,𝑛=0,1,2,…. In this section, we will derive the formulaes determining the direction, stability, and period of these periodic solutions bifurcating from 𝑃_∗ at 𝜏. The basic techniques used here are the normal form and the center manifold theory developed by Hassard et al. [3].

Assuming the system (1.6) undergoes Hopf bifurcations at the positive equilibrium 𝑃_∗ at 𝜏=𝜏_𝑛, we let 𝑖𝜔 be the corresponding purely imaginary root of the characteristic equation at the positive equilibrium. In order to investigate the direction of Hopf bifurcation and the stability of the bifurcating periodic solutions of system (1.6), we consider the parametric system (2.5) of system (2.4).

Let 𝑡=𝑡𝜏,𝜏=𝜏_𝑛+𝜇; then the parametric system (2.5) is equivalent to the following functional differential equation system: 𝑧̇₁=𝜏_𝑛+𝜇−𝑎₁₁𝑥_∗𝑧₁−𝑎₁₂𝑥_∗𝑧₂,𝑧̇₂=𝜏_𝑛+𝜇𝑎₂₁𝑦_∗𝑧₁+𝛼𝐸_∗𝑦_∗𝛼𝑦_∗𝑧−𝑐²−𝑎₂₂𝑦_∗𝑧₂(𝑡−1), (2.17) where 𝑧₁=𝑥₁(𝑡𝜏),𝑧₂=𝑥₂(𝑡𝜏).

For Φ∈𝐶([−1,0],𝑅²), define a family of operators 𝐿_𝜇Φ∶=𝜏_𝑛+𝜇−𝑎₁₁𝑥_∗−𝑎₁₂𝑥_∗𝑎₂₁𝑦_∗𝛼𝐸_∗𝑦_∗𝛼𝑦_∗Φ(0)+𝜏−𝑐^𝑛+𝜇000−𝑎₂₂𝑦_∗Φ(),𝑓𝜏−1𝜇,Φ∶=_𝑛+𝜇×−𝑎₁₁Φ²₁(0)−𝑎₁₂Φ₁(0)Φ₂(0)𝑎+⋯₂₁Φ₁(0)Φ₂(0)−𝑎₂₂Φ₂(0)Φ₂(−1)−𝛼𝐸_∗𝑐𝛼𝑦_∗−𝑐²Φ²₂(0),+⋯ (2.18) where Φ=(Φ₁,Φ₂)∈𝐶. By the Riesz representation theorem, there exists a function of bounded variation for 𝜃∈[−1,0], such that 𝐿_𝜇Φ=⁰₋₁𝑑𝜂𝜃,𝜇Ψ(𝜃),Ψ∈𝐶. (2.19)

In fact, we can choose 𝜂𝜃,𝜇=𝜏_𝑛+𝜇−𝑎₁₁𝑥_∗−𝑎₁₂𝑥_∗𝑎₂₁𝑦_∗𝛼𝐸_∗𝑦_∗𝛼𝑦_∗𝛿(𝜃)+𝜏−𝑐^𝑛+𝜇000𝑎₂₂𝑦_∗𝛿(𝜃+1). (2.20)

For Φ∈𝐶, define 𝐴𝜇(𝜃)Φ∶=𝑑Φ[),∫𝑑𝜃,𝜃∈−1,0₀₋₁𝑅𝜇[),𝑓()𝑑𝜂𝑡,𝜇Φ,𝜃=0,∶=0,𝜃∈−1,0𝑡,Φ,𝜃=0. (2.21)

Hence system (2.17) can be rewritten as 𝑍̇_𝑡=𝐴𝜇𝑍_𝑡+𝑅𝜇𝑍_𝑡. (2.22)

For Ψ∈𝐶([−1,0],(𝑅²)^∗), define 𝐴^∗Ψ(𝑠)(𝑠)∶=𝑑Ψ[),∫𝑑𝑠,𝑠∈−1,0₀₋₁𝑑𝜂^𝑇(𝑠,0)Ψ(−𝑠),𝑠=0, (2.23) and a bilinear inner product ⟨Ψ,Φ⟩=Ψ^𝑇(0)Φ(0)−⁰₋₁^𝜃₀Ψ^𝑇(𝜉−𝜃)𝑑𝜂(𝜃)Φ(𝜉)𝑑𝜉, (2.24) where 𝜂(𝜃)=𝜂(𝜃,0). Then 𝐴^∗ and 𝐴(0) are adjoint operators. By the discussion in Section 2.1, we know that ±𝑖𝜔 are eigenvalues of 𝐴(0). Thus they are also eigenvalues of 𝐴^∗. Suppose that 𝑞(𝜃)=(1,𝑞₁)^𝑇𝑒^{𝑖𝜔𝜏𝜃} and 𝑞^∗(𝑠)=𝐷(𝑞^∗₁,1)𝑒^{𝑖𝜔𝜏𝑠} are eigenvectors of 𝐴(0) and 𝐴^∗ corresponding to 𝑖𝜔𝜏 and −𝑖𝜔𝜏, respectively. It is easy to obtain 𝑞₁=−𝑖𝜔+𝑎₁₁𝑥_∗𝑎₁₂𝑥_∗,𝑞^∗₁=𝑖𝜔−𝑎₂₂𝑦_∗𝑒^𝑖𝜔𝜏𝑎₁₂𝑥_∗,𝐷=𝑞₁+𝑞^∗₁+𝜏𝑞₁𝑎₂₂𝑦_∗𝑒^{−𝑖𝜔𝜏−1}. (2.25) Moreover, ⟨𝑞^∗(𝑠),𝑞(𝜃)⟩=1,⟨𝑞^∗(𝑠),𝑞(𝜃)⟩=0.

Now, we calculate the coordinates to describe the center manifold 𝐶₀ at 𝜇=0. Using the same notations as in [3], we define 𝑧=𝑞^∗,𝑍_𝑡,𝑊(𝑡,𝜃)=𝑍_𝑡(𝜃)−2Re𝑧𝑞(𝜃). (2.26)

On the center manifold 𝐶₀, we have 𝑊(𝑡,𝜃)=𝑊(𝑧,𝑧,𝜃), where 𝑊(𝑧,𝑧,𝜃)=𝑊₂₀(𝜃)𝑧²2+𝑊₁₁(𝜃)𝑧𝑧+𝑊₀₂(𝜃)𝑧²2+⋯, (2.27) and 𝑧 and 𝑧 are local coordinates for 𝐶₀ in the direction of 𝑞^∗ and 𝑞^∗. Note that 𝑊 is real if 𝑍_𝑡 is real; we only consider real solutions. For 𝑍_𝑡∈𝐶₀, we have 𝑧̇=𝑖𝜔𝜏𝑧+𝑞^∗(0)𝑓₀(𝑧,𝑧)∶=𝑖𝜔𝜏𝑧+𝑔(𝑧,𝑧), (2.28) where 𝑔(𝑧,𝑧)=𝑞^∗(0)𝑓₀(𝑧,𝑧)=𝑞^∗(0)𝑓0,𝑍_𝑡=𝑔₂₀(𝜃)𝑧²2+𝑔₁₁(𝜃)𝑧𝑧+𝑔₀₂(𝜃)𝑧²2+⋯. (2.29)

Moreover, we have 𝑔₂₀=2𝐷𝜏𝑎₂₁𝑞₁−𝑎₁₁𝑞^∗₁−𝑎₁₂𝑞^∗₁𝑞₁−𝛼𝐸_∗𝑐𝛼𝑦_∗−𝑐²𝑞²₁−𝑎₂₂𝑞²₁𝑒^{−𝑖𝜔𝜏𝜃},𝑔₁₁=2𝐷𝜏𝑎₂₁Re𝑞₁−𝑎₁₁𝑞^∗₁−𝑎₁₂𝑞^∗₁Re𝑞₁−𝑎₂₂Re𝑞₁𝑞₁𝑒^{𝑖𝜔𝜏𝜃}−𝛼𝐸_∗𝑐𝛼𝑦_∗−𝑐²𝑞₁𝑞₁,𝑔₀₂=2𝐷𝜏𝑎₂₁𝑞₁−𝑎₁₁𝑞^∗₁−𝑎₁₂𝑞^∗₁𝑞₁−𝛼𝐸_∗𝑐𝛼𝑦_∗−𝑐²𝑞²₁−𝑎₂₂𝑞²₁𝑒^{𝑖𝜔𝜏𝜃},𝑔₂₁=2𝐷𝜏𝑎₂₁𝑞₁−2𝑎₁₁𝑞^∗₁−𝑎₁₂𝑞^∗₁𝑞₁𝑊⁽¹⁾₁₁(0)+12𝑎₂₁𝑞₁−𝑎₁₂𝑞^∗₁𝑞₁−2𝑎₁₁𝑞^∗₁𝑊⁽¹⁾₂₀(0)+𝑎₂₁−𝑎₂₂𝑞₁𝑒^{−𝑖𝜔𝜏𝜃}−𝑎₁₂𝑞^∗₁−2𝛼𝐸_∗𝑐𝛼𝑦_∗−𝑐²𝑞₁𝑊⁽²⁾₁₁(0)+12𝑎₂₁−𝑎₂₂𝑞₁𝑒^{𝑖𝜔𝜏𝜃}−𝑎₁₂𝑞^∗₁−2𝛼𝐸_∗𝑐𝛼𝑦_∗−𝑐²𝑞₁𝑊⁽²⁾₂₀(0)−𝑎₂₂𝑞₁𝑊⁽²⁾₁₁()−1−12𝑎₂₂𝑞₁𝑊⁽²⁾₂₀(−1). (2.30)

Next, we will calculate 𝑊₁₁(𝜃) and 𝑊₂₀(𝜃). From (2.22) and (2.28), we have 𝑊̇=𝑍̇_𝑡−𝑧𝑞̇−𝑧̇𝑞=AW−2Re𝑞^∗(0)𝑓(𝑧,𝑧)𝑞(𝜃),𝜃∈[−1,0),AW−2Re𝑞^∗(0)𝑓(𝑧,𝑧)𝑞(𝜃()+𝑓,𝜃=0,∶=AW+𝐻𝑧,𝑧,𝜃), (2.31) where 𝐻(𝑧,𝑧,𝜃)=𝐻₂₀(𝜃)(𝑧²/2)+𝐻₁₁(𝜃)𝑧𝑧+𝐻₀₂(𝜃)(𝑧²/2)+⋯.

For 𝜃∈[−1,0), we can get (𝐴−2𝑖𝜔𝜏)𝑊₂₀(𝜃)=−𝐻₂₀(𝜃),AW₁₁(𝜃)=−𝐻₁₁(𝜃). (2.32)

On the other hand, 𝐻(𝑧,𝑧,𝜃)=AW−2Re{𝑞∗(0)𝑓(𝑧,𝑧)𝑞(𝜃)}; hence we can obtain 𝐻₂₀(𝜃)=−𝑔₂₀𝑞(𝜃)−𝑔₀₂𝑞(𝜃),𝐻₁₁(𝜃)=−𝑔₁₁𝑞(𝜃)−𝑔₁₁𝑞(𝜃). (2.33)

It follows from (2.32) that 𝑊̇₂₀(𝜃)=2𝑖𝜔𝑊₂₀(𝜃)+𝑔₂₀𝑞(𝜃)+𝑔₀₂𝑞(𝜃),𝑊̇₁₁(𝜃)=𝑔₁₁𝑞(𝜃)+𝑔₁₁𝑞(𝜃). (2.34)

Solving it, we have 𝑊₂₀(𝜃)=𝑖𝑔₂₀𝑞(0)𝑒𝜔𝜏^{𝑖𝜔𝜏𝜃}+𝑖𝑔₀₂3𝜔𝜏𝑞(0)𝑒^{−𝑖𝜔𝜏𝜃}+𝑀𝑒^{2𝑖𝜔𝜏𝜃},𝑊₁₁(𝜃)=−𝑖𝑔₁₁𝑞(0)𝑒𝜔𝜏^{𝑖𝜔𝜏𝜃}+𝑖𝑔₁₁𝜔𝜏𝑞(0)𝑒^{−𝑖𝜔𝜏𝜃}+𝑁. (2.35)

Now we will seek appropriate 𝑀 and 𝑁. From (2.29) and (2.31), we have 𝐻₂₀(𝜃)=−𝑔₂₀𝑞(0)−𝑔₀₂𝑞(0)+2𝜏𝐻₁,𝐻₁₁(𝜃)=−𝑔₁₁𝑞(0)−𝑔₁₁𝑞(0)+2𝜏𝐻₂, (2.36) where 𝐻₁=−𝑎₁₁−𝑎₁₂𝑞₁𝑎₂₁𝑞₁−𝛼𝐸_∗𝑐𝛼𝑦_∗−𝑐²𝑞²₁−𝑎₂₂𝑞²₁𝑒^{−𝑖𝜔𝜏𝜃},𝐻₂=−𝑎₁₁−𝑎₁₂Re𝑞₁𝑎₂₁Re𝑞₁−𝑎₂₂Re𝑞₁𝑞₁𝑒^{𝑖𝜔𝜏𝜃}−𝛼𝐸_∗𝑐𝛼𝑦_∗−𝑐²𝑞₁𝑞₁. (2.37)

Noting that 𝑖𝜔𝜏𝐼−⁰₋₁𝑒^{𝑖𝜔𝜏𝜃}𝑑𝜂(𝜃)𝑞(0)=0,−𝑖𝜔𝜏𝐼−⁰₋₁𝑒^{−𝑖𝜔𝜏𝜃}𝑑𝜂(𝜃)𝑞(0)=0, (2.38) we have the following linear equations: 2𝑖𝜔+𝑎₁₁𝑥_∗𝑎₁₂𝑥_∗−𝑎₂₁𝑦_∗2𝑖𝜔−𝑎₂₂𝑦_∗−𝛼𝐸_∗𝑦_∗𝛼𝑦_∗−𝑐𝑀=2𝐻₁,𝑎₁₁𝑥_∗𝑎₁₂𝑥_∗−𝑎₂₁𝑦_∗−𝑎₂₂𝑦_∗−𝛼𝐸_∗𝑦_∗𝛼𝑦_∗−𝑐𝑁=2𝐻₂. (2.39)

It is easy to get 𝑀 and 𝑁. Furthermore, we can get the following values: 𝐶₁(0)=𝑖𝑔2𝜔𝜏¹¹𝑔₂₀−2𝑔₁₁²−𝑔₀₂²3+𝑔₂₁2,𝜇₂𝐶=−Re₁(0)Re𝜆^′(𝜏),𝛽₂=2Re𝐶₁(0),𝑇₂𝐶=−Im₁(0)+𝜇₂Im𝜆^′(𝜏).𝜔𝜏 (2.40)

It is well known that, at the critical value 𝜏, the sign of 𝜇₂ determines the direction of the Hopf bifurcation, 𝛽₂ determines the stability of bifurcated periodic solutions, while 𝑇₂ determines the period of bifurcated periodic solutions, respectively. In fact, we have the theorem as follows.

Theorem 2.2. (i) The Hopf bifurcation is supercritical if 𝜇₂>0 and subcritical if 𝜇₂<0.
(ii) The bifurcated periodic solutions are unstable if 𝛽₂>0 and stable if 𝛽₂<0.
(iii) The period of bifurcated periodic solutions increases if 𝑇₂>0 and decreases if 𝑇₂<0.

3. Numerical Examples

We now perform some simulations for better understanding of our analytical treatment.

Example 3.1. Consider the following differential-algebraic system: ̇(𝑡),̇()𝑥=𝑥6−2𝑥−1.5𝑦𝑦=𝑦−0.5+2𝑥−𝑦𝑡−𝜏−𝐸𝑦,0=𝐸2𝑦−1−1.5. (3.1)
The only positive equilibrium point of the system (3.1) is 𝑃⁺=(1.5,2.0,0.5). It is easy to get 𝑝=2.3333,𝑞=6,𝑟=7,𝑠=2.
When 𝜏=0, from (2.9), the characteristic equation of the linearized system of parametric system of (3.1) at 𝑃⁺ is given by 𝜆²+𝑝+𝑠𝜆+𝑟+𝑞=0. (3.2) By simple calculating, 𝑝+𝑠=4.3333>0,𝑞+𝑟=13>0. According to Routh-Hurwitz stability criterion, 𝑃⁺ is stable.
When 𝜏>0, we can further get 𝑟²−𝑞²=13>0,𝑠²−𝑝²+2𝑟=12.5556>0, and (𝑠²−𝑝²+2𝑟)²−4(𝑟²−𝑞²)=105.6420>0 which mean (H1), (H2), and (H5) hold. Furthermore, we have 𝜔₋=1.0671,𝜔₊=3.3789, and hence 𝜏⁺₀=0.5638<𝜏⁺₁=2.4233<𝜏⁻₀=2.8880<𝜏⁺₂=4.2829<𝜏⁺₃=6.1424<𝜏⁻₁=8.7762<𝜏⁻₂𝜏=14.6645<₋₃=20.5527<⋯. According to conclusion (iii) of Theorem 2.1, the equilibrium point 𝑃⁺ is stable only for 𝜏<𝜏⁺₀, unstable for any 𝜏>𝜏⁺₀, and the Hopf bifurcation occurs at 𝜏=𝜏⁺₀.
By the aid of Mathematica, we can obtain the following values according to equalities in (2.40): 𝐶₁(0)=3.9725−4.7012𝑖,𝜆^′𝜏⁺₀=1.3427−2.9776𝑖. (3.3) So we have 𝜇₂=−2.9586<0,𝛽₂=7.9451>0, and 𝑇₂=−2.1567<0. Thus from Theorem 2.2, we can conclude that the Hopf bifurcation of system (3.1) is subcritical, the bifurcating periodic solution exists when 𝜏 crosses 𝜏⁺₀ to the left, and the bifurcating periodic solution is unstable and decreases.
Figure 1 shows the positive equilibrium point 𝑃⁺ of system (3.1) is locally asymptotically stable when 𝜏=0.55<𝜏⁺₀=0.5638. The periodic solutions occur from 𝑃⁺ when 𝜏=0.563<𝜏⁺₀=0.5638 as is illustrated in Figure 2. When 𝜏=0.57>𝜏⁺₀=0.5638, the positive equilibrium point 𝑃⁺ becomes unstable as shown in Figure 3.

Figure 1 

: Dynamic behavior of the differential-algebraic system (3.1) with 𝜏=0.55<𝜏+0, and the positive equilibrium point 𝑃+ is stable.

Figure 2 

Dynamic behavior of the differential-algebraic system (3.1) with 𝜏=0.563<𝜏+0, periodic solutions bifurcating from the positive equilibrium point 𝑃+.

Figure 3 

Dynamic behavior of the differential-algebraic system (3.1) with 𝜏=0.57>𝜏+0, and the positive equilibrium point 𝑃+ is unstable.

Acknowledgments

The authors would like to thank the reviewers for their valuable comments which improved the paper. This work is supported by the Foundation for the Authors of the National Excellent Doctoral Thesis Award of China (200720), the “333 Project” Foundation of Jiangsu Province, and the Foundation for Innovative Program of Jiangsu Province(CXZZ11_0870, CXZZ12_0383).