Bifurcation Analysis in a Kind of Fourth-Order Delay Differential Equation

Abstract

A kind of fourth-order delay differential equation is considered. Firstly, the linear stability is investigated by analyzing the associated characteristic equation. It is found that there are stability switches for time delay and Hopf bifurcations when time delay cross through some critical values. Then the direction and stability of the Hopf bifurcation are determined, using the normal form method and the center manifold theorem. Finally, some numerical simulations are carried out to illustrate the analytic results.

1. Introduction

Sadek [1] has considered the following fourth-order delay differential equation:𝑥⁽⁴⁾(𝑡)+𝛼₁𝑥⃛(𝑡)+𝛼₂𝑥̈(𝑡)+𝜙(𝑥̇(𝑡−𝜏))+𝑓(𝑥(𝑡))=0. (1.1)By constructing Lyapunov functionals, it was given a group of conditions to ensure that the zero solution of (1.1) is globally asymptotically stable when the delay 𝜏 is suitable small, but if the sufficient conditions are not satisfied, what are the behaviors of the solutions? This is a interesting question in mathematics. The purpose of the present paper is to study the dynamics of (1.1) from bifurcation. We will give a detailed analysis on the above mentioned question. By regarding the delay 𝜏 as a bifurcation parameter, we analyze the distribution of the roots of the characteristic equation of (1.1) and obtain the existence of stability switches and Hopf bifurcation when the delay varies. Then by using the center manifold theory and normal form method, we derive an explicit algorithm for determining the direction of the Hopf bifurcation and the stability of the bifurcating periodic solutions.

We would like to mention that there are several articles on the stability of fourth-order delay differential equations, we refer the readers to [1–8] and the references cited therein.

The rest of this paper is organized as follows. In Section 2, we firstly focus mainly on the local stability of the zero solution. This analysis is performed through the study of a characteristic equation, which takes the form of a fourth-degree exponential polynomial. Using the approach of Ruan and Wei [9], we show that the stability of the zero solution can be destroyed through a Hopf bifurcation. In Section 3, we investigate the stability and direction of bifurcating periodic solutions by using the normal form theory and center manifold theorem presented in Hassard et al. [10]. In Section 4, we illustrate our results by numerical simulations. Section 5 with conclusion completes the paper.

2. Stability and Hopf Bifurcation

In this section, we will study the stability of the zero solution and the existence of Hopf bifurcation by analyzing the distribution of the eigenvalues. For convenience, we give the following assumptions:𝜏>0,𝛼₁>0,𝛼₂>0,𝜙(0)=0,𝑓(0)=0, (𝐇_𝟏)with 𝜙 and 𝑓 are both continuous functions and those three-order differential quotients at origin are existent. We rewrite (1.1) as the following form:̇̇̇̇𝑥=𝑦,𝑦=𝑢,𝑢=𝑣,𝑣=−𝛼₂𝑢−𝛼₁𝑣−𝑓(𝑥)−𝜙(𝑦(𝑡−𝜏)). (2.1)It is easy to see that (0,0,0,0) is the only trivial solution to the system (2.1) and the linearization around (0,0,0,0) is given bẏ̇̇̇𝑥=𝑦,𝑦=𝑢,𝑢=𝑣,𝑣=−𝑓^′(0)𝑥−𝛼₂𝑢−𝛼₁𝑣−𝜙^′(0)𝑦(𝑡−𝜏). (2.2)Its characteristic equation is𝜆⁴+𝛼₁𝜆³+𝛼₂𝜆²+𝜙^′(0)𝜆𝑒^−𝜆𝜏+𝑓^′(0)=0. (2.3)

Lemma 2.1. Suppose (𝐇_𝟏) and𝛼₁𝛼₂−𝜙^′(0)>0,𝑓^′(0)>0,𝜙^′(0)(𝛼₁𝛼₂−𝜙^′(0))−𝛼²₁𝑓^′(0)>0 (𝐇_𝟐)are satisfied. Then the trivial solution (0,0,0,0) is asymptotically stable when 𝜏=0.

Proof. When 𝜏=0, (2.3) becomes𝜆⁴+𝛼₁𝜆³+𝛼₂𝜆²+𝜙^′(0)𝜆+𝑓^′(0)=0. (2.4)By Routh-Hurwitz criterion, all roots of (2.4) have negative real parts if and only if𝛼₁>0,𝛼₁𝛼₂−𝜙^′(0)>0,𝑓^′(0)>0,𝜙^′(0)(𝛼₁𝛼₂−𝜙^′(0))−𝛼²₁𝑓^′(0)>0. (2.5)The conclusion follows from (𝐇_𝟏) and (𝐇_𝟐).
Let 𝑖𝜔(𝜔>0) be a root of (2.3), then we have𝜔⁴−𝑖𝛼₁𝜔³−𝛼₂𝜔²+𝜙^′(0)(𝑖𝜔)𝑒^{−𝑖𝜔𝜏}+𝑓^′(0)=0. (2.6)Separating the real and imaginary parts gives−𝜔⁴+𝛼₂𝜔²−𝑓^′(0)=𝜙^′𝛼(0)𝜔sin𝜔𝜏,₁𝜔³=𝜙^′(0)𝜔cos𝜔𝜏. (2.7)Adding up the squares of both equations yields𝜔⁸+(𝛼²₁−2𝛼₂)𝜔⁶+(𝛼²₂+2𝑓^′(0))𝜔⁴−(𝜙^′2(0)+2𝛼₂𝑓^′(0))𝜔²+𝑓^′2(0)=0. (2.8)Let 𝑉=𝜔², and denote𝑃=𝛼²₁−2𝛼₂,𝑄=𝛼²₂+2𝑓^′(0),𝐾=−𝜙^′2(0)−2𝛼₂𝑓^′(0). (2.9)Then (2.8) becomes𝑉⁴+𝑃𝑉³+𝑄𝑉²+𝐾𝑉+𝑓^′2(0)=0. (2.10)Setℎ(𝑉)=𝑉⁴+𝑃𝑉³+𝑄𝑉²+𝐾𝑉+𝑓^′2(0). (2.11)Then we haveℎ^′(𝑉)=4𝑉³+3𝑃𝑉²+2𝑄𝑉+𝐾. (2.12)Consider4𝑉³+3𝑃𝑉²+2𝑄𝑉+𝐾=0. (2.13)Let 𝑈=𝑉+(3/4)𝑃. Then (2.13) becomes𝑈³+𝑃₁𝑈+𝑄₁=0, (2.14)where𝑃₁=𝑄2−3𝑃16²,𝑄₁=1𝑃32³−18𝑃𝑄+𝐾. (2.15)Define𝑄𝑀=(¹2)²𝑃+(¹3)³,√𝜎=−1+3𝑖2,𝑈₁=³−𝑄₁2+√𝑀+³−𝑄₁2−√𝑈𝑀,₂=³−𝑄₁2+√𝑀𝜎+³−𝑄₁2−√𝑀𝜎²,𝑈₂=³−𝑄₁2+√𝑀𝜎²+³−𝑄₁2−√𝑉𝑀𝜎,_𝑖=𝑈_𝑖−34𝑃,𝑖=1,2,3. (2.16)Then by Lemma 2.2 in Li and Wei [11], we have the following results on the distribution of the roots of (2.10).

