# Hopf Bifurcation Analysis for the van der Pol Equation with Discrete and Distributed Delays

Authors: Xiaobing Zhou; Murong Jiang; Xiaomei Cai
Journal: Discrete Dynamics in Nature and Society
 (2011)

---

## Abstract

We consider the van der Pol equation with discrete and distributed delays. Linear
stability of this equation is investigated by analyzing the transcendental characteristic equation of
its linearized equation. It is found that this equation undergoes a sequence of Hopf bifurcations
by choosing the discrete time delay as a bifurcation parameter. In addition, the properties of Hopf
bifurcation were analyzed in detail by applying the center manifold theorem and the normal form theory. Finally, some numerical simulations are performed to illustrate and verify the theoretical analysis.


---

## Body

## 1. Introduction

Since its introduction in 1927, the van der Pol equation [1] has served as a basic model of self-excited oscillations in physics, electronics, biology, neurology and other disciplines [2–15]. The intensively studied van der Pol equation is governed by the following second-order nonlinear damped oscillatory system:(1.1)ẋ(t)=y(t)-f(x(t)),ẏ(t)=-x(t),
where f(x)=ax+bx3, a and b are real constants.In 1999, Murakami [16] introduced a discrete time delay into system (1.1), and obtained the following pair of delay differential equations(1.2)ẋ(t)=y(t-τ)-f(x(t-τ)),ẏ(t)=-x(t-τ).
By using the center manifold theorem, he [16] found that periodic solutions existed in system (1.2). The stability of bifurcating periodic solutions was discussed in detail by Yu and Cao [17].It is well known that dynamical systems with distributed delay are more general than those with discrete delay. So Liao et al. [18] proposed the following van der Pol equation with distributed delay:(1.3)ẋ(t)=∫0∞F(t-s)y(s)ds-f[∫0∞F(t-s)x(s)ds],ẏ(t)=-∫0∞F(t-s)x(s)ds.
The existence of Hopf bifurcation and the stability of the bifurcating periodic solutions of system (1.3) were analyzed in [18–25] for the weak and strong kernels, respectively.In [26], Liao considered the following system with two discrete time delays:(1.4)ẋ(t)=y(t-τ2)-f(x(t-τ1)),ẏ(t)=-x(t-τ1).
By choosing one of the delays as a bifurcation parameter, system (1.4) was found to undergo a sequence of Hopf bifurcations. The author had also found that resonant codimension two bifurcation occurred in this system.In this paper, we consider the following van der Pol equation with discrete and distributed time delays:(1.5)ẋ(t)=∫-∞tF(t-s)y(s)ds-f(x(t-τ)),ẏ(t)=-x(t-τ),
with initial conditions x(θ1)=φ1(θ1), y(θ2)=φ2(θ2), (-τ≤θ1≤0,-∞<θ2≤0), τ≥0, where φ1(θ1) and φ2(θ2) are bounded and are continuous functions. The weight function F(s) is a nonnegative bounded function, which describes the influence of the past states on the current dynamics. It is assumed that the presence of the distributed time delay does not affect the system equilibrium. Hence, O(0,0) is the unique equilibrium of system (1.5).We normalize the kernel in the following way:(1.6)∫0∞F(s)ds=1.
Usually, the following form:(1.7)F(s)=αp+1spe-αsp!,p=0,1,2,…,
is taken as the kernel. The kernel is called “weak” when p=0 and “strong” when p=1, respectively. The analysis of weak and strong kernels is similar, so we only consider the weak kernel in this paper, that is,(1.8)F(s)=αe-αs,α>0,
where α reflects the mean delay of the weak kernel.The purpose of this paper is to discuss the stability and bifurcation of system (1.5), which is an extension of the aforementioned systems. By taking the discrete delay τ as the bifurcation parameter, we will show that the equilibrium of system (1.5) loses its stability and Hopf bifurcation occurs when τ passes through a certain critical value.The remainder of this paper is organized as follows. In Section2, the linear stability of system (1.5) is discussed and some sufficient conditions for the existence of Hopf bifurcations are derived. The properties of Hopf bifurcation are analyzed in detail by using the center manifold theorem and the normal form theory in Section 3. In Section 4, some numerical simulations are performed to illustrate and verify the theoretical analysis. Finally, conclusions are drawn in Section 5.


