Zhuang and Jia SpringerPlus (2016)5:1303 DOI 10.1186/s40064-016-2958-y
Open Access
RESEARCH
Delay‑induced periodic phenomenon in a diffusive regulated logistic model Kejun Zhuang1,2* and Gao Jia3 *Correspondence: [email protected] 1 Business School, University of Shanghai for Science and Technology, Shanghai 200093, China Full list of author information is available at the end of the article
Abstract The diffusive logistic growth model with time delay and feedback control is considered. First, the well-posedness and permanence of solutions are discussed by using some comparison techniques. Then, the sufficient conditions for stability of nonnegative constant steady states are established, and the occurrence of Hopf bifurcation at positive steady state is performed. Next, the bifurcation properties are derived by computing the normal form on center manifold. Our results not only supplement but also generalized some existing ones. Finally, some numerical simulations show the feasibility of our theoretical analyses. Keywords: Logistic model, Feedback control, Hopf bifurcation, Reaction–diffusion system, Delay Mathematics Subject Classification: 35K57, 35B10, 35B32
Background The classic logistic model N (t) dN (t) = rN (t) 1 − , dt K
r, k ∈ (0, +∞)
(1)
was first proposed by Verhulst in 1838. It can be utilized to describe the single–species growth and has been the basis of varieties of models in population ecology and epidemiology. For system (1) and its generalized forms, the significant results involve the asymptotic properties (Berezansky et al. 2004; Röst 2011), permanence and stability (Fan and Wang 2010; Chen et al. 2006), periodicity (Sun and Chen 2007) and almost periodicity (Yang and Yuan 2008) of solutions, Hopf bifurcation (Sun et al. 2007; Song and Yuan 2007; Song and Peng 2006; Chen and Shi 2012), traveling wave front (Zhang and Sun 2014), free boundary problem (Gu and Lin 2014), and so on. In addition, the Hopf bifurcation analyses for some diffusive predator–prey systems were also done (see Yang 2015; Yang and Zhang 2016a, b). In particular, Gopalsamy (1993) considered the controlled delay system in the following form
dN (t) a1 N (t)+a2 N (t−τ ) K dt = rN (t) 1 − du(t) = bN (t − τ ) − au(t), dt
− cu(t) ,
(2)
© 2016 The Author(s). This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.
Zhuang and Jia SpringerPlus (2016)5:1303
Page 2 of 15
where all the coefficients and time delay τ are positive constants, N(t) is the number of individuals at time t, and variable u(t) denotes an indirect control variable (see Aizerman and Gantmacher 1964; Lefschetz 1965). They have derived the sufficient conditions to guarantee that the positive equilibrium solution is globally asymptotical stable. Strictly speaking, spatial diffusion can not be ignored in studying the natural biological system (Murray 2003; Ghergu and Radulescu 2012). In the real world, most populations are moving and the densities are dependent of time and space. Therefore, diffusion should be taken into account in studying the basic logistic equation. However, there have been very few results on the influence of time delay on the reaction–diffusion logistic model with feedback control. Inspired by the previous discussions, we mainly consider the reaction–diffusion system as follows:
∂N (x,t) = d1 ∆N (x, t) + rN (x, t) 1 − a1 N (x,t)+aK2 N (x,t−τ ) ∂t ∂u(x,t) = d2 ∆u(x, t) + bN (x, t − τ ) − au(x, t), ∂t
− cu(x, t) ,
(3)
where (x, t) ∈ Ω × [0, +∞), Ω = (0, lπ ). The model (3) is considered with the initial value conditions as follows
N (x, t) = η1 (x, t) ≥ 0,
u(x, t) = η2 (x, t) ≥ 0,
x ∈ [0, lπ] × [−τ , 0].
(4)
We also assume that the model (3) is closed and there is no emigration or immigration across the boundary. Hence, the boundary conditions are considered as
∂u ∂N = = 0, ∂ν ∂ν
(x, t) ∈ ∂Ω × [0, +∞),
(5)
where ∂/∂ν represents the outward normal derivative on the boundary ∂Ω. In this paper, we develop a reaction–diffusion logistic model with time delay and diffusion, which makes up perfectly for the deficiencies of the previous literatures. The main objective is to explore the dynamics of system (3) by regarding τ as the bifurcation parameter. The structure of this paper is arranged as follows. In section “Preliminaries”, we derive the well–posedness of solutions and the permanence of the system. In section “Occurrence of the Hopf bifurcation”, we establish the existence of Hopf bifurcation. In section “Bifurcation properties”, we get the formulae for determining the Hopf bifurcation properties. In section “Numerical simulations”, we illustrate our theoretical results by some numerical simulations. Finally, we give some discussions and conclusions.
Preliminaries As we know, spatial diffusion and time delay do not change the number and locations of constant equilibria because of no-flux boundary conditions. Then system (3) has two nonnegative equlibria E0 = (0, 0) and E ∗ = (N ∗ , u∗ ), where N∗ =
aK , a(a1 + a2 ) + bcK
u∗ =
b ∗ bK . N = a a(a1 + a2 ) + bcK
Zhuang and Jia SpringerPlus (2016)5:1303
Page 3 of 15
Well–posedness of solutions
Here, for problem (3)–(5), we devote ourselves to the existence, uniqueness, nonnegativity and boundedness of solutions. Theorem 1 For any given initial data satisfying the conditions (4) and boundary conditions (5), system (3) has a unique global solution of system and the solution maintains nonnegative and uniformly bounded for all t ≥ 0. Proof Using the similar methods in Hattaf and Yousfi (2015), Hattaf and Yousfi (2015), we can get the local existence and uniqueness of solution (N(x, t), u(x, t)) with x ∈ Ω¯ and t ∈ [0, T ), where T is the maximal existence time of solution. It is easy to find that 0 = (0, 0) and M = (M1 , M2 ) are a pair of coupled lower–upper solutions to problem (3)–(5), where
K M1 = max , sup �ϕ1 (·, s)�C(Ω, ¯ R) , a1 −τ ≤s≤0 bM1 , sup �ϕ2 (·, s)�C(Ω, M2 = max ¯ R) . a −τ ≤s≤0 By means of the comparison theorem, we can obtain that 0 ≤ N (x, t) ≤ M1 and 0 ≤ u(x, t) ≤ M2 for x ∈ Ω¯ and t ∈ [0, T ). It is obvious that the upper bound of solution is independent of the maximal existence interval [0, T). It follows from the standard theory for semilinear parabolic systems (Wu 1996; Henry 1993) that the solution glob ally exists. The proof is complete Dissipativeness and permanence
In the following, we will show that system (3) is permanent, which means that any nonnegative solution of (3) is bounded as t → +∞ for all x ∈ Ω. Theorem 2 (Dissipativeness) The nonnegative solution (N, u) of system (3) satisfies
lim sup N (x, t) ≤ t→+∞
K bK , lim sup u(x, t) ≤ . a1 t→+∞ aa1
Proof Based on the first equation in system (3), we get
a1 ∂N (x, t) − d1 ∆N (x, t) ≤ rN (x, t) 1 − N (x, t) ∂t K
for
(x, t) ∈ Ω × [0, +∞).
Then from the standard comparison principle of parabolic equations, we can easily get
lim sup N (x, t) ≤ t→+∞
K . a1
For an arbitrary ε1 > 0, we could get a positive constant T1 such that for any t ≥ T1,
N (x, t) ≤
K + ε1 . a1
Zhuang and Jia SpringerPlus (2016)5:1303
Page 4 of 15
Thus, for any T ∈ [T1 + τ , +∞), we have K ∂u(x, t) − d2 ∆u(x, t) ≤ b + ε1 − au(x, t). ∂t a1 This implies
lim sup u(x, t) ≤ t→+∞
bK aa1
by comparison principle of parabolic equations and the arbitrariness of ε1.
Theorem 3 If aa1 > aa2 + bcK , then system (3) is permanent. Proof From Theorem 2, for an arbitrary ε2 > 0, we can find a constant T > T1 + T2, such that
u(x, t) ≤
bK + ε2 aa1
in Ω × [T2 , +∞). Moreover, we can obtain
∂N (x, t) − d1 ∆N (x, t) ∂t bK a2 K a1 + ε1 − c + ε2 − N (x, t) , ≥ rN (x, t) 1 − K a1 aa1 K the comparison principle shows that
lim inf N (x, t) ≥ t→+∞
K aa1 − aa2 − bcK >0 a1 aa1
due to the continuity as ε1 → 0 and ε2 → 0. Similarly, we can also have
lim inf u(x, t) ≥ t→+∞
bK aa1 − aa2 − bcK > 0. aa1 aa1
Combining the results in Theorem 2, we can easily conclude that system (3) is permanent.
Occurrence of the Hopf bifurcation For system (3), we shall study the local stability of two constant steady states and the occurrence of Hopf bifurcation phenomenon through discussing the distribution of characteristic values. Denote u1 (t) = N (x, t),
u2 (t) = u(x, t),
U (t) = (u1 (t), u2 (t))T .
By defining the phase space C = C([−τ , 0], X), we can rewritten system (3) as the semilinear functional differential equation:
Zhuang and Jia SpringerPlus (2016)5:1303
Page 5 of 15
U˙ (t) = D∆U (t) + G(Ut ),
(6)
where X = {(u, v) ∈ H 2 (0, lπ ) × H 2 (0, lπ )|ux = vx = 0, x = 0, lπ}, Ut (·) = U (t + ·), D = diag {d1 , d2 }, ∆ = diag {∂ 2 /∂x2 , ∂ 2 /∂x2 }, and G(Ut ) : C → X is defined by
G(Ut ) =
ru1 (t) 1 − a1 u1 (t)+aK2 u1 (t−τ ) − cu2 (t) . bu1 (t − τ ) − au2 (t)
The linear system of (6) at E0 (0, 0) is
U˙ (t) = D∆U (t) + LE0 (Ut ),
(7)
where
LE0 (ϕ) =
rϕ1 (0) bϕ1 (−τ )
0 − aϕ2 (0)
for ϕ(θ) = Ut (θ), ϕ = (ϕ1 , ϕ2 )T ∈ C. The characteristic equation of (7) is
y − D∆y − LE0 (e ·y ) = 0,
(8)
where y ∈ dom (∆)\{0}, dom ∆ ⊂ X and e · (θ)y = e θ y for θ ∈ [−τ , 0]. We know that the operator ∆ in Ω with homogeneous Neumann boundary condition has the eigenvalues −n2 /l 2 and the corresponding eigenfunctions cos(nx/l), n ∈ N0 = {0, 1, 2, . . .}. By using the Fourier expansion in (8),
∞ αn y= cos(nx/l), γn n=0
where αn, γn ∈ C. Therefore, the characteristic equation (8) can be transferred into + d n2 − r 0 1 l2 = 0, n ∈ N0 . 2 −be− τ + d2 nl 2 + a
We then obtain the characteristic values as follows
1,n = −d1
n2 + r, l2
2,n = −d2
n2 − a, l2
n ∈ N0 .
It is obvious that 1,0 = r > 0, and we can establish the instability of E0. Theorem 4 The trivial equilibrium E0 of system (3) is always unstable. Next, we will focus on the occurrence of Hopf bifurcation phenomenon. Linearizing system (3) at E ∗ = (N ∗ , u∗ ), we get
U˙ (t) = D∆U (t) + L(Ut ), where L : C → X is given by ra ∗ − K1 N ϕ1 (0) − raK2 N ∗ ϕ1 (−τ ) L(ϕ) = bϕ1 (−τ )
(9)
− cN ∗ ϕ2 (0) − aϕ2 (0)
Zhuang and Jia SpringerPlus (2016)5:1303
Page 6 of 15
with ϕ(θ) = Ut (θ), ϕ = (ϕ1 , ϕ2 )T ∈ C. Similar to the previous discussion, we can obtain the characteristic equation
2 + An + Bn + e− τ (C + Dn ) = 0,
n ∈ N0 ,
(10)
where
ra1 ∗ n2 N > 0, +a+ 2 l K n4 ra1 ∗ n2 raa1 ∗ Bn = d1 d2 4 + a + N d1 2 + N > 0, l K l K ra2 ∗ N > 0, C= K ra2 ∗ n2 raa2 ∗ Dn = N d1 2 + N + bcN ∗ > 0. K l K An = (d1 + d2 )
For τ = 0, Eq. (10) can be reduced to
2 + (An + C) + Bn + Dn = 0 with An + C > 0 and Bn + D > 0. On the basis of Routh–Hurwitz stability criterion, we can obtain the local stability of E ∗ when τ = 0. Lemma 1 The positive equilibrium is always locally asymptotically stable without time delay. Remark 1 From Lemma 1, we can find that there is no Turing instability without time delay. For τ � = 0, let us suppose that = iω(ω > 0) satisfies Eq. (10). First, plugging = iω into Eq. (10) and then segregating the real and imaginary components with the help of Euler’s formula, we can get the following two equations of ω ω2 − Bn = Dn cos ωτ + Cω sin ωτ , −ωAn = Cω cos ωτ − Dn sin ωτ . Second, solving these equations, we can obtain 2 −B D n n cos ωτ = (Dn −A2n C)ω , C ω2 +D2 n
sin ωτ =
Cω3 +(An Dn −Bn C)ω . C 2 ω2 +Dn2
(11)
Third, squaring both sides of those two equations and then adding them up, we get the following equation
ω4 + A2n − 2Bn − C 2 ω2 + Bn2 − Dn2 = 0,
(12)
Zhuang and Jia SpringerPlus (2016)5:1303
Page 7 of 15
where
ra n2 n4 r 2 a21 − a22 ∗2 1 ∗ N d1 + ad2 2 + N , A2n − 2Bn − C 2 = d12 + d22 4 + 2 l K l K2 n4 ra1 ∗ n2 raa1 ∗ ra2 ∗ n2 N d1 2 + N − N d1 2 Bn2 − Dn2 = (Bn + Dn ) d1 d2 4 + a + l K l K K l raa2 ∗ ∗ − N − bcN . K Lemma 2 For τ > 0, we have (i) If a1 > a2 + bcK ar , then Eq. (10) does not have purely imaginary root. (ii) If a2 < a1 < a2 + bcK ar , then there exists N0 ∈ N0, such that Eq. (10) does not have purely imaginary root when n > N0, and has a pair of conjugate purely imaginary eigenvalues when 0 ≤ n ≤ N0. Proof We can easily verify that A2n − 2Bn − C 2 > 0 and Bn2 − Dn2 > 0 when a1 > a2 + bcK ar . This means that Eq. (12) has no positive root. In other words, there could be no purely imaginary root in Eq. (10) for any τ > 0. 2 2 On the contrary, if a2 < a1 < a2 + bcK ar , then B0 − D0 < 0 and there exists N ∈ N0 such that 2 Bn − Dn2 < 0, n = 0, 1, 2, . . . , N0 , Bn2 − Dn2 ≥ 0, n = N0 + 1, N0 + 2, . . . .
That is to say, Eq. (12) has no positive root when n > N0 and has the unique positive root ωn when 0 ≤ n ≤ N0, where
ωn =
−(A2n − 2Bn − C 2 ) +
1 (A2n − 2Bn − C 2 )2 − 4(Bn2 − Dn2 ) 2 . 2
By direct computation, we have
2 ra2 ∗ 2 n4 n 2ra2 ∗ ∗ ∗ N d2 4 + bcN d1 + N d2 + bcN d2 2 An Dn − Bn C = K l K l +
ra1 bc ∗2 ra2 a2 ∗ N + abcN ∗ + N K K
> 0. Moreover, Eq. (10) has characteristic values ±iωn with (n)
τj
(n)
= τ0 +
2jπ , ωn
0 ≤ n ≤ N0 ,
j = 0, 1, 2, . . . ,
where (n)
τ0
=
(Dn − An C)ωn2 − Bn Dn 1 arccos . ωn Dn2 + C 2 ωn2
Zhuang and Jia SpringerPlus (2016)5:1303
Page 8 of 15
This completes the proof. We now check the transversality condition. Lemma 3 If a2 < a1 < a2 + n ∈ {0, 1, 2, . . . , N0 } .
bcK ar ,
then
dRe ( ) dτ τ =τ (n)
for
>0
j
j ∈ N0
and
Proof By taking the derivatives on both sides of (10) with respect to τ, we can get d d d d + An + Ce− τ + (C + Dn )e− τ − − τ = 0, 2 dτ dτ dτ dτ and
d dτ
−1
=
2 + An + Ce− τ − τ e− τ (C + Dn ) e− τ (C + Dn )
=
τ (2 + An )e τ + C − . (C + Dn )
On the basis of (11) and (12), we get
d dτ
−1
(n)
(n) (n) +C (2iωn + An ) cos ωn τj + i sin ωn τj
=
iωn (iCωn + Dn )
τ =τj
(n)
C + An cos ωn τj
=
(n)
− 2ωn sin ωn τj
−Cωn2 + iDn ωn +
(n)
i 2ωn cos ωn τj
(n)
+ An sin ωn τj
−Cωn2 + iDn ωn
(n)
−
τj
iωn
(n)
τj − iωn .
Further simplification will lead to
Re
d dτ
−1
(n)
τ =τj
(n) (n) Dn ωn 2ωn cos ωn τj + An sin ωn τj = 2 Cωn2 + (Dn ωn )2 (n) (n) Cωn2 C + An cos ωn τj − 2ωn sin ωn τj − 2 Cωn2 + (Dn ωn )2
=
ωn4 + Dn2 − Bn2 2 Cωn2 + (Dn ωn )2
> 0.
The proof is complete. According to Lemmas 1–3 and the Hopf bifurcation theory developed by Wu (1996), the following conclusions can be drawn. Theorem 5 Define
τ0 =
min
n∈{0,1,2,...,N0 },j∈N0
(n)
τj
.
Zhuang and Jia SpringerPlus (2016)5:1303
Page 9 of 15
(i) If a1 > a2 + bcK , then for any τ > 0, the positive equilibrium E ∗ is always locally ar asymptotically stable. ∗ (ii) If a2 < a1 < a2 + bcK ar , then E is locally asymptotically stable when τ ∈ [0, τ0 ), and is unstable when τ ∈ (τ0 , +∞). (n) (iii) System (3) has a Hopf bifurcation from E ∗ at τj with n ∈ {0, 1, 2, . . . , N0 } and j ∈ N0. If n = 0, the periodic solutions bifurcating positive equilibrium are all spatially homogeneous. Otherwise, these bifurcating periodic solutions are spatially inhomogeneous.
Bifurcation properties In Theorem 5, we have demonstrated that there exist some spatially homogeneous or inhomogeneous periodic solutions when time delay crosses through some particular values. We are now in the position to investigate the bifurcation properties. (n) with j ∈ N0 and In general, we use τ ∗ to denote an arbitrary value of τj ∗ n ∈ {0, 1, 2, . . . , N0 }. And we also use ±iω to denote the corresponding simply purely imaginary roots ±iωn. ˜ t)), and τ = τ ∗ + α with ˜ t) = u(·, τ t), U˜ (t) = (N˜ (·, t), u(·, Set N˜ (·, t) = N (·, τ t), u(·, α ∈ R. For simplicity we drop the tilde and rewrite system (3) as follows, dU (t) = τ D∆U (t) + L(α)(Ut ) + f (Ut , α), dt
(13)
where ϕ = (ϕ1 , ϕ2 )T ∈ C, L(α)(·) : C → X and f : C × R → X are respectively denoted by ra ∗ − K1 N ϕ1 (0) + raK2 N ∗ ϕ1 (−1) − cN ∗ ϕ2 (0) L(α)(ϕ) = (τ ∗ + α) bϕ1 (−1) − aϕ2 (0) and
2ra − K 1 ϕ12 (0) − f (ϕ, α) = (τ ∗ + α)
ra2 K ϕ1 (0)ϕ1 (−1) − rcϕ1 (0)ϕ2 (0)
0
.
Note that α = τ − τ ∗, we can find that system (13) may causes a Hopf bifurcation when α = 0. For the following linear differential equation:
U˙ (t) = τ D∆U (t) + L(α)(Ut ),
(14)
we can easily deduce that the corresponding characteristic equation has characteristic values ± i ω∗ τ ∗ when α = 0. Next, we discuss the following differential equation:
Y˙ (t) = −τ Dn2 Y (t) + L(α)(Yt ).
(15)
We can use Riesz representation theorem here, which tells us that there is a 2 × 2 matrix function η(θ, α) (−1 ≤ θ ≤ 0) with bounded variation elements satisfying 0 n2 −τ D 2 ϕ(0) + L(α)(ϕ) = d[η(θ, α)]ϕ(θ), l −1
Zhuang and Jia SpringerPlus (2016)5:1303
Page 10 of 15
where
� � 2 −d1 nl 2 − raK1 N ∗ − cN ∗ ∗ , (τ + α) 2 0 − d2 nl 2 − a η(θ, α) = 0, � ra ∗ � − K2 N 0 ∗ , (τ + α) b 0
θ = 0, θ ∈ (−1, 0), θ = −1.
For Φ ∈ C 1 ([−1, 0], R2 ), Ψ ∈ C 1 ([0, 1], R2 ), we define ) dΦ(θ θ ∈ [−1, 0), dθ , A1 (Φ(θ)) = � 0 [dη0 (θ)]Φ(θ), θ = 0, −1 dΨ − ds(s) , s ∈ (0, 1], ∗ A1 (Ψ (s)) = � 0 s = 0. −1 [dη0 (θ)]Ψ (−θ ),
Then the formal adjoint, A∗1, of A1 is given by
(Ψ , Φ)0 = Ψ (0)Φ(0) −
0
−1
= Ψ (0)Φ(0) + τ
∗
θ
Ψ (ζ − θ )d[η(θ, 0)]Φ(ζ )dζ
ζ =0 0
ra ∗ − K2 N 0 Ψ (ζ + 1) Φ(ζ )dζ . b 0 −1
∗
∗
By calculation, we can find that q(θ) = (1, ξ )T eiω θ τ and q ∗ (s) = M(1, η)eiω eigenvectors of A1 and A∗1 associated with iω∗ τ ∗, respectively, where
θ ∈ [−1, 0],
∗ sτ ∗
are
s ∈ [0, 1],
and ∗ ∗
be−iω τ ra2 ∗ iω∗ τ ∗ −1 cN ∗ ∗ ξ= ξ η + τ bη − N e . , η = , M = 1 + iω∗ + a iω∗ + a K Then P = span{q(θ), q(θ)}, P ∗ = span{q ∗ (s), q ∗ (s)} are the center subspace of system (3). T 2 nx T . Define h · fn = h1 βn1 + h2 βn2, fn = βn1 , βn2 and βn1 = cos nx l , 0 , βn = 0, cos l The complex-valued L2 inner product on Hilbert space XC are lπ 1 �U1 , U2 � = (u1 v1 + u2 v2 )dx, (16) lπ 0 for U1 = (u1 , u2 ), U2 = (v1 , v2 ) ∈ XC . And �β0i , β0i � = 1, �βni , βni � = 12, i = 1, 2, n = 1, 2, . . .,
�Φ, fn � = �Φ, βn1 �, �Φ, βn2 � ,
(17)
where Φ ∈ C([−1, 0], X). We can establish the center subspace of system (14) at α = 0 as follows PCN L = (q(θ)z + q(θ)¯z ) · fn , z ∈ C . (18)
Zhuang and Jia SpringerPlus (2016)5:1303
Page 11 of 15
Based on the conclusions drawn by Wu (1996) and Hassard et al. (1981), the solutions of (13) are
Ut = (q(θ)z(t) + q(θ)¯z (t)) · fn + W (z(t), z¯ (t), θ ), where
W (z, z¯ , θ ) = W20
z¯ 2 z2 + W11 z z¯ + W02 + · · · . 2 2
(19)
Moreover, for Ut ∈ C0 of (13) at τ = τ ∗, we have z˙ = iω∗ τ ∗ z + g(z, z¯ ), where
g(z, z¯ ) = q ∗ (0)�f (Ut , 0), fn � = g20
z2 z¯ 2 z 2 z¯ + g11 z z¯ + g02 + g21 + ··· . 2 2 2
(20)
By (16)–(20), we can compute
g20 =
g11 =
0, −2τ ∗ M raK1 +
+ rcξ ,
0, −2τ ∗ M raK1 + rc Re {ξ } +
g02 = g20 , g21
ra2 −iω∗ τ ∗ K e
ra2 K
Re
n = 1, 2, . . . , n = 0,
∗ ∗ {eiω τ }
,
n = 1, 2, . . . , n = 0,
lπ 2ra Mτ ∗ nx 1 (1) (1) (W11 dx = −2 (0) + W20 (0)) cos2 lπ K l 0 lπ nx 1 1 (2) (1) (1) (2) ¯ dx rc W11 (0) + W20 (0) + ξ W20 (0) + ξ W11 (0) cos2 + 2 2 l 0 lπ nx ra2 −iω∗ τ ∗ (2) 1 ∗ ∗ (2) 1 (1) + dx e W11 (0) + eiω τ W20 (0) + ξ¯ W20 (−1) cos2 K 2 2 l 0 lπ ra2 (11) 2 nx ξ W11 (−1) cos dx . + K l 0
Then we should compute W20 (θ) and W11 (θ) to determine g21. Following the formulas in Wu (1996), We can obtain that ig20 ig02 ∗ ∗ W20 (θ) = q(θ) · fn + E1 e2iω τ θ , q(θ) + ω∗ τ ∗ 3ω∗ τ ∗ ig11 q(θ) ig11 · fn + E2 , W11 (θ) = − ∗ ∗ 1(θ) + ω τ ω∗ τ ∗ 2ra ∗ ∗ 2 −iω τ nx − K 1 − 2ra e − 2rcξ ′ K cos2 , E1 = E1 × 0 l −1 2 ∗ ∗ 2iω∗ + d1 nl 2 + raK1 N ∗ + raK2 N ∗ e−2iω τ cN ∗ ′ E1 = , 2 ∗ ∗ −be−2iω τ 2iω∗ + a + d2 nl 2 and
Zhuang and Jia SpringerPlus (2016)5:1303
Page 12 of 15
1 − 2ra K − 2rc Re {ξ } − 0
E2 =
E2′
×
E2′
d1 nl 2 +
=
2
ra1 ∗ K N
+
ra2 ∗ K N
2ra2 K
Re {eiω
cN ∗ 2
a + d2 nl 2
−b
−1
∗τ ∗
}
cos2
nx , l
.
From the previous expressions of g20, g11, g02 and g21, we can further compute
g21 1 2 g20 g11 − 2|g11 | − |g02 | + , c1 (0) = 2ω∗ τ ∗ 3 2 Re (c1 (0)) , µ2 = − Re ( ′ (τ ∗ )) β2 = 2 Re (c1 (0)), 1 T2 = − ∗ ∗ ( Im (c1 (0)) + µ2 Im ( ′ (τ ∗ ))). ω τ i
2
On account of preceding calculations, we arrive at the following conclusion on the bifurcation properties. Theorem 6 The bifurcation direction is supercritical if µ2 > 0, which means that the periodic solution exists for τ > τ0. On the contrary, the bifurcation direction is subcritical if µ2 < 0, which means that the periodic solution exists for τ < τ0. Moreover, the periodic solution is orbitally asymptotically stable if β2 < 0, or unstable if β2 > 0. The period of periodic solution is monotonically increasing at the time delay τ when T2 > 0, or is monotonically decreasing at the time delay τ when T2 < 0.
Numerical simulations In this section, we give some numerical examples to test the preceding results with assistance of MATLAB. For system (3), let Ω = (0, 2π ) and choose d1 = 1,
d2 = 0.5,
r = 0.6,
a = b = c = 1,
a1 = a2 = 2,
K = 1,
and the initial values N (x, 0) = 0.5 and u(x, 0) = 0.9. Then we can get the positive equilibrium E ∗ = (0.2, 0.2). By direct computation, we have N0 = 0, ω0 ≈ 0.348266, and τ0(0) ≈ 5.81966, then the Hopf bifurcation values are given by (n)
τj
(0)
= τj
(0)
= τ0 +
2jπ , ω0
j = 0, 1, 2, . . .
(0) (0) (0) Concretely, τ0 = τ0 ≈ 5.81966, τ1 ≈ 23.861, τ2 ≈ 41.9024, ... From Fig. 1, we can see the asymptotical stability of positive equilibrium E ∗ when time delay is slightly smaller than the first bifurcation value τ0. Moreover, we can obtain c1 (0) ≈ −1.4328 + 1.53343i. From Theorem 6, the Hopf bifurcation is supercritical, that is, the periodic solutions exist for τ > τ0, and they are orbitally asymptotically stable (see Fig. 2).
Zhuang and Jia SpringerPlus (2016)5:1303
Fig. 1 The equilibrium E ∗ is stable when τ = 2 < τ0
Fig. 2 Spatially periodic solution exists when τ = 10
Fig. 3 The spatially periodic solution still exists when τ = 50
Fig. 4 The spatially periodic solution still exists even when τ = 130
Page 13 of 15
Zhuang and Jia SpringerPlus (2016)5:1303
In the light of these simulations, we can find that spatially periodic solutions still exist (0) (0) (0) (0) even when τ = 50 ∈ (τ2 , τ3 ) and τ = 130 ∈ (τ6 , τ7 ) (see Figs. 3, 4).
Discussions and conclusions In this paper, we considered the reaction–diffusion regulated logistic growth model. We have investigated the basic properties and Hopf bifurcation under the Neumann boundary conditions. It is shown that the logistic model may undergo Hopf bifurcation when time delay varies. We further give the formulae for determining the bifurcation properties, such as the direction of bifurcation, the stability of periodic solution and the monotonicity of period of periodic solution. Here, we only discussed the single–species diffusive model with feedback control. In fact, how spatial diffusion and time delay affect the dynamic behaviors of multi–species controlled model remains unclear. We will focus on these novel and interesting models in the future. Furthermore, from the numerical simulations in section “Numerical simulations”, we conjecture that the Hopf bifurcation induced by time delay is global. This means that the periodic solutions due to Hopf bifurcation still exist even if the time delay is sufficiently large. Authors’ contributions KZ carried out the genetic studies and drafted the manusctipt. GJ designed the structure of this paper and helped to draft the manuscript. Both authors read and approved the final manuscript. Author details 1 Business School, University of Shanghai for Science and Technology, Shanghai 200093, China. 2 School of Statistics and Applied Mathematics, Anhui University of Finance and Economics, Bengbu 233030, China. 3 College of Science, University of Shanghai for Science and Technology, Shanghai 200093, China. Acknowledgements This work is supported by the National Natural Science Foundation of China (11171220) and the Hujiang Foundation of China (B14005). This work is also sponsored by Key Project for Excellent Young Talents Fund Program of Higher Education Institutions of Anhui Province (gxyqZD2016100). Competing interests The authors declare that they have no competing interests. Received: 10 May 2016 Accepted: 29 July 2016
References Aizerman MA, Gantmacher FR (1964) Absolute stability of refgulator systems. Holden Day, San Francisco Berezansky L, Braverman E, Idels L (2004) Delay differential logistic equations with harvesting. Math Comput Model 39:1243–1259 Chen F, Xie X, Chen X (2006) Permanence and global attractivity of a delayed periodic logistic equation. Appl Math Comput 177:118–127 Chen S, Shi J (2012) Stability and Hopf bifurcation in a diffusive logistic population model with nonlocal delay effect. J Diff Eq 253:3440–3470 Fan Y, Wang L (2010) Global asymptotical stability of a Logistic model with feedback control. Nonlinear Anal Real World Appl 11:2686–2697 Ghergu M, Radulescu VD (2012) Nonlinear PDEs: mathematical models in biology. Chemistry and population genetics. Springer, Heidelberg Gopalsamy K (1993) Peixuan weng, feedback regulation of logistic growth. Intern J Math Math Sci 16:177–192 Hassard B, Kazarinoff N, Yan Y (1981) Theory and applications of Hopf bifurcation. Cambridge University Press, Cambridge Hattaf K, Yousfi N (2015) A generalized HBV model with diffusion and two delays. Comput Math Appl 69:31–40 Hattaf K, Yousfi N (2015) Global dynamics of a delay reaction-diffusion model for viral infection with specific functional response. Comput Appl Math 34:807–818 Henry D (1993) Geometric theory of semilinear parabolic equations. Springer, Berlin Gu H, Lin Z (2014) Long time behavior of solutions of a diffusion-advection logistic model with free boundaries. Appl Math Lett 37:49–53
Page 14 of 15
Zhuang and Jia SpringerPlus (2016)5:1303
Lefschetz S (1965) Stability of nonlinear control systems. Academic Press, New York Murray JD (2003) Mathematical biology II: spatial models and biomedical applications. Springer, New York Röst G (2011) On an approximate method for the delay logistic equation. Commun Nonlinear Sci Numer Simul 16:3470–3474 Song Y, Peng Y (2006) Stability and bifurcation analysis on a Logistic model with discrete and distributed delays. Appl Math Comput 181:1745–1757 Song Y, Yuan S (2007) Bifurcation analysis for a regulated logistic growth model. Appl Math Model 31:1729–1738 Sun C, Han M, Lin Y (2007) Analysis of stability and Hopf bifurcation for a delayed logistic equation. Chaos Solitons Fractals 31:672–682 Sun S, Chen L (2007) Existence of positive periodic solution of an impulsive delay Logistic model. Appl Math Comput 184:617–623 Wu J (1996) Theory and applications of partial functional differrential equations. Springer, New York Yang R (2015) Hopf bifurcation analysis of a delayed diffusive predator-prey system with nonconstant death rate. Chaos Solitons and Fractals 81:224–232 Yang X, Yuan R (2008) Global attractivity and positive almost periodic solution for delay logistic differential equation. Nonlinear Anal 68:54–72 Yang R, Zhang C (2016a) The effect of prey refuge and time delay on a diffusive predator-prey system with hyperbolic mortality. Complexity 50:105–113 Yang R, Zhang C (2016b) Dynamics in a diffusive predator-prey system with a constant prey refuge and delay. Nonlinear Anal Real World Appl 31:1–22 Zhang J, Sun Y (2014) Dynamical analysis of a logistic equation with spatio-temporal delay. Appl Math Comput 247:996–1002
Page 15 of 15
Open Access
RESEARCH
Delay‑induced periodic phenomenon in a diffusive regulated logistic model Kejun Zhuang1,2* and Gao Jia3 *Correspondence: [email protected] 1 Business School, University of Shanghai for Science and Technology, Shanghai 200093, China Full list of author information is available at the end of the article
Abstract The diffusive logistic growth model with time delay and feedback control is considered. First, the well-posedness and permanence of solutions are discussed by using some comparison techniques. Then, the sufficient conditions for stability of nonnegative constant steady states are established, and the occurrence of Hopf bifurcation at positive steady state is performed. Next, the bifurcation properties are derived by computing the normal form on center manifold. Our results not only supplement but also generalized some existing ones. Finally, some numerical simulations show the feasibility of our theoretical analyses. Keywords: Logistic model, Feedback control, Hopf bifurcation, Reaction–diffusion system, Delay Mathematics Subject Classification: 35K57, 35B10, 35B32
Background The classic logistic model N (t) dN (t) = rN (t) 1 − , dt K
r, k ∈ (0, +∞)
(1)
was first proposed by Verhulst in 1838. It can be utilized to describe the single–species growth and has been the basis of varieties of models in population ecology and epidemiology. For system (1) and its generalized forms, the significant results involve the asymptotic properties (Berezansky et al. 2004; Röst 2011), permanence and stability (Fan and Wang 2010; Chen et al. 2006), periodicity (Sun and Chen 2007) and almost periodicity (Yang and Yuan 2008) of solutions, Hopf bifurcation (Sun et al. 2007; Song and Yuan 2007; Song and Peng 2006; Chen and Shi 2012), traveling wave front (Zhang and Sun 2014), free boundary problem (Gu and Lin 2014), and so on. In addition, the Hopf bifurcation analyses for some diffusive predator–prey systems were also done (see Yang 2015; Yang and Zhang 2016a, b). In particular, Gopalsamy (1993) considered the controlled delay system in the following form
dN (t) a1 N (t)+a2 N (t−τ ) K dt = rN (t) 1 − du(t) = bN (t − τ ) − au(t), dt
− cu(t) ,
(2)
© 2016 The Author(s). This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.
Zhuang and Jia SpringerPlus (2016)5:1303
Page 2 of 15
where all the coefficients and time delay τ are positive constants, N(t) is the number of individuals at time t, and variable u(t) denotes an indirect control variable (see Aizerman and Gantmacher 1964; Lefschetz 1965). They have derived the sufficient conditions to guarantee that the positive equilibrium solution is globally asymptotical stable. Strictly speaking, spatial diffusion can not be ignored in studying the natural biological system (Murray 2003; Ghergu and Radulescu 2012). In the real world, most populations are moving and the densities are dependent of time and space. Therefore, diffusion should be taken into account in studying the basic logistic equation. However, there have been very few results on the influence of time delay on the reaction–diffusion logistic model with feedback control. Inspired by the previous discussions, we mainly consider the reaction–diffusion system as follows:
∂N (x,t) = d1 ∆N (x, t) + rN (x, t) 1 − a1 N (x,t)+aK2 N (x,t−τ ) ∂t ∂u(x,t) = d2 ∆u(x, t) + bN (x, t − τ ) − au(x, t), ∂t
− cu(x, t) ,
(3)
where (x, t) ∈ Ω × [0, +∞), Ω = (0, lπ ). The model (3) is considered with the initial value conditions as follows
N (x, t) = η1 (x, t) ≥ 0,
u(x, t) = η2 (x, t) ≥ 0,
x ∈ [0, lπ] × [−τ , 0].
(4)
We also assume that the model (3) is closed and there is no emigration or immigration across the boundary. Hence, the boundary conditions are considered as
∂u ∂N = = 0, ∂ν ∂ν
(x, t) ∈ ∂Ω × [0, +∞),
(5)
where ∂/∂ν represents the outward normal derivative on the boundary ∂Ω. In this paper, we develop a reaction–diffusion logistic model with time delay and diffusion, which makes up perfectly for the deficiencies of the previous literatures. The main objective is to explore the dynamics of system (3) by regarding τ as the bifurcation parameter. The structure of this paper is arranged as follows. In section “Preliminaries”, we derive the well–posedness of solutions and the permanence of the system. In section “Occurrence of the Hopf bifurcation”, we establish the existence of Hopf bifurcation. In section “Bifurcation properties”, we get the formulae for determining the Hopf bifurcation properties. In section “Numerical simulations”, we illustrate our theoretical results by some numerical simulations. Finally, we give some discussions and conclusions.
Preliminaries As we know, spatial diffusion and time delay do not change the number and locations of constant equilibria because of no-flux boundary conditions. Then system (3) has two nonnegative equlibria E0 = (0, 0) and E ∗ = (N ∗ , u∗ ), where N∗ =
aK , a(a1 + a2 ) + bcK
u∗ =
b ∗ bK . N = a a(a1 + a2 ) + bcK
Zhuang and Jia SpringerPlus (2016)5:1303
Page 3 of 15
Well–posedness of solutions
Here, for problem (3)–(5), we devote ourselves to the existence, uniqueness, nonnegativity and boundedness of solutions. Theorem 1 For any given initial data satisfying the conditions (4) and boundary conditions (5), system (3) has a unique global solution of system and the solution maintains nonnegative and uniformly bounded for all t ≥ 0. Proof Using the similar methods in Hattaf and Yousfi (2015), Hattaf and Yousfi (2015), we can get the local existence and uniqueness of solution (N(x, t), u(x, t)) with x ∈ Ω¯ and t ∈ [0, T ), where T is the maximal existence time of solution. It is easy to find that 0 = (0, 0) and M = (M1 , M2 ) are a pair of coupled lower–upper solutions to problem (3)–(5), where
K M1 = max , sup �ϕ1 (·, s)�C(Ω, ¯ R) , a1 −τ ≤s≤0 bM1 , sup �ϕ2 (·, s)�C(Ω, M2 = max ¯ R) . a −τ ≤s≤0 By means of the comparison theorem, we can obtain that 0 ≤ N (x, t) ≤ M1 and 0 ≤ u(x, t) ≤ M2 for x ∈ Ω¯ and t ∈ [0, T ). It is obvious that the upper bound of solution is independent of the maximal existence interval [0, T). It follows from the standard theory for semilinear parabolic systems (Wu 1996; Henry 1993) that the solution glob ally exists. The proof is complete Dissipativeness and permanence
In the following, we will show that system (3) is permanent, which means that any nonnegative solution of (3) is bounded as t → +∞ for all x ∈ Ω. Theorem 2 (Dissipativeness) The nonnegative solution (N, u) of system (3) satisfies
lim sup N (x, t) ≤ t→+∞
K bK , lim sup u(x, t) ≤ . a1 t→+∞ aa1
Proof Based on the first equation in system (3), we get
a1 ∂N (x, t) − d1 ∆N (x, t) ≤ rN (x, t) 1 − N (x, t) ∂t K
for
(x, t) ∈ Ω × [0, +∞).
Then from the standard comparison principle of parabolic equations, we can easily get
lim sup N (x, t) ≤ t→+∞
K . a1
For an arbitrary ε1 > 0, we could get a positive constant T1 such that for any t ≥ T1,
N (x, t) ≤
K + ε1 . a1
Zhuang and Jia SpringerPlus (2016)5:1303
Page 4 of 15
Thus, for any T ∈ [T1 + τ , +∞), we have K ∂u(x, t) − d2 ∆u(x, t) ≤ b + ε1 − au(x, t). ∂t a1 This implies
lim sup u(x, t) ≤ t→+∞
bK aa1
by comparison principle of parabolic equations and the arbitrariness of ε1.
Theorem 3 If aa1 > aa2 + bcK , then system (3) is permanent. Proof From Theorem 2, for an arbitrary ε2 > 0, we can find a constant T > T1 + T2, such that
u(x, t) ≤
bK + ε2 aa1
in Ω × [T2 , +∞). Moreover, we can obtain
∂N (x, t) − d1 ∆N (x, t) ∂t bK a2 K a1 + ε1 − c + ε2 − N (x, t) , ≥ rN (x, t) 1 − K a1 aa1 K the comparison principle shows that
lim inf N (x, t) ≥ t→+∞
K aa1 − aa2 − bcK >0 a1 aa1
due to the continuity as ε1 → 0 and ε2 → 0. Similarly, we can also have
lim inf u(x, t) ≥ t→+∞
bK aa1 − aa2 − bcK > 0. aa1 aa1
Combining the results in Theorem 2, we can easily conclude that system (3) is permanent.
Occurrence of the Hopf bifurcation For system (3), we shall study the local stability of two constant steady states and the occurrence of Hopf bifurcation phenomenon through discussing the distribution of characteristic values. Denote u1 (t) = N (x, t),
u2 (t) = u(x, t),
U (t) = (u1 (t), u2 (t))T .
By defining the phase space C = C([−τ , 0], X), we can rewritten system (3) as the semilinear functional differential equation:
Zhuang and Jia SpringerPlus (2016)5:1303
Page 5 of 15
U˙ (t) = D∆U (t) + G(Ut ),
(6)
where X = {(u, v) ∈ H 2 (0, lπ ) × H 2 (0, lπ )|ux = vx = 0, x = 0, lπ}, Ut (·) = U (t + ·), D = diag {d1 , d2 }, ∆ = diag {∂ 2 /∂x2 , ∂ 2 /∂x2 }, and G(Ut ) : C → X is defined by
G(Ut ) =
ru1 (t) 1 − a1 u1 (t)+aK2 u1 (t−τ ) − cu2 (t) . bu1 (t − τ ) − au2 (t)
The linear system of (6) at E0 (0, 0) is
U˙ (t) = D∆U (t) + LE0 (Ut ),
(7)
where
LE0 (ϕ) =
rϕ1 (0) bϕ1 (−τ )
0 − aϕ2 (0)
for ϕ(θ) = Ut (θ), ϕ = (ϕ1 , ϕ2 )T ∈ C. The characteristic equation of (7) is
y − D∆y − LE0 (e ·y ) = 0,
(8)
where y ∈ dom (∆)\{0}, dom ∆ ⊂ X and e · (θ)y = e θ y for θ ∈ [−τ , 0]. We know that the operator ∆ in Ω with homogeneous Neumann boundary condition has the eigenvalues −n2 /l 2 and the corresponding eigenfunctions cos(nx/l), n ∈ N0 = {0, 1, 2, . . .}. By using the Fourier expansion in (8),
∞ αn y= cos(nx/l), γn n=0
where αn, γn ∈ C. Therefore, the characteristic equation (8) can be transferred into + d n2 − r 0 1 l2 = 0, n ∈ N0 . 2 −be− τ + d2 nl 2 + a
We then obtain the characteristic values as follows
1,n = −d1
n2 + r, l2
2,n = −d2
n2 − a, l2
n ∈ N0 .
It is obvious that 1,0 = r > 0, and we can establish the instability of E0. Theorem 4 The trivial equilibrium E0 of system (3) is always unstable. Next, we will focus on the occurrence of Hopf bifurcation phenomenon. Linearizing system (3) at E ∗ = (N ∗ , u∗ ), we get
U˙ (t) = D∆U (t) + L(Ut ), where L : C → X is given by ra ∗ − K1 N ϕ1 (0) − raK2 N ∗ ϕ1 (−τ ) L(ϕ) = bϕ1 (−τ )
(9)
− cN ∗ ϕ2 (0) − aϕ2 (0)
Zhuang and Jia SpringerPlus (2016)5:1303
Page 6 of 15
with ϕ(θ) = Ut (θ), ϕ = (ϕ1 , ϕ2 )T ∈ C. Similar to the previous discussion, we can obtain the characteristic equation
2 + An + Bn + e− τ (C + Dn ) = 0,
n ∈ N0 ,
(10)
where
ra1 ∗ n2 N > 0, +a+ 2 l K n4 ra1 ∗ n2 raa1 ∗ Bn = d1 d2 4 + a + N d1 2 + N > 0, l K l K ra2 ∗ N > 0, C= K ra2 ∗ n2 raa2 ∗ Dn = N d1 2 + N + bcN ∗ > 0. K l K An = (d1 + d2 )
For τ = 0, Eq. (10) can be reduced to
2 + (An + C) + Bn + Dn = 0 with An + C > 0 and Bn + D > 0. On the basis of Routh–Hurwitz stability criterion, we can obtain the local stability of E ∗ when τ = 0. Lemma 1 The positive equilibrium is always locally asymptotically stable without time delay. Remark 1 From Lemma 1, we can find that there is no Turing instability without time delay. For τ � = 0, let us suppose that = iω(ω > 0) satisfies Eq. (10). First, plugging = iω into Eq. (10) and then segregating the real and imaginary components with the help of Euler’s formula, we can get the following two equations of ω ω2 − Bn = Dn cos ωτ + Cω sin ωτ , −ωAn = Cω cos ωτ − Dn sin ωτ . Second, solving these equations, we can obtain 2 −B D n n cos ωτ = (Dn −A2n C)ω , C ω2 +D2 n
sin ωτ =
Cω3 +(An Dn −Bn C)ω . C 2 ω2 +Dn2
(11)
Third, squaring both sides of those two equations and then adding them up, we get the following equation
ω4 + A2n − 2Bn − C 2 ω2 + Bn2 − Dn2 = 0,
(12)
Zhuang and Jia SpringerPlus (2016)5:1303
Page 7 of 15
where
ra n2 n4 r 2 a21 − a22 ∗2 1 ∗ N d1 + ad2 2 + N , A2n − 2Bn − C 2 = d12 + d22 4 + 2 l K l K2 n4 ra1 ∗ n2 raa1 ∗ ra2 ∗ n2 N d1 2 + N − N d1 2 Bn2 − Dn2 = (Bn + Dn ) d1 d2 4 + a + l K l K K l raa2 ∗ ∗ − N − bcN . K Lemma 2 For τ > 0, we have (i) If a1 > a2 + bcK ar , then Eq. (10) does not have purely imaginary root. (ii) If a2 < a1 < a2 + bcK ar , then there exists N0 ∈ N0, such that Eq. (10) does not have purely imaginary root when n > N0, and has a pair of conjugate purely imaginary eigenvalues when 0 ≤ n ≤ N0. Proof We can easily verify that A2n − 2Bn − C 2 > 0 and Bn2 − Dn2 > 0 when a1 > a2 + bcK ar . This means that Eq. (12) has no positive root. In other words, there could be no purely imaginary root in Eq. (10) for any τ > 0. 2 2 On the contrary, if a2 < a1 < a2 + bcK ar , then B0 − D0 < 0 and there exists N ∈ N0 such that 2 Bn − Dn2 < 0, n = 0, 1, 2, . . . , N0 , Bn2 − Dn2 ≥ 0, n = N0 + 1, N0 + 2, . . . .
That is to say, Eq. (12) has no positive root when n > N0 and has the unique positive root ωn when 0 ≤ n ≤ N0, where
ωn =
−(A2n − 2Bn − C 2 ) +
1 (A2n − 2Bn − C 2 )2 − 4(Bn2 − Dn2 ) 2 . 2
By direct computation, we have
2 ra2 ∗ 2 n4 n 2ra2 ∗ ∗ ∗ N d2 4 + bcN d1 + N d2 + bcN d2 2 An Dn − Bn C = K l K l +
ra1 bc ∗2 ra2 a2 ∗ N + abcN ∗ + N K K
> 0. Moreover, Eq. (10) has characteristic values ±iωn with (n)
τj
(n)
= τ0 +
2jπ , ωn
0 ≤ n ≤ N0 ,
j = 0, 1, 2, . . . ,
where (n)
τ0
=
(Dn − An C)ωn2 − Bn Dn 1 arccos . ωn Dn2 + C 2 ωn2
Zhuang and Jia SpringerPlus (2016)5:1303
Page 8 of 15
This completes the proof. We now check the transversality condition. Lemma 3 If a2 < a1 < a2 + n ∈ {0, 1, 2, . . . , N0 } .
bcK ar ,
then
dRe ( ) dτ τ =τ (n)
for
>0
j
j ∈ N0
and
Proof By taking the derivatives on both sides of (10) with respect to τ, we can get d d d d + An + Ce− τ + (C + Dn )e− τ − − τ = 0, 2 dτ dτ dτ dτ and
d dτ
−1
=
2 + An + Ce− τ − τ e− τ (C + Dn ) e− τ (C + Dn )
=
τ (2 + An )e τ + C − . (C + Dn )
On the basis of (11) and (12), we get
d dτ
−1
(n)
(n) (n) +C (2iωn + An ) cos ωn τj + i sin ωn τj
=
iωn (iCωn + Dn )
τ =τj
(n)
C + An cos ωn τj
=
(n)
− 2ωn sin ωn τj
−Cωn2 + iDn ωn +
(n)
i 2ωn cos ωn τj
(n)
+ An sin ωn τj
−Cωn2 + iDn ωn
(n)
−
τj
iωn
(n)
τj − iωn .
Further simplification will lead to
Re
d dτ
−1
(n)
τ =τj
(n) (n) Dn ωn 2ωn cos ωn τj + An sin ωn τj = 2 Cωn2 + (Dn ωn )2 (n) (n) Cωn2 C + An cos ωn τj − 2ωn sin ωn τj − 2 Cωn2 + (Dn ωn )2
=
ωn4 + Dn2 − Bn2 2 Cωn2 + (Dn ωn )2
> 0.
The proof is complete. According to Lemmas 1–3 and the Hopf bifurcation theory developed by Wu (1996), the following conclusions can be drawn. Theorem 5 Define
τ0 =
min
n∈{0,1,2,...,N0 },j∈N0
(n)
τj
.
Zhuang and Jia SpringerPlus (2016)5:1303
Page 9 of 15
(i) If a1 > a2 + bcK , then for any τ > 0, the positive equilibrium E ∗ is always locally ar asymptotically stable. ∗ (ii) If a2 < a1 < a2 + bcK ar , then E is locally asymptotically stable when τ ∈ [0, τ0 ), and is unstable when τ ∈ (τ0 , +∞). (n) (iii) System (3) has a Hopf bifurcation from E ∗ at τj with n ∈ {0, 1, 2, . . . , N0 } and j ∈ N0. If n = 0, the periodic solutions bifurcating positive equilibrium are all spatially homogeneous. Otherwise, these bifurcating periodic solutions are spatially inhomogeneous.
Bifurcation properties In Theorem 5, we have demonstrated that there exist some spatially homogeneous or inhomogeneous periodic solutions when time delay crosses through some particular values. We are now in the position to investigate the bifurcation properties. (n) with j ∈ N0 and In general, we use τ ∗ to denote an arbitrary value of τj ∗ n ∈ {0, 1, 2, . . . , N0 }. And we also use ±iω to denote the corresponding simply purely imaginary roots ±iωn. ˜ t)), and τ = τ ∗ + α with ˜ t) = u(·, τ t), U˜ (t) = (N˜ (·, t), u(·, Set N˜ (·, t) = N (·, τ t), u(·, α ∈ R. For simplicity we drop the tilde and rewrite system (3) as follows, dU (t) = τ D∆U (t) + L(α)(Ut ) + f (Ut , α), dt
(13)
where ϕ = (ϕ1 , ϕ2 )T ∈ C, L(α)(·) : C → X and f : C × R → X are respectively denoted by ra ∗ − K1 N ϕ1 (0) + raK2 N ∗ ϕ1 (−1) − cN ∗ ϕ2 (0) L(α)(ϕ) = (τ ∗ + α) bϕ1 (−1) − aϕ2 (0) and
2ra − K 1 ϕ12 (0) − f (ϕ, α) = (τ ∗ + α)
ra2 K ϕ1 (0)ϕ1 (−1) − rcϕ1 (0)ϕ2 (0)
0
.
Note that α = τ − τ ∗, we can find that system (13) may causes a Hopf bifurcation when α = 0. For the following linear differential equation:
U˙ (t) = τ D∆U (t) + L(α)(Ut ),
(14)
we can easily deduce that the corresponding characteristic equation has characteristic values ± i ω∗ τ ∗ when α = 0. Next, we discuss the following differential equation:
Y˙ (t) = −τ Dn2 Y (t) + L(α)(Yt ).
(15)
We can use Riesz representation theorem here, which tells us that there is a 2 × 2 matrix function η(θ, α) (−1 ≤ θ ≤ 0) with bounded variation elements satisfying 0 n2 −τ D 2 ϕ(0) + L(α)(ϕ) = d[η(θ, α)]ϕ(θ), l −1
Zhuang and Jia SpringerPlus (2016)5:1303
Page 10 of 15
where
� � 2 −d1 nl 2 − raK1 N ∗ − cN ∗ ∗ , (τ + α) 2 0 − d2 nl 2 − a η(θ, α) = 0, � ra ∗ � − K2 N 0 ∗ , (τ + α) b 0
θ = 0, θ ∈ (−1, 0), θ = −1.
For Φ ∈ C 1 ([−1, 0], R2 ), Ψ ∈ C 1 ([0, 1], R2 ), we define ) dΦ(θ θ ∈ [−1, 0), dθ , A1 (Φ(θ)) = � 0 [dη0 (θ)]Φ(θ), θ = 0, −1 dΨ − ds(s) , s ∈ (0, 1], ∗ A1 (Ψ (s)) = � 0 s = 0. −1 [dη0 (θ)]Ψ (−θ ),
Then the formal adjoint, A∗1, of A1 is given by
(Ψ , Φ)0 = Ψ (0)Φ(0) −
0
−1
= Ψ (0)Φ(0) + τ
∗
θ
Ψ (ζ − θ )d[η(θ, 0)]Φ(ζ )dζ
ζ =0 0
ra ∗ − K2 N 0 Ψ (ζ + 1) Φ(ζ )dζ . b 0 −1
∗
∗
By calculation, we can find that q(θ) = (1, ξ )T eiω θ τ and q ∗ (s) = M(1, η)eiω eigenvectors of A1 and A∗1 associated with iω∗ τ ∗, respectively, where
θ ∈ [−1, 0],
∗ sτ ∗
are
s ∈ [0, 1],
and ∗ ∗
be−iω τ ra2 ∗ iω∗ τ ∗ −1 cN ∗ ∗ ξ= ξ η + τ bη − N e . , η = , M = 1 + iω∗ + a iω∗ + a K Then P = span{q(θ), q(θ)}, P ∗ = span{q ∗ (s), q ∗ (s)} are the center subspace of system (3). T 2 nx T . Define h · fn = h1 βn1 + h2 βn2, fn = βn1 , βn2 and βn1 = cos nx l , 0 , βn = 0, cos l The complex-valued L2 inner product on Hilbert space XC are lπ 1 �U1 , U2 � = (u1 v1 + u2 v2 )dx, (16) lπ 0 for U1 = (u1 , u2 ), U2 = (v1 , v2 ) ∈ XC . And �β0i , β0i � = 1, �βni , βni � = 12, i = 1, 2, n = 1, 2, . . .,
�Φ, fn � = �Φ, βn1 �, �Φ, βn2 � ,
(17)
where Φ ∈ C([−1, 0], X). We can establish the center subspace of system (14) at α = 0 as follows PCN L = (q(θ)z + q(θ)¯z ) · fn , z ∈ C . (18)
Zhuang and Jia SpringerPlus (2016)5:1303
Page 11 of 15
Based on the conclusions drawn by Wu (1996) and Hassard et al. (1981), the solutions of (13) are
Ut = (q(θ)z(t) + q(θ)¯z (t)) · fn + W (z(t), z¯ (t), θ ), where
W (z, z¯ , θ ) = W20
z¯ 2 z2 + W11 z z¯ + W02 + · · · . 2 2
(19)
Moreover, for Ut ∈ C0 of (13) at τ = τ ∗, we have z˙ = iω∗ τ ∗ z + g(z, z¯ ), where
g(z, z¯ ) = q ∗ (0)�f (Ut , 0), fn � = g20
z2 z¯ 2 z 2 z¯ + g11 z z¯ + g02 + g21 + ··· . 2 2 2
(20)
By (16)–(20), we can compute
g20 =
g11 =
0, −2τ ∗ M raK1 +
+ rcξ ,
0, −2τ ∗ M raK1 + rc Re {ξ } +
g02 = g20 , g21
ra2 −iω∗ τ ∗ K e
ra2 K
Re
n = 1, 2, . . . , n = 0,
∗ ∗ {eiω τ }
,
n = 1, 2, . . . , n = 0,
lπ 2ra Mτ ∗ nx 1 (1) (1) (W11 dx = −2 (0) + W20 (0)) cos2 lπ K l 0 lπ nx 1 1 (2) (1) (1) (2) ¯ dx rc W11 (0) + W20 (0) + ξ W20 (0) + ξ W11 (0) cos2 + 2 2 l 0 lπ nx ra2 −iω∗ τ ∗ (2) 1 ∗ ∗ (2) 1 (1) + dx e W11 (0) + eiω τ W20 (0) + ξ¯ W20 (−1) cos2 K 2 2 l 0 lπ ra2 (11) 2 nx ξ W11 (−1) cos dx . + K l 0
Then we should compute W20 (θ) and W11 (θ) to determine g21. Following the formulas in Wu (1996), We can obtain that ig20 ig02 ∗ ∗ W20 (θ) = q(θ) · fn + E1 e2iω τ θ , q(θ) + ω∗ τ ∗ 3ω∗ τ ∗ ig11 q(θ) ig11 · fn + E2 , W11 (θ) = − ∗ ∗ 1(θ) + ω τ ω∗ τ ∗ 2ra ∗ ∗ 2 −iω τ nx − K 1 − 2ra e − 2rcξ ′ K cos2 , E1 = E1 × 0 l −1 2 ∗ ∗ 2iω∗ + d1 nl 2 + raK1 N ∗ + raK2 N ∗ e−2iω τ cN ∗ ′ E1 = , 2 ∗ ∗ −be−2iω τ 2iω∗ + a + d2 nl 2 and
Zhuang and Jia SpringerPlus (2016)5:1303
Page 12 of 15
1 − 2ra K − 2rc Re {ξ } − 0
E2 =
E2′
×
E2′
d1 nl 2 +
=
2
ra1 ∗ K N
+
ra2 ∗ K N
2ra2 K
Re {eiω
cN ∗ 2
a + d2 nl 2
−b
−1
∗τ ∗
}
cos2
nx , l
.
From the previous expressions of g20, g11, g02 and g21, we can further compute
g21 1 2 g20 g11 − 2|g11 | − |g02 | + , c1 (0) = 2ω∗ τ ∗ 3 2 Re (c1 (0)) , µ2 = − Re ( ′ (τ ∗ )) β2 = 2 Re (c1 (0)), 1 T2 = − ∗ ∗ ( Im (c1 (0)) + µ2 Im ( ′ (τ ∗ ))). ω τ i
2
On account of preceding calculations, we arrive at the following conclusion on the bifurcation properties. Theorem 6 The bifurcation direction is supercritical if µ2 > 0, which means that the periodic solution exists for τ > τ0. On the contrary, the bifurcation direction is subcritical if µ2 < 0, which means that the periodic solution exists for τ < τ0. Moreover, the periodic solution is orbitally asymptotically stable if β2 < 0, or unstable if β2 > 0. The period of periodic solution is monotonically increasing at the time delay τ when T2 > 0, or is monotonically decreasing at the time delay τ when T2 < 0.
Numerical simulations In this section, we give some numerical examples to test the preceding results with assistance of MATLAB. For system (3), let Ω = (0, 2π ) and choose d1 = 1,
d2 = 0.5,
r = 0.6,
a = b = c = 1,
a1 = a2 = 2,
K = 1,
and the initial values N (x, 0) = 0.5 and u(x, 0) = 0.9. Then we can get the positive equilibrium E ∗ = (0.2, 0.2). By direct computation, we have N0 = 0, ω0 ≈ 0.348266, and τ0(0) ≈ 5.81966, then the Hopf bifurcation values are given by (n)
τj
(0)
= τj
(0)
= τ0 +
2jπ , ω0
j = 0, 1, 2, . . .
(0) (0) (0) Concretely, τ0 = τ0 ≈ 5.81966, τ1 ≈ 23.861, τ2 ≈ 41.9024, ... From Fig. 1, we can see the asymptotical stability of positive equilibrium E ∗ when time delay is slightly smaller than the first bifurcation value τ0. Moreover, we can obtain c1 (0) ≈ −1.4328 + 1.53343i. From Theorem 6, the Hopf bifurcation is supercritical, that is, the periodic solutions exist for τ > τ0, and they are orbitally asymptotically stable (see Fig. 2).
Zhuang and Jia SpringerPlus (2016)5:1303
Fig. 1 The equilibrium E ∗ is stable when τ = 2 < τ0
Fig. 2 Spatially periodic solution exists when τ = 10
Fig. 3 The spatially periodic solution still exists when τ = 50
Fig. 4 The spatially periodic solution still exists even when τ = 130
Page 13 of 15
Zhuang and Jia SpringerPlus (2016)5:1303
In the light of these simulations, we can find that spatially periodic solutions still exist (0) (0) (0) (0) even when τ = 50 ∈ (τ2 , τ3 ) and τ = 130 ∈ (τ6 , τ7 ) (see Figs. 3, 4).
Discussions and conclusions In this paper, we considered the reaction–diffusion regulated logistic growth model. We have investigated the basic properties and Hopf bifurcation under the Neumann boundary conditions. It is shown that the logistic model may undergo Hopf bifurcation when time delay varies. We further give the formulae for determining the bifurcation properties, such as the direction of bifurcation, the stability of periodic solution and the monotonicity of period of periodic solution. Here, we only discussed the single–species diffusive model with feedback control. In fact, how spatial diffusion and time delay affect the dynamic behaviors of multi–species controlled model remains unclear. We will focus on these novel and interesting models in the future. Furthermore, from the numerical simulations in section “Numerical simulations”, we conjecture that the Hopf bifurcation induced by time delay is global. This means that the periodic solutions due to Hopf bifurcation still exist even if the time delay is sufficiently large. Authors’ contributions KZ carried out the genetic studies and drafted the manusctipt. GJ designed the structure of this paper and helped to draft the manuscript. Both authors read and approved the final manuscript. Author details 1 Business School, University of Shanghai for Science and Technology, Shanghai 200093, China. 2 School of Statistics and Applied Mathematics, Anhui University of Finance and Economics, Bengbu 233030, China. 3 College of Science, University of Shanghai for Science and Technology, Shanghai 200093, China. Acknowledgements This work is supported by the National Natural Science Foundation of China (11171220) and the Hujiang Foundation of China (B14005). This work is also sponsored by Key Project for Excellent Young Talents Fund Program of Higher Education Institutions of Anhui Province (gxyqZD2016100). Competing interests The authors declare that they have no competing interests. Received: 10 May 2016 Accepted: 29 July 2016
References Aizerman MA, Gantmacher FR (1964) Absolute stability of refgulator systems. Holden Day, San Francisco Berezansky L, Braverman E, Idels L (2004) Delay differential logistic equations with harvesting. Math Comput Model 39:1243–1259 Chen F, Xie X, Chen X (2006) Permanence and global attractivity of a delayed periodic logistic equation. Appl Math Comput 177:118–127 Chen S, Shi J (2012) Stability and Hopf bifurcation in a diffusive logistic population model with nonlocal delay effect. J Diff Eq 253:3440–3470 Fan Y, Wang L (2010) Global asymptotical stability of a Logistic model with feedback control. Nonlinear Anal Real World Appl 11:2686–2697 Ghergu M, Radulescu VD (2012) Nonlinear PDEs: mathematical models in biology. Chemistry and population genetics. Springer, Heidelberg Gopalsamy K (1993) Peixuan weng, feedback regulation of logistic growth. Intern J Math Math Sci 16:177–192 Hassard B, Kazarinoff N, Yan Y (1981) Theory and applications of Hopf bifurcation. Cambridge University Press, Cambridge Hattaf K, Yousfi N (2015) A generalized HBV model with diffusion and two delays. Comput Math Appl 69:31–40 Hattaf K, Yousfi N (2015) Global dynamics of a delay reaction-diffusion model for viral infection with specific functional response. Comput Appl Math 34:807–818 Henry D (1993) Geometric theory of semilinear parabolic equations. Springer, Berlin Gu H, Lin Z (2014) Long time behavior of solutions of a diffusion-advection logistic model with free boundaries. Appl Math Lett 37:49–53
Page 14 of 15
Zhuang and Jia SpringerPlus (2016)5:1303
Lefschetz S (1965) Stability of nonlinear control systems. Academic Press, New York Murray JD (2003) Mathematical biology II: spatial models and biomedical applications. Springer, New York Röst G (2011) On an approximate method for the delay logistic equation. Commun Nonlinear Sci Numer Simul 16:3470–3474 Song Y, Peng Y (2006) Stability and bifurcation analysis on a Logistic model with discrete and distributed delays. Appl Math Comput 181:1745–1757 Song Y, Yuan S (2007) Bifurcation analysis for a regulated logistic growth model. Appl Math Model 31:1729–1738 Sun C, Han M, Lin Y (2007) Analysis of stability and Hopf bifurcation for a delayed logistic equation. Chaos Solitons Fractals 31:672–682 Sun S, Chen L (2007) Existence of positive periodic solution of an impulsive delay Logistic model. Appl Math Comput 184:617–623 Wu J (1996) Theory and applications of partial functional differrential equations. Springer, New York Yang R (2015) Hopf bifurcation analysis of a delayed diffusive predator-prey system with nonconstant death rate. Chaos Solitons and Fractals 81:224–232 Yang X, Yuan R (2008) Global attractivity and positive almost periodic solution for delay logistic differential equation. Nonlinear Anal 68:54–72 Yang R, Zhang C (2016a) The effect of prey refuge and time delay on a diffusive predator-prey system with hyperbolic mortality. Complexity 50:105–113 Yang R, Zhang C (2016b) Dynamics in a diffusive predator-prey system with a constant prey refuge and delay. Nonlinear Anal Real World Appl 31:1–22 Zhang J, Sun Y (2014) Dynamical analysis of a logistic equation with spatio-temporal delay. Appl Math Comput 247:996–1002
Page 15 of 15