Lemma 2.2. (i) If 𝑀≥0, then (2.10) has positive roots if and only if 𝑉₁>0 and ℎ(𝑉₁)<0.
(ii) If 𝑀<0, then (2.10) has positive roots if and only if there exists at least one 𝑉^∗∈{𝑉₁,𝑉₂,𝑉₃}, such that 𝑉^∗>0 and ℎ(𝑉^∗)≤0.

Without loss of generality, we assume that equation ℎ(𝑉)=0 has four positive roots denoted by 𝑉₁,𝑉₂,𝑉₃, and 𝑉₄, respectively. Then (2.8) also has four positive roots, say 𝜔_𝑖=𝑉_𝑖,𝑖=1,2,3,4.

From (2.7), and conditions (𝐇_𝟏) and (𝐇_𝟐), we have that𝛼cos𝜔𝜏=¹𝜔²𝜙^′(0)>0. (2.17)Hence, we define𝜏^𝑗_𝑘=1𝜔_𝑘𝛼[arccos¹𝜔²_𝑘𝜙^′(0)+2𝑗𝜋],𝑘=1,2,3,4,𝑗=0,1,…, (2.18)when−𝜔⁴_𝑘+𝛼₂𝜔²_𝑘−𝑓^′(0)𝜙^′(0)𝜔_𝑘𝜏>0,_𝑗𝑘=1𝜔_𝑘𝛼[−arccos¹𝜔²_𝑘𝜙^′(0)+2(𝑗+1)𝜋],𝑘=1,2,3,4,𝑗=0,1,…, (2.19)when−𝜔⁴_𝑘+𝛼₂𝜔²_𝑘−𝑓^′(0)𝜙^′(0)𝜔_𝑘<0. (2.20)Let𝜆(𝜏)=𝛼(𝜏)+𝑖𝛽(𝜏) (2.21)be the root of (2.3) satisfying 𝛼(𝜏^𝑗_𝑘)=0,𝛽(𝜏^𝑗_𝑘)=𝜔_𝑘.

Lemma 2.3. Suppose ℎ^′(𝑉_𝑖)≠0(𝑖=1,2,3,4). If 𝜏=𝜏^𝑗_𝑘, then ±𝑖𝜔_𝑘 is a pair of simple purely imaginary roots of (2.3); and Re(𝑑𝜆(𝜏^𝑗_𝑘)/𝑑𝜏)>0 when 𝑘=2,4; and Re(𝑑𝜆(𝜏^𝑗_𝑘)/𝑑𝜏)<0 when 𝑘=1,3.

Proof. Substituting 𝜆(𝜏) into (2.3) and differentiating with respet to 𝜏 gives𝑑𝜆=𝜙𝑑𝜏^′(0)𝜆²𝑒^−𝜆𝜏4𝜆³+3𝛼₁𝜆²+2𝛼₂𝜆+𝜙^′(0)𝑒^−𝜆𝜏−𝜏𝜙^′(0)𝜆𝑒^−𝜆𝜏=−𝜆²(𝜆⁴+𝛼₁𝜆³+𝛼₂𝜆²+𝑓^′(0))𝜏𝜆⁵+(3+𝜏𝛼₁)𝜆⁴+(3𝛼₁+𝜏𝛼₂)𝜆³+𝛼₂𝜆²+𝑓^′(0)𝜏𝜆−𝑓^′.(0) (2.22)ThenRe𝑑𝜆(𝜏^𝑗)=𝑑𝜏(𝜔⁶_𝑘−𝛼₂𝜔⁴_𝑘+𝑓^′(0)𝜔²_𝑘)((3+𝜏𝛼₁)𝜔⁴_𝑘−𝛼₂𝜔²_𝑘−𝑓^′(0))−𝛼₁𝜔⁵_𝑘(𝜏𝜔⁵_𝑘−(3𝛼₁+𝜏𝛼₂)𝜔³_𝑘+𝑓^′(0)𝜏𝜔_𝑘)((3+𝜏𝛼₁)𝜔⁴_𝑘−𝛼₂𝜔²_𝑘−𝑓^′(0))²+(𝜏𝜔⁵_𝑘−(3𝛼₁+𝜏𝛼₂)𝜔³_𝑘+𝑓^′(0)𝜏𝜔_𝑘)²=𝜔²_𝑘Δℎ^′(𝑉_𝑘), (2.23)whereΔ=((3+𝜏𝛼₁)𝜔⁴_𝑘−𝛼₂𝜔²_𝑘−𝑓^′(0))²+(𝜏𝜔⁵_𝑘−(3𝛼₁+𝜏𝛼₂)𝜔³_𝑘+𝑓^′(0)𝜏𝜔_𝑘)²; (2.24)and for ℎ^′(𝑉_𝑘)≠0(𝑘=1,2,3,4), ℎ(0)=𝑓^′2(0)>0 and lim_𝑉→±∞ℎ(𝑉)=∞, we can know that Re(𝑑𝜆(𝜏^𝑗_𝑘)/𝑑𝜏)>0 when 𝑘=2,4; and Re(𝑑𝜆(𝜏^𝑗_𝑘)/𝑑𝜏)<0 when 𝑘=1,3. This completes the proof.