## 2. Linear Stability and Existence of Hopf Bifurcation

In this section, we discuss the linear stability of the equilibriumO(0,0) of system (1.5) and the existence of Hopf bifurcations. For analysis convenience, we define the following variable:(2.1)z(t)=∫-∞tαe-α(t-s)y(s)ds.
Then by the linear chain trick technique, system (1.5) can be transformed into the following system with only discrete time delay:(2.2)ẋ(t)=z(t)-f(x(t-τ)),ẏ(t)=-x(t-τ),ż(t)=αy(t)-αz(t).
It is obvious that system (2.2) has a unique equilibrium O(0,0,0).The linearization of system (2.2) at the equilibrium O(0,0,0) is(2.3)ẋ(t)=z(t)-ax(t-τ),ẏ(t)=-x(t-τ),ż(t)=αy(t)-αz(t).The associated characteristic equation of (2.3) is(2.4)det⁡(λ+ae-λτ0-1e-λτλ00-αλ+α)=0,
which is equivalent to(2.5)λ3+αλ2+(aλ2+aαλ+α)e-λτ=0.In the following, we investigate the distribution of roots of (2.5) and obtain the conditions under which system (2.2) undergoes Hopf bifurcation.We know thatiω(ω>0) is a root of (2.5) if and only if ω satisfies(2.6)-ω3i-αω2+(-aω2+aαωi+α)(cos⁡⁡ωτ-isin⁡ωτ)=0.
Separating the real and imaginary parts, yields(2.7)αω2=(-aω2+α)cos⁡⁡ωτ+aαωsin⁡ωτ,-ω3=(-aω2+α)sin⁡ωτ-aαωcos⁡⁡ωτ.
Taking square on both sides of the equations in system (2.7) and adding them up yield(2.8)ω6+(α2-a2)ω4+(2aα-a2α2)ω2-α2=0.
Let z=ω2, p=α2-a2, q=2aα-a2α2, and r=-α2. Then (2.8) becomes(2.9)z3+pz2+qz+r=0.Denote(2.10)h(z)=z3+pz2+qz+r.
Since α>0, then r=-α2<0, and limz→∞h(z)=∞, we can conclude that (2.9) has at least one positive root.Without loss of generality, we assume that (2.9) has three positive roots, defined by z1, z2, and z3, respectively. Then (2.8) has three positive roots(2.11)ω1=z1,ω2=z2,ω3=z3.
From (2.7), we have that(2.12)cos⁡⁡ωτ=α2ω2(aω2-α)2+a2α2ω2.
Thus, if we denote(2.13)τk(j)=1ωk{cos⁡-1(α2ωk2(aωk2-α)2+a2α2ωk2)+2jπ},
where k=1,2,3, and j=0,1,…, then ±iωk is a pair of purely imaginary roots of (2.5) with τk(j). Define(2.14)τ0=τk0(0)=min⁡k∈{1,2,3}{τk(0)},ω0=ωk0.In order to further investigate (2.5), we need to introduce a result proposed by Ruan and Wei [27], which is stated as follows.Lemma 2.1 (see [27]).
Consider the exponential polynomial(2.15)P(λ,e-λτ1,…,e-λτm)=λn+p1(0)λn-1+⋯+pn-1(0)λ+pn(0)+[p1(1)λn-1+⋯+pn-1(1)λ+pn(1)]e-λτ1+⋯+[p1(m)λn-1+⋯+pn-1(m)λ+pn(m)]e-λτm,
where τi≥0(i=1,2,…,m) and pj(i)(i=0,1,…,m;j=1,2,…,n) are constants. As (τ1,τ2,…,τm) vary, the sum of the order of the zeros of P(λ,e-λτ1,…,e-λτm) on the right half plane can change only if a zero appears on or crosses the imaginary axis.By applying Lemma2.1, one can easily obtain the following result on the distribution of roots of (2.5).Lemma 2.2.
For the third-degree transcendental (2.5), if r<0, then all roots with positive real parts of (2.5) have the same sum to those of the polynomial (2.5) for τ∈[0,τ0).Let(2.16)λ(τ)=ξ(τ)+iω(τ)
be the root of (2.5) near τ=τk(j), satisfying(2.17)ξ(τk(j))=0,ω(τk(j))=ωk.
Then the following transversality condition holds.Lemma 2.3.
Suppose thatzk=wk2 and  h′(zk)≠0, where  h(z) is defined by (2.10); then
(2.18)d(Reλ(τk(j)))dτ≠0,
and the sign of d(Reλ(τk(j)))/dτ is consistent with that of h′(zk).Proof.
The proof is similar to those in [28–32], so we omit it here.
Whenτ=0, (2.5) becomes
(2.19)λ3+(a+α)λ2+aαλ+α=0.
According to the Routh-Hurwitz criterion, if the following conditions:
(2.20)a+α>0,a(a+α)>1,
hold, then all roots of (2.19) have negative real parts, which means that the equilibrium O(0,0,0) of system (2.19) is stable.By applying Lemmas2.2 and 2.3 to (2.5), we have the following theorem.Theorem 2.4.
Letτk(j),  ω0,  τ0 and be defined by (2.13) and (2.14), respectively. Suppose that  a+α>0 and  a(a+α)>1. Then one has the following.  (i)
Ifr<0, then the equilibrium O(0,0,0) of system (2.2) is asymptotically stable for τ∈[0,τ0).(ii)
Ifr<0 and h′(zk)≠0, then system (2.2) undergoes a Hopf bifurcation at its equilibrium O(0,0,0) when τ=τk(j).