From ℎ(0)=𝑓^′2(0)>0 and lim_𝑉→±∞ℎ(𝑉)=∞, it is easy to know that: if ℎ(𝑉) satisfies ℎ^′(𝑉_𝑖)≠0(𝑖=1,2,3,4), if the equation ℎ(𝑉)=0 has positive roots, then the number of the roots must be even; and from Lemma 2.3, we have that the sign of 𝛼^′(𝜏^𝑗_𝑘) changes as 𝜏^𝑗_𝑘 varies, and then the stability switches may happen.

From Lemmas 2.1–2.3 and the theory in [9], we have the following.

Lemma 2.4. Suppose that (𝐇_𝟏), (𝐇_𝟐) and ℎ^′(𝑉_𝑖)≠0(𝑖=1,2,3,4) are satisfied.
(i)If conditions (i) and (ii) in Lemma 2.2 are not satisfied, then all the roots of (2.3) have negative real parts for all 𝜏≥0.(ii)If one of conditions (i) and (ii) in Lemma 2.2 is satisfied, let 𝜏^∗=min(𝜏⁰₂,𝜏⁰₄), then all roots of (2.3) have negative real parts when 𝜏∈[0,𝜏^∗); and there may exist an integer 𝑚≥0 such that 0<𝜏₁<𝜏₂<⋯<𝜏_𝑚−1<𝜏_𝑚<𝜏_𝑚+1<⋯, and all the roots of (2.3) have negative real parts when 𝜏∈[0,𝜏₁)∪(𝜏₂,𝜏₃)∪⋯∪(𝜏_𝑚−1,𝜏_𝑚), and (2.3) has at least a pair of roots with positive real parts when 𝜏∈(𝜏₁,𝜏₂)∪(𝜏₃,𝜏₄)∪⋯∪(𝜏_𝑚,∞), where 𝜏_𝑚∈{𝜏^𝑗_𝑘}.

From Lemma 2.4 and applying the Hopf bifurcation theorem for functional differential equations [12, Chapter 11, Theorem 1.1], we have the following results.

Theorem 2.5. Suppose (𝐇_𝟏), (𝐇_𝟐), and ℎ^′(𝑉_𝑖)≠0(𝑖=1,2,3,4) are satisfied.
(i)If conditions (i) and (ii) in Lemma 2.2 are not satisfied, then the trivial solution (0,0,0,0) of system (2.1) is asymptotically stable when 𝜏>0.(ii)If one of conditions (i) and (ii) in Lemma 2.2 is satisfied, let 𝜏^∗=min{𝜏⁰₂,𝜏⁰₄}, then the trivial solution (0,0,0,0) of system (2.1) is asymptotically stable when 𝜏∈[0,𝜏^∗); and there may exist an integer 𝑚≥0 such that 0<𝜏₁<𝜏₂<⋯<𝜏_𝑚−1<𝜏_𝑚<𝜏_𝑚+1<⋯, and the trivial solution (0,0,0,0) of system (2.1) is asymptotically stable when 𝜏∈[0,𝜏₁)∪(𝜏₂,𝜏₃)∪⋯∪(𝜏_𝑚−1,𝜏_𝑚), and is unstable when 𝜏∈(𝜏₁,𝜏₂)∪(𝜏₃,𝜏₄)∪⋯∪(𝜏_𝑚,∞), where 𝜏_𝑚∈{𝜏^𝑗_𝑘}.(iii)The system (2.1) undergoes a Hopf bifurcation at the origin when 𝜏=𝜏^𝑗_𝑘, with 𝑘=1,2,3,4;𝑗=0,1,2,….

3. Direction and Stability of the Hopf Bifurcation

In this section, we will study the direction, stability, and the period of the bifurcating periodic solution. The method we used is based on the normal form method and the center manifold theory presented by Hassard et al. [10].

We first rescale the time by 𝑡→𝑡/𝜏 to normalize the delay so that system (2.1) can be written as the forṁ̇̇̇𝑥=𝜏𝑦,𝑦=𝜏𝑢,𝑢=𝜏𝑣,𝑣=−𝛼₂𝜏𝑢−𝛼₁𝜏𝑣−𝜏𝑓(𝑥)−𝜏𝜙(𝑦(𝑡−1)). (3.1)The linearization around (0,0,0,0) is given bẏ̇̇̇𝑥=𝜏𝑦,𝑦=𝜏𝑢,𝑢=𝜏𝑣,𝑣=−𝜏𝑓^′(0)𝑥−𝛼₂𝜏𝑢−𝛼₁𝜏𝑣−𝜏𝜙^′(0)𝑦(𝑡−1); (3.2)and the nonlinear term is000𝐹=−𝑓^′′(0)2𝑥²−𝑓^′′′(0)𝑥³6−𝜏𝜙^′′(0)𝑦²(𝑡−1)2−𝜏𝜙^′′′(0)𝑦³(𝑡−1)6−⋯. (3.3)The characteristic equation associated with (3.2) is𝛾⁴+𝛼₁𝜏𝛾³+𝛼₂𝜏²𝛾²+𝜙^′(0)𝜏³𝛾𝑒^−𝛾+𝜏⁴𝑓^′(0)=0. (3.4)Comparing (3.4) with (2.3), one can find out that 𝛾=𝜏𝜆, and hence, (3.4) has a pair of imaginary roots ±𝑖𝜏^𝑗_𝑘𝜔_𝑘, when 𝜏=𝜏^𝑗_𝑘 for 𝑘=1,2,3,4, 𝑗=0,1,2,…, and the transversal condition holds.

Let 𝜏=𝜏₀+𝜇, 𝜇∈𝑅 where 𝜏₀∈{𝜏^𝑗_𝑘}, 𝜔₀∈{𝜔_𝑘}, 𝑘=1,2,3,4, 𝑗=0,1,2,…. Then 𝜇=0 is the Hopf bifurcation value for (3.1). Let 𝑖𝜏₀𝜔₀ be the root of (3.4).

For 𝜑∈𝐶([−1,0],𝑅⁴), let𝐿_𝜇𝜑=(𝜏₀+𝜇)010000100001−𝑓^′(0)0−𝛼₂−𝛼₁𝜑(0)−(𝜏₀+𝜇)𝜙^′=000(0)0000000000000100𝜑(−1),𝐹(𝜇,𝜑)−(𝜏₀+𝜇)𝑓^′′(0)𝜑²₁(0)2−(𝜏₀+𝜇)𝑓^′′′(0)𝜑³₁(0)6−(𝜏₀+𝜇)𝜙^′′(0)𝜑²₂(−1)2−(𝜏₀+𝜇)𝜙^′′′(0)𝜑³₂(−1)6.−⋯ (3.5)By the Riesz representation theorem, there exists a matrix whose components are bounded variation functions 𝜂(𝜃,𝜇) in 𝜃∈[−1,0] such that𝐿_𝜇𝜑=⁰₋₁𝑑𝜂(𝜃,𝜇)𝜑(𝜃),𝜑∈𝐶([−1,0],𝑅⁴). (3.6)In fact, we choose𝜂(𝜃,𝜇)=(𝜏₀+𝜇)010000100001−𝑓^′(0)0−𝛼₂−𝛼₁𝛿(𝜃)+(𝜏₀+𝜇)𝜙^′(0)0000000000000100𝛿(𝜃+1), (3.7)where 𝛿(𝜃)=1,𝜃=0,0,𝜃≠0. (3.8)

For 𝜑∈𝐶¹([−1,0],𝐶⁴), define𝐴(𝜇)𝜑=𝑑𝜑(𝜃)𝑑𝜃,𝜃∈[−1,0),₀₋₁𝑑𝜂(𝑡,𝜇)𝜑(𝑡),𝜃=0,𝑅(𝜇)𝜑=0,𝜃∈[−1,0),𝐹(𝜇,𝜑),𝜃=0. (3.9)Hence, we can rewrite (3.1) in the following form:𝑤̇_𝑡=𝐴(𝜇)𝑤_𝑡+𝑅(𝜇)𝑤_𝑡, (3.10)where 𝑤=(𝑥,𝑦,𝑢,𝑣)^𝑇, 𝑤_𝑡=𝑤(𝑡+𝜃) for 𝜃∈[−1,0].

For 𝜓∈𝐶¹([0,1],𝐶⁴), define𝐴^∗−𝜓(𝑠)=𝑑𝜓(𝑠)𝑑𝑠,𝑠∈(0,1],₀₋₁𝜓(−𝑡)𝑑𝜂(𝑡,0),𝑠=0. (3.11)For 𝜑∈𝐶([−1,0],𝐶⁴) and 𝜓∈𝐶([0,1],𝐶⁴), define the bilinear form⟨𝜓,𝜑⟩=𝜓(0)𝜑(0)−⁰₋₁^𝜃_𝜉=0𝜓^𝑇(𝜉−𝜃)𝑑𝜂(𝜃)𝜑(𝜉)𝑑𝜉, (3.12)where 𝜂(𝜃)=𝜂(𝜃,0). Then 𝐴^∗ and 𝐴(0) are adjoint operators, and ±𝑖𝜏₀𝜔₀ are eigenvalues of 𝐴(0). Thus, they are also eigenvalues of 𝐴^∗.

By direct computation, we obtain that1𝑞(𝜃)=𝑖𝜔₀−𝜔²₀−𝑖𝜔³₀𝑒^{𝑖𝜏0𝜔0𝜃} (3.13)is the eigenvector of 𝐴(0) corresponding to 𝑖𝜏₀𝜔₀, and𝑞^∗(𝜃)=𝐷−𝑖𝜔³₀−𝜔²₀𝛼₁+𝑖𝜔₀𝛼₂+𝜙^′(0)𝑒^{−𝑖𝜏0𝜔0}−𝜔²₀+𝑖𝜔₀𝛼₁+𝛼₂𝑖𝜔₀+𝛼₁1^𝑇𝑒^{𝑖𝜏0𝜔0𝑠} (3.14)is the eigenvector of 𝐴^∗ corresponding to −𝑖𝜏⁰𝜔₀. Moreover,⟨𝑞^∗,𝑞⟩=1,⟨𝑞^∗,𝑞⟩=0, (3.15)where1𝐷=−𝜔²₀𝛼₁+𝜙^′(0)𝑒^{−𝑖𝜏0𝜔0}−𝑖𝜏₀𝜔₀𝜙^′(0)𝑒^{𝑖𝜏0𝜔0}. (3.16)

Using the same notation as in Hassard et al. [10], we first compute the coordinates to describe the center manifold 𝒞₀ at 𝜇=0. Let 𝑤_𝑡 be the solution of (3.1) when 𝜇=0.

Define𝑧_𝑡=⟨𝑞^∗,𝑤_𝑡⟩,𝑊(𝑡,𝜃)=𝑤_𝑡(𝜃)−2Re{𝑧(𝑡)𝑞(𝜃)}. (3.17)On the center manifold 𝒞₀, we have𝑊(𝑡,𝜃)=𝑊(𝑧(𝑡),𝑧(𝑡),𝜃), (3.18)where𝑊(𝑧,𝑧,𝜃)=𝑊₂₀𝑧(𝜃)²2+𝑊₁₁(𝜃)𝑧𝑧+𝑊₀₂(𝜃)𝑧²2+𝑊₃₀𝑧(𝜃)³6+⋯, (3.19)𝑧 and 𝑧 are local coordinates for center manifold 𝒞₀ in the direction of 𝑞^∗ and 𝑞^∗. Note that 𝑊 is real if 𝑤_𝑡 is real. We consider only real solutions.