## 3. The Properties of Hopf Bifurcation

In the previous section, we obtain some conditions for Hopf bifurcations to occur at the critical valueτk(j). In this section, we analyze the properties of the Hopf bifurcation by virtue of the method proposed by Hassard et al. [33], namely, to determine the direction of Hopf bifurcation and the stability of bifurcating periodic solutions bifurcating from the equilibrium O(0,0,0) for system (2.2) by applying the normal form theory and the center manifold theorem.For analysis convenience, lett=sτ, x1(s)=x(sτ), y1(s)=y(sτ), z1(s)=z(sτ), and τ=τk(j)+μ, μ∈R. Denote t=s; then system (2.2) is transformed into the following form:(3.1)ẋ1(t)=(τk(j)+μ)(z1(t)-f(x1(t-1))),ẏ1(t)=(τk(j)+μ)(-x1(t-1)),ż1(t)=(τk(j)+μ)(αy1(t)-αz1(t)),
where f(x1(t-1))=ax1(t-1)+bx13(t-1).Its linear part is(3.2)ẋ1(t)=(τk(j)+μ)(z1(t)-ax1(t-1)),ẏ1(t)=(τk(j)+μ)(-x1(t-1)),ż1(t)=(τk(j)+μ)(αy1(t)-αz1(t)),
and the nonlinear part is as follows:(3.3)f(μ,xt)=(τk(j)+μ)(-bx13(t-1)00).
Denote Ck[-1,0]={ϕ|ϕ:[-1,0]→R3, each component of ϕ has k-order continuous derivative}. For ϕ(θ)=(ϕ1(θ),ϕ2(θ),ϕ3(θ))T∈C0[-1,0], we define(3.4)Lμϕ=(τk(j)+μ)(0010000α-α)ϕ(0)+(τk(j)+μ)(-a00-100000)ϕ(-1),
where Lμ is a one-parameter family of bounded linear operators in C0[-1,0]→R3. By the Riesz representation theorem, there exists a 3×3 matrix whose components are bounded variation functions η(θ,μ) in [-1,0]→R3 such that(3.5)Lμϕ=∫-10dη(θ,μ)ϕ(θ).
In fact, we can choose(3.6)η(θ,μ)=(τk(j)+μ)(0010000α-α)δ(θ)-(τk(j)+μ)(-a00-100000)δ(θ+1),
where δ is the Dirac delta function, which is defined by δ(θ)={0,θ≠0,1,θ=0.Forϕ∈C1([-1,0],R3), define(3.7)A(μ)ϕ={dϕ(θ)dθ,θ∈[-1,0),∫-10dη(μ,s)ϕ(s),θ=0,R(μ)ϕ={0,θ∈[-1,0),f(μ,ϕ),θ=0.
Then, we can transform system (2.2) into an operator equation as the following form:(3.8)ẋt=A(μ)xt+R(μ)xt,
where xt(θ)=(x1(t+θ),y1(t+θ),z1(t+θ))T for θ∈[-1,0].Forψ∈C1([0,1],(R3)*), define(3.9)A*ψ(s)={-dψ(s)ds,s∈(0,1],∫-10ψ(-t)dη(t,0),s=0
and a bilinear form(3.10)〈ψ(s),ϕ(θ)〉=ψ¯(0)ϕ(0)-∫-10∫ξ=0θψ¯(ξ-θ)dη(θ)ϕ(ξ)dξ,
where η(θ)=η(θ,0). Then A(0) and A* are adjoint operators. By the discussion in Section 2, we know that ±iωkτk(j) are eigenvalues of A(0), and they are also eigenvalues of A*.Suppose thatq(θ)=(1,β,γ)Teiθωkτk(j)is the eigenvector of A(0) corresponding to iτk(j)ωk; then A(0)q(θ)=iτk(j)ωkq(θ). It follows from the definition of A(0) and (3.5) and (3.6) that(3.11)τk(j)(0010000α-α)(1βγ)+τk(j)(-a00-100000)(e-iωkτk(j)βe-iωkτk(j)γe-iωkτk(j))=iωkτk(j)(1βγ).
Thus, we can easily obtain(3.12)β=(α+iωk)(iωk+ae-iωkτk(j))α,γ=iωk+ae-iωkτk(j).
On the other hand, suppose that q*(s)=B(1,β*,γ*)Teisωkτk(j)is the eigenvector of A* corresponding to -iωkτk(j). By the definition of A* and (3.5) and (3.6), we have(3.13)τk(j)(00000α10-α)(1β*γ*)+τk(j)(-a-10000000)(eiωkτk(j)β*eiωkτk(j)γ*eiωkτk(j))=-iωkτk(j)(1β*γ*).
Therefore, we have(3.14)β*=αiωk(iωk-α),γ*=1α-iωk.
In order to assure 〈q*(s),q(θ)〉=1, we need to determine the value of B. From (3.10), we have that(3.15)〈q*(s),q(θ)〉=q*¯(0)q(0)-∫-10∫ξ=0θq*¯(ξ-θ)dη(θ)q(ξ)dξ=B¯(1,β*¯,γ*¯)(1,β,γ)T-∫-10∫ξ=0θB¯(1,β*¯,γ*¯)e-i(ξ-θ)ωkτk(j)dη(θ)(1,β,γ)Teiξωkτk(j)dξ=B¯{1+ββ*¯+γγ*¯-∫-10(1,β*¯,γ*¯)θeiθωkτk(j)dη(θ)(1,β,γ)T}=B¯(1+ββ*¯+γγ*¯-τk(j)(a+β)e-iωkτk(j)).
Thus, we can choose B as(3.16)B=11+β¯β*+γ¯γ*-τk(j)(a+β¯)eiωkτk(j).
Similarly, we can get 〈q*(s),q¯(θ)〉=0.Using the same notation as in Hassard et al. [33], we compute the coordinates to describe the center manifold C0 at μ=0. Let xt be the solution of (3.8) when μ=0.Define(3.17)z(t)=〈q*,xt〉,W(t,θ)=xt(θ)-2Re{z(t)q(θ)}.
On the center manifold C0 we have that(3.18)W(t,θ)=W(z,z¯,θ),
where(3.19)W(z,z¯,θ)=W20(θ)z22+W11(θ)zz¯+W02(θ)z¯22+⋅⋅⋅,
in which z and z¯ are local coordinates for center manifold C0 in the direction of q* and q¯*, respectively. Note that W is real if xt is real. We consider only real solutions.For the solutionxt∈C0 of (3.8), since μ=0, we have that(3.20)ż(t)=iτk(j)ωkz+〈q¯*(θ),f(0,W(z,z¯,θ)+2Re{zq(θ)})〉=iτk(j)ωkz+q¯*(0)f(0,W(z,z¯,0)+2Re{zq(0)})=iτk(j)ωkz+q¯*(0)f0(z,z¯).
We rewrite this equation as(3.21)ż(t)=iτk(j)ωkz(t)+g(z,z¯),
where(3.22)g(z,z¯)=q¯*(0)f0(z,z¯)=g20z22+g11zz¯+g02z¯22+g21z2z¯2+⋅⋅⋅.
Noting that xt(θ)=W(t,θ)+zq(θ)+z¯q¯(θ) and q(θ)=(1,β,γ)Teiθωkτk, we have that(3.23)x1(-1)=z+z¯+W20(1)(0)z22+W11(1)(0)zz¯+W02(1)(0)z¯22+O(|(z,z¯)|3).According to (3.20) and (3.21), we know that(3.24)g(z,z¯)=q¯*(0)f(z,z¯)=τk(j)B¯(1,β*¯,γ*¯)(-bx13(t-1)00),
where(3.25)x1(t+θ)=W(1)(t,θ)+z(t)q(1)(θ)+z¯(t)q¯(1)(θ).
From (3.20) and (3.24), we have that(3.26)g(z,z¯)=-τk(j)bB¯x13(t-1)=-τk(j)bB¯[W(1)(t,θ)+z(t)q(1)(θ)+z¯(t)q¯(1)(θ)]3=-τk(j)bB¯[W20(1)(-1)z22+W11(1)(-1)zz¯+W02(1)(-1)z¯22+z(t)q(1)(-1)+z¯(t)q¯(1)(-1)]3.
Comparing the coefficients in (3.26) with those in (3.22), we have that(3.27)g20=0,g11=0,g02=0,g21=-6τk(j)bB¯(q(1)(-1))2q¯(1)(-1).
Thus, we can calculate the following values:(3.28)c1(0)=i2τk(j)ωk(g11g20-2|g11|2-|g02|23)+g212,μ2=-Re{c1(0)}Re{λ′(τk(j))},β2=2Re{c1(0)},t2=-Im⁡{c1(0)}+μ2Im⁡{λ′(τk(j))}ωkτk(j),
which we need to investigate the properties of Hopf bifurcation. According to [33], we know that μ2 determines the direction of the Hopf bifurcation: if μ2>0(μ2<0), then the Hopf bifurcation is supercritical (subcritical) and the bifurcating periodic solutions exist for τ>τk(j)(τ<τk(j)); β2 determines the stability of the bifurcating periodic solutions: the bifurcating periodic solutions are stable (unstable) if β2<0(β2>0); t2 determines the period of the bifurcating periodic solutions: the period increases (decreases) if t2>0(t2<0).


## 4. A Numerical Example

In this section, we use the formulae obtained in Sections2 and 3 to verify the existence of a Hopf bifurcation and calculate the Hopf bifurcation value and the direction of the Hopf bifurcation of system (2.2) with α=1, a=0.9, and b=2.By the results in Section2, we can determine that(4.1)z1=0.6506,ω0=0.8066,τ0=0.4628.
It follows from (3.28) that(4.2)c1(0)=-1.2260-0.4739i,μ2=5.8875,β2=-2.452,t2=-1.9408.In light of Theorem2.4, the equilibrium O(0,0,0) of system (2.2) is stable when τ<τ0. This is illustrated in Figure 1 with τ=0.37. Since μ2>0, when τ passes through the critical value τ0=0.4628, the equilibrium O(0,0,0) loses its stability and a Hopf bifurcation occurs, that is, periodic solutions bifurcate from the equilibrium O(0,0,0). The individual periodic orbits are stable since β2<0. Figure 2 shows that there are stable limit cycles for system (2.2) with τ=0.5. Since t2<0, the period of the periodic solutions decreases as τ increases. For τ=0.55, the phase portrait and the waveform portraits are shown in Figure 3. We can tell from Figures 2 and 3 that the period of τ=0.55 is slightly smaller than that of τ=0.5.Phase portrait and waveform portraits of system (2.2) with τ=0.37.
(a)(b)(c)(d)Phase portrait and waveform portraits of system (2.2) with τ=0.5.
(a)(b)(c)(d)Phase portrait and waveform portraits of system (2.2) with τ=0.55.
(a)(b)(c)(d)


## 5. Conclusions

The van der Pol equation with discrete and distributed delays is investigated in this paper. Sufficient conditions on the linear stability of this van der Pol equation have been obtained by analyzing the associated transcendental characteristic equation. By choosing the discrete time delay as a bifurcation parameter, we have shown that this equation undergoes a sequence of Hopf bifurcations. In addition, formulae for determining the direction of Hopf bifurcation and the stability of bifurcating periodic solutions are derived. Simulation results have verified and demonstrated the correctness of the theoretical analysis.



---

*Source: 569141-2011-05-12.xml*