For solution 𝑤_𝑡 in 𝒞₀ of (3.1), since 𝜇=0,𝑧(̇𝑡)=𝑖𝜏₀𝜔₀𝑧+⟨𝑞^∗(𝜃),𝐹(𝑊(𝑧,𝑧,𝜃)+2Re{𝑧(𝑡)𝑞(𝜃)})⟩=𝑖𝜏₀𝜔₀𝑧+𝑞^∗(0)𝐹(𝑊(𝑧,𝑧,0)+2Re{𝑧(𝑡)𝑞(0)})_def=𝑖𝜏₀𝜔₀𝑧+𝑞^∗(0)𝐹₀(𝑧,𝑧). (3.20)We rewrite this as𝑧(̇𝑡)=𝑖𝜏₀𝜔₀𝑧(𝑡)+𝑔(𝑧,𝑧), (3.21)where𝐹₀(𝑧,𝑧)=𝐹_𝑧2𝑧²2+𝐹_𝑧²𝑧²2+𝐹_𝑧𝑧𝑧𝑧+𝐹_𝑧2𝑧𝑧²𝑧2+⋯, (3.22) 𝑔(𝑧,𝑧)=𝑞^∗(0)𝐹(𝑊(𝑧,𝑧,0)+2Re{𝑧(𝑡)𝑞(0)})=𝑔₂₀𝑧²2+𝑔₁₁𝑧𝑧+𝑔₀₂𝑧²2+𝑔₂₁𝑧²𝑧2+⋯. (3.23)Compare the coefficients of (3.20) and (3.21), noticing (3.23), we have𝑔₂₀=𝑞^∗(0)𝐹_𝑧2,𝑔₁₁=𝑞^∗(0)𝐹_𝑧𝑧,𝑔₀₂=𝑞^∗(0)𝐹_𝑧²,𝑔₂₁=𝑞^∗(0)𝐹_𝑧2𝑧. (3.24)By (3.10) and (3.21), it follows that𝑊̇=𝑤̇_𝑡−𝑧𝑞̇−𝑧̇𝑞=𝐴𝑊−2Re{𝑞^∗(0)𝐹₀𝑞(𝜃)}−1≤𝜃<0𝐴𝑊−2Re{𝑞^∗(0)𝐹₀𝑞(0)}+𝐹₀𝜃=0_def=𝐴𝑊+𝐻(𝑧,𝑧,𝜃), (3.25)where𝐻(𝑧,𝑧,𝜃)=𝐻₂₀𝑧²2+𝐻₁₁(𝜃)𝑧𝑧+𝐻₀₂(𝜃)𝑧²2+⋯. (3.26)Expanding the above series and comparing the coefficients, we obtain(𝐴−2𝑖𝜏₀𝜔₀𝐼)𝑊₂₀(𝜃)=−𝐻₂₀(𝜃),𝐴𝑊₁₁(𝜃)=−𝐻₁₁(𝜃),(𝐴−2𝑖𝜏₀𝜔₀)𝑊₀₂(𝜃)=−𝐻₀₂(𝜃). (3.27)Notice that𝑤_𝑡(𝜃)=𝑊(𝑧,𝑧,𝜃)+𝑧𝑞(𝜃)+𝑧𝑞(𝜃), (3.28)that is,𝑥(𝑡)=𝑧+𝑧+𝑊⁽¹⁾(𝑧,𝑧,0)=𝑧+𝑧+𝑊⁽¹⁾₂₀𝑧(0)²2+𝑊⁽¹⁾₁₁(0)𝑧𝑧+𝑊⁽¹⁾₀₂(0)𝑧²2+⋯,𝑦(𝑡−1)=𝑖𝜔₀𝑒^{−𝑖𝜏0𝜔0}𝑧−𝑖𝜔₀𝑒^{𝑖𝜏0𝜔0}𝑧+𝑊⁽²⁾(𝑧,𝑧,−1)=𝑖𝜔₀𝑒^{−𝑖𝜏0𝜔0}𝑧−𝑖𝜔₀𝑒^{𝑖𝜏0𝜔0}𝑧+𝑊⁽²⁾₂₀𝑧(−1)²2+𝑊⁽²⁾₁₁(−1)𝑧𝑧+𝑊⁽²⁾₀₂(−1)𝑧²2+⋯. (3.29)Thus𝑥²𝑧(𝑥)=2²2+2𝑧𝑧+2𝑧²2+(4𝑊¹₁₁(0)+2𝑊¹₂₀𝑧(0))²𝑧2𝑥+⋯,₃𝑧(𝑥)=6²𝑧2𝑦+⋯,₂(𝑡−1)=−2𝜔²₀𝑒^{−2𝑖𝜏0𝜔0}𝑧²2+2𝜔²₀𝑧𝑧−2𝜔²₀𝑒^{2𝑖𝜏0𝜔0}𝑧²2+(4𝑖𝜔₀𝑊⁽²⁾₁₁(−1)𝑒^{−𝑖𝜏0𝜔0}−𝑖𝜔₀2𝑊⁽²⁾₂₀(−1)𝑒^{𝑖𝜏0𝜔0})𝑧²𝑧2𝑦+⋯,₃(𝑡−1)=6𝑖𝜔³₀𝑒^{−𝑖𝜏0𝜔0}𝑧²𝑧2+⋯; (3.30)and we have𝐹(0,𝜔_𝑡000−𝜏)=₀𝑓^′′(0)𝑥²(𝑡)2−𝜏₀𝑓^′′′(0)𝑥³(𝑡)6−𝜏₀𝜙^′′(0)𝑦²(𝑡−1)2−𝜏₀𝜙^′′(0)𝑦³(𝑡−1)6,𝑞−⋯_∗(0)=𝐷−𝑖𝜔³₀−𝜔²₀𝛼₁+𝑖𝜔₀𝛼₂+𝜙^′(0)𝑒^{−𝑖𝜏0𝜔0}−𝜔²₀+𝑖𝜔₀𝛼₁+𝛼₂𝑖𝜔₀+𝛼₁1^𝑇,𝑔(𝑧,𝑧)=𝑞^∗(0)𝐹₀(𝑧,𝑧). (3.31)Then we have𝑔₂₀=−𝜏₀𝐷𝑓^′′(0)+𝜏₀𝐷𝜙^′′(0)𝜔²𝑒^{−2𝑖𝜏0𝜔0},𝑔₁₁=−𝜏₀𝐷𝑓^′′(0)−𝜏₀𝐷𝜙^′′(0)𝜔²,𝑔₀₂=−𝜏₀𝐷𝑓^′′(0)+𝜏₀𝐷𝜙^′′(0)𝜔²𝑒^{2𝑖𝜏0𝜔0},𝑔₂₁=−𝜏₀𝐷𝑓^′′(0)[2𝑊⁽¹⁾₁₁(0)+𝑊⁽¹⁾₂₀(0)]−𝜏₀𝐷𝑓^′′′(0)−𝜏₀𝐷𝜙^′′(0)[2𝑖𝜔₀𝑊⁽²⁾₁₁(−1)𝑒^{−𝑖𝜏0𝜔0}−𝑖𝜔₀𝑊⁽²⁾₂₀(−1)𝑒^{𝑖𝜏0𝜔0}]−𝜏₀𝐷𝜙^′′′(0). (3.32)So we only need to find out 𝑊⁽¹⁾₁₁(0), 𝑊⁽¹⁾₂₀(0), 𝑊⁽²⁾₁₁(−1), and 𝑊⁽²⁾₂₀(−1) to obtain 𝑔₂₁.

When 𝜃∈[−1,0), we have𝐻(𝑧,𝑧,𝜃)=−2Re{𝑞^∗(0)𝐹₀𝑞(𝜃)}=−𝑞^∗(0)𝐹₀𝑞(𝜃)−𝑞^∗(0)𝐹₀𝑞(𝜃)=−𝑔𝑞(𝜃)−𝑔𝑞(𝜃). (3.33)Comparing the coefficients with (3.26), we get𝐻₂₀(𝜃)=−𝑔₂₀𝑞(𝜃)−𝑔₀₂𝐻𝑞(𝜃),₁₁(𝜃)=−𝑔₁₁𝑞(𝜃)−𝑔₁₁𝑞(𝜃). (3.34)From (3.27), (3.32), (3.33), and (3.34), we derive𝑊̇₂₀(𝜃)=2𝑖𝜏₀𝜔₀𝑊₂₀(𝜃)+𝑔₂₀𝑞(𝜃)+𝑔₀₂𝑊̇𝑞(𝜃),₁₁(𝜃)=𝑔₁₁𝑞(𝜃)+𝑔₁₁𝑞(𝜃). (3.35)Then we can get𝑊₂₀(𝜃)=𝑖𝑔₂₀𝜏₀𝜔₀𝑞(0)𝑒^{𝑖𝜏0𝜔0𝜃}+𝑖𝑔₀₂3𝜏₀𝜔₀𝑞(0)𝑒^{−𝑖𝜏0𝜔0𝜃}+𝐸₁𝑒^{2𝑖𝜏0𝜔0𝜃},𝑊₁₁𝑔(𝜃)=¹¹𝑖𝜏₀𝜔₀𝑞(0)𝑒^{𝑖𝜏0𝜔0𝜃}−𝑔₁₁𝑖𝜏₀𝜔₀𝑞(0)𝑒^{−𝑖𝜏0𝜔0𝜃}+𝐸₂. (3.36)Notice that(2𝑖𝜏₀𝜔₀𝐼−∫⁰₋₁𝑒^{2𝑖𝜏0𝜔0𝜃}𝑑𝜂(𝜃))𝐸₁=𝐹_𝑧2,(∫⁰₋₁𝑑𝜂(𝜃))𝐸₂=−𝐹_𝑧𝑧. (3.37)We obtain𝐸₁=𝐸^∗2𝑖𝜔₀𝐸^∗−4𝜔²₀𝐸^∗−8𝑖𝜔³₀𝐸^∗,𝐸₂=−𝑓^′′(0)+𝜙^′′(0)𝑓^′000(0), (3.38)where 𝐸^∗=−𝑓^′′(0)+𝜙^′′(0)𝜔²₀𝑒^{−2𝑖𝜏0𝜔0}𝑓^′(0)+2𝑖𝜔₀𝑒^{−2𝑖𝜏0𝜔0}−4𝜔²₀𝛼₂−8𝑖𝜔³₀(2𝑖𝜔₀+𝛼₁). (1)Hence𝑊⁽¹⁾₂₀(0)=𝑖𝑔₂₀𝜏₀𝜔₀𝑞₁𝑖(0)+𝑔₀₂3𝜏₀𝜔₀𝑞₁(0)+𝐸^∗𝑒^{−2𝑖𝜏0𝜔0}=𝑖(−𝐷𝑓^′′(0)+𝐷𝜙^′′(0)𝜔²₀𝑒^{−2𝑖𝜏0𝜔0})𝜔₀+𝑖(−𝐷𝑓^′′(0)+𝐷𝜙^′′(0)𝜔²₀𝑒^{−2𝑖𝜏0𝜔0})3𝜔₀+𝐸^∗𝑒^{−2𝑖𝜏0𝜔0},𝑊⁽²⁾₂₀(−1)=𝑖𝑔₂₀𝜏₀𝜔₀𝑞₂(0)𝑒^{−𝑖𝜏0𝜔0}+𝑖𝑔₀₂3𝜏₀𝜔₀𝑞₂(0)𝑒^{𝑖𝜏0𝜔0}+(2𝑖𝜔₀)𝐸^∗𝑒^{−2𝑖𝜏0𝜔0}=−(−𝐷𝑓^′′(0)+𝐷𝜙^′′(0)𝜔²₀𝑒^{−2𝑖𝜏0𝜔0})𝑒^{−𝑖𝜏0𝜔0}+−𝐷𝑓^′′(0)+𝐷𝜙^′′(0)𝜔²₀𝑒^{−2𝑖𝜏0𝜔0}3𝑒^{𝑖𝜏0𝜔0}+(2𝑖𝜔₀)𝐸^∗𝑒^{−2𝑖𝜏0𝜔0},𝑊⁽¹⁾₁₁𝑔(0)=¹¹𝑖𝜏₀𝜔₀𝑞₁(0)−𝑔₁₁𝑖𝜏₀𝜔₀𝑞₁𝑓(0)−^′′(0)+𝜔²₀𝜙^′(0)𝑓^′=−(0)𝐷𝑓^′′(0)−𝐷𝜙^′′(0)𝜔²₀𝑖𝜔₀+𝐷𝑓^′′(0)+𝐷𝜙^′′(0)𝜔²₀𝑖𝜔₀−𝑓^′′(0)+𝜔²₀𝜙^′(0)𝑓^′,𝑊(0)₍₂₎₁₁𝑔(−1)=¹¹𝑖𝜏₀𝜔₀𝑞₂(0)𝑒^{−𝑖𝜏0𝜔0}−𝑔₁₁𝑖𝜏₀𝜔₀𝑞₂(0)𝑒^{𝑖𝜏0𝜔0}=(−𝐷𝑓^′′(0)−𝐷𝜙^′′(0)𝜔²₀)𝑒^{−𝑖𝜏0𝜔0}+(𝐷𝑓^′′(0)+𝐷𝜙^′′(0)𝜔²₀)𝑒^{𝑖𝜏0𝜔0}. (3.39)Consequently, from (3.32),𝑔₂₁=−𝜏₀𝐷𝑓^′′(0)[2𝑊⁽¹⁾₁₁(0)+𝑊⁽¹⁾₂₀(0)]−𝜏₀𝐷𝑓^′′′(0)−𝜏₀𝐷𝜙^′′(0)[2𝑖𝜔₀𝑊⁽²⁾₁₁(−1)𝑒^{−𝑖𝜏0𝜔0}−𝑖𝜔₀𝑊⁽²⁾₂₀(−1)𝑒^{𝑖𝜏0𝜔0}]−𝜏₀𝐷𝜙^′′′(0)=−2𝜏₀𝐷𝑓^′′−(0)[𝐷𝑓^′′(0)−𝐷𝜙^′′(0)𝜔²₀𝑖𝜔₀+𝐷𝑓^′′(0)+𝐷𝜙^′′(0)𝜔²₀𝑖𝜔₀−𝑓^′′(0)+𝜔²₀𝜙^′′(0)𝑓^′(0)]−𝜏₀𝐷𝑓^′′(0)×[𝑖(−𝐷𝑓^′′(0)+𝐷𝜙^′′(0)𝜔²₀𝑒^{−2𝑖𝜏0𝜔0})𝜔₀+𝑖(−𝐷𝑓^′′(0)+𝐷𝜙^′′(0)𝜔²₀𝑒^{−2𝑖𝜏0𝜔0})3𝜔₀+𝐸^∗𝑒^{−2𝑖𝜏0𝜔0}]−𝜏₀𝐷𝑓^′′′(0)−2𝜏₀𝐷𝜙^′′(0)(𝑖𝜔₀)[(−𝐷𝑓^′′(0)−𝐷𝜙^′′(0)𝜔²₀)𝑒^{−𝑖𝜏0𝜔0}+(𝐷𝑓^′′(0)+𝐷𝜙^′′(0)𝜔²₀)𝑒^{𝑖𝜏0𝜔0}]𝑒^{−𝑖𝜏0𝜔0}+𝜏₀𝐷𝜙^′′(0)(𝑖𝜔₀)[−(−𝐷𝑓^′′(0)+𝐷𝜙^′′(0)𝜔²₀𝑒^{−2𝑖𝜏0𝜔0})𝑒^{−𝑖𝜏0𝜔0}]𝑒^{𝑖𝜏0𝜔0}+𝜏₀𝐷𝜙^′′(0)(𝑖𝜔₀)[−𝐷𝑓^′′(0)+𝐷𝜙^′′(0)𝜔²₀𝑒^{−2𝑖𝜏0𝜔0}3𝑒^{𝑖𝜏0𝜔0}+(2𝑖𝜔₀)𝐸^∗𝑒^{−2𝑖𝜏0𝜔0}]𝑒^{𝑖𝜏0𝜔0}−𝜏₀𝐷𝜙^′′′(0). (3.40)Substituting 𝑔₂₀, 𝑔₁₁, 𝑔₀₂, and 𝑔₂₁ into𝐶₁𝑖(0)=2𝜏₀𝜔₀(𝑔₂₀𝑔₁₁−2|𝑔₁₁|²−13|𝑔₀₂|²𝑔)+²¹2, (3.41)we can obtain Re𝐶₁(0). Then we obtain the sign of𝛽₂=2Re𝐶₁𝜇(0),₂=−Re𝐶₁(0)𝛼^′(𝜏₀). (3.42)

By the general theory due to Hassard et al. [10], we know that the quantity of 𝛽₂ determines the stability of the bifurcating periodic solutions on the center manifold, and 𝜇₂ determines the direction of the bifurcation; and we have the following.

Theorem 3.1. (i) If 𝜇₂>0(<0), then the Hopf bifurcation at the origin of system (1.1) is supercritical (subcritical).
(ii) If 𝛽₂<0(>0), then the bifurcating periodic solutions of system (1.1) are asymptotically stable (unstable).

4. An Example and Numerical Simulations

In this section, we give an example and present some numerical simulations to illustrate the analytic results.

Example 4.1. Consider the following equation:𝑥⁽⁴⁾(𝑡)+𝑥⃛(𝑡)+14𝑥̈(𝑡)+12sin𝑥̇(𝑡−𝜏)+𝑥³(𝑡)+18𝑥(𝑡)=0. (4.1)Clearly,𝛼₁=1,𝛼₂=14,𝜙(𝑥̇(𝑡−𝜏))=12sin𝑥̇(𝑡−𝜏),𝑓(𝑥)=𝑥³+18𝑥,𝜙(0)=𝜙^′′(0)=0,𝜙^′(0)=12,𝜙^′′′(0)=−12,𝑓(0)=𝑓^′′(0)=0,𝑓^′(0)=18,𝑓^′′′(0)=6. (4.2)

By direct computation, we know (𝐇_𝟏) and (𝐇_𝟐) are satisfied. That is, the data satisfy the conditions of Lemma 2.1. The characteristic equation is𝜆⁴+𝜆³+14𝜆²+12𝜆𝑒^−𝜆𝜏+18=0; (4.3)and we can obtainℎ(𝑉)=𝑉⁴−27𝑉³+232𝑉²−648𝑉+324. (4.4) As shown in Figure 1, the equation ℎ(𝑉)=0 has four roots as𝑉₁=0.633,𝑉₂=4.166,𝑉₃=10.464,𝑉₄=11.737; (4.5)and ℎ^′(𝑉₁)<0,ℎ^′(𝑉₂)>0,ℎ^′(𝑉₃)<0,ℎ^′(𝑉₄)>0. (4.6)Hence𝜔₁=0.796,𝜔₂=2.041,𝜔₃=3.235,𝜔₄𝜏=3.426,₀₁=5.988,𝜏⁰₂=0.596,𝜏⁰₃=0.158,𝜏⁰₄=0.061. (4.7)For 𝜏⁰₄<𝜏⁰₃<𝜏⁰₂<𝜏⁰₁, we obtain that the zero solution of system (4.1) is asymptotically stable when 𝜏∈[0,0.061)∪(0.158,0.596).

Figure 1 

The curve of function $h(V)=V^4-27V^3+232V^2-648V+324$.

According to the formula given in Section 3, we can obtain that𝐷₄=−1.898+0.133𝑖,𝐷₃=−0.323+0.056𝑖,𝐷₂𝑔=−0.061+0.027𝑖,₍₄₎₂₁=−0.697+0.049𝑖,𝑔⁽³⁾₂₁=−0.306+0.053𝑖,𝑔⁽²⁾₂₁𝑔=−0.217+0.098𝑖,₀₂=𝑔₁₁=𝑔₀₂=𝐸^∗=0. (4.8)Then we have𝐶₁𝑔(0)=²¹2. (4.9)Hence, when 𝜏∈{𝜏⁰₁,𝜏⁰₂,𝜏⁰₃}, we have𝛽₂=2Re𝐶₁𝜇(0)<0,₂=−Re𝐶₁(0)𝛼^′(𝜏₀)>0. (4.10)

Conclusion of (4.1)
The zero solution of system (4.1) is asymptotically stable when 𝜏∈[0,0.061)∪(0.158,0.596). The Hopf bifurcation at the origin when 𝜏₀=𝜏⁰_𝑘 is supercritical, and the bifurcating periodic solutions are asymptotically stable.
The following is the results of numerical simulations to system (4.1).
(i)We choose 𝜏=0.4∈(0.158,0.596), then the zero solution of system (4.1) is asymptotically stable, as shown in Figure 2.(ii) We choose 𝜏=0.64 being near to 𝜏⁰₂=0.596, a periodic solution bifurcates from the origin and is asymptotically stable, as shown in Figure 3.

(a) y-x phase plots

(b) y waveform

(c) x waveform

(d) v-u phase plots

(e) v waveform

(f) u waveform

Figure 2 

The zero solution of system (4.1) is asymptotically stable when $\tau =0.4$.

(a) y-x phase plots

(b) y waveform

(c) x waveform

(d) v-u phase plots

(e) v waveform

(f) u waveform

Figure 3 

For (4.1) with $\tau =0.64>{\tau }_2^0$ and sufficiently near ${\tau }_2^0=0.596$, the bifurcating periodic solution from zero solution occurs and is asymptotically stable.

5. Conclusion

In this paper, we consider a certain fourth-order delay differential equation. The linear stability is investigated by analyzing the associated characteristic equation. It is found that there may exist the stability switches when delay varies, and the Hopf bifurcation occurs when the delay passes through a sequence of critical values. Then the direction and the stability of the Hopf bifurcation are determined using the normal form method and the center manifold theorem. Finally, an example is given and numerical simulations are carried out to illustrate the results. By using Lyapunov's second method, Sadek [1] investigated the stability of system (1.1). The main result is as the following.

Theorem 5.1. Suppose that the following hold.
(i)There are constants 𝛼₁>0,𝛼₂>0,𝛼₃>0,𝛼₄>0, and Δ>0 such that(𝛼₁𝛼₂−𝜙^′(𝑦))𝛼₃−𝛼²₁𝛼₄≥Δ (5.1)for all 𝑦.(ii)𝑓(0)=0,𝑥𝑓(𝑥)>0(𝑥≠0),𝐹(𝑥)=∫^𝑥₀𝑓(𝜉)𝑑𝜉→∞ as |𝑥|→∞, and0≤𝛼₄−𝑓^′(𝑥)≤𝜀𝑑₀𝛼⁰₁ (5.2)for all 𝑥, where 𝜀 is a positive constant such that1𝜀=min(𝛼₁,Δ4𝛼₁𝛼₃𝑑₀,𝛼₃4𝛼₄𝑑₀(2Δ𝛼₄𝛼₁𝛼²₃−𝛿₁)) (5.3)with 𝑑₀=𝛼₁𝛼₂+𝛼₂𝛼₃𝛼₄⁻¹.(iii)𝜙(0)=0 and 𝜙^′(𝑦)≥𝛼₃>0 for all 𝑦, and 0≤𝜙^′(𝑦)−𝜙(𝑦)/𝑦≤𝛿₁<2Δ𝛼₄/𝛼₁𝛼²₃ for all 𝑦≠0. Then the zero solution of (1.1) is asymptotically stable, provided that𝜀𝜏<min{𝑑₂𝛼₂,Δ2𝛼₁𝛼₃(𝛼₂𝛼₃,𝛼+2𝜇)¹𝜀𝑑₁𝛼₂𝛼₃}, (5.4)with 𝜇=(𝛼₂𝛼₃/2)(𝑑₁+𝑑₂+1)>0.

Comparing Theorem 5.1 with Theorem 2.5 obtained in Section 2, one can find out that if the sufficient conditions to ensure the globally asymptotical stability of system (1.1) given in [10] are not satisfied, we can also get the stability of system (1.1), but here the stability means local stability, and the system undergoes a Hopf bifurcation at the origin. Otherwise, here we just need to give the condition on the origin of 𝑓(𝑥) and 𝜙(𝑥), the condition is relatively weak.

Acknowledgment

This work was supported by the National Natural Science Foundation of China (No10771045), and Program of Excellent Team in Harbin Institute of Technology.