# Stability and Hopf Bifurcation Analysis of a Plant Virus Propagation Model with Two Delays

Authors: Junli Liu; Tailei Zhang
Journal: Discrete Dynamics in Nature and Society
 (2018)

---

## Abstract

To understand the interaction between the insects and the plants, a system of delay differential equations is proposed and studied. We prove that ifR0≤1, the disease-free equilibrium is globally asymptotically stable for any length of time delays by constructing a Lyapunov functional, and the system admits a unique endemic equilibrium if R0>1. We establish the sufficient conditions for the stability of the endemic equilibrium and existence of Hopf bifurcation. Using the normal form theory and center manifold theorem, the explicit formulae which determine the stability, direction, and other properties of bifurcating periodic solutions are derived. Some numerical simulations are given to confirm our analytic results.


---

## Body

## 1. Introduction

Plants are very important not only to man’s survival but to every species on Earth; however, plants may contract a disease by many different ways. Tremendous crop losses and global threat to food security have been caused by plant diseases [1, 2]. In recent years, plant diseases have attracted the interest of many mathematical modeling researchers and epidemiologists [3–7].Mathematical models provide powerful tools for investigating how an infection propagates within a population of plants. Shi et al. [8] have proposed an epidemic model which describes vector-borne plant diseases, and the global dynamics of the system have been analyzed in terms of the basic reproduction number. Luo et al. [9] studied a discrete plant virus disease model with roguing and replanting; they proved that the basic reproduction number serves as a threshold parameter in determining the global dynamics of the model. The plant diseases epidemic models have been extensively studied by many authors (see [10–14]).In [15], a delay differential equations was proposed to model the interaction between plants, a plant virus, and the insect-vector that transfers the virus from one plant to another. Since insects can only bite a limited number of plants, then the interaction between vector and plant is of predator-prey Holling type II [16]. In order to consider the time it takes for the virus to spread throughout the plant or insect-vector, a couple of delays were introduced to the model (see [15] for more details). They obtained the following model:(1)dSdt=μK-S-βYt-τ11+αYt-τ1St-τ1+dI,dIdt=βYt-τ11+αYt-τ1St-τ1-d+μ+γI,dRdt=γI-μR,dXdt=Λ-β1It-τ21+α1It-τ2Xt-τ2-mX,dYdt=β1It-τ21+α1It-τ2Xt-τ2-mY,where the state variables S(t), I(t), and R(t) represent the number of susceptible, infected, and recovered plants at time t, respectively. Because when a plant dies by the virus or natural death in farms, it is replaced with a new healthy plant. The new plant shares the same characteristics of the plant it replaced, before it was infected. Then it is supposed that the total number of plants stabilizes at K, K=S+I+R. μ is the natural death rate of plants; β is the infection rate of plants due to vectors; α is the saturation constant of plants due to vectors; d is the death rate of infected plants due to the disease; γ is the recovery rate of plants. The insect-vectors are divided into two populations: susceptible and infective denoted by X(t) and Y(t), respectively. The total number of insects is denoted by N(t), and then N(t)=X(t)+Y(t). Λ is the replenishing rate of vectors (birth and/or immigration); β1 is the infection rate of vectors due to plants; α1 is the saturation constant of vectors due to plants; m is the natural death rate of vectors. τ1 is the time it takes a plant to become infected after contagion, and τ2 is the time it takes a vector to become infected after contagion.Notice that(2)dNdt=Λ-mN,and then N(t)→Λ/m as t→∞.Thus, one can consider the following reduced system:(3)dSdt=μK-S-βYt-τ11+αYt-τ1St-τ1+dI,dIdt=βYt-τ11+αYt-τ1St-τ1-aI,dYdt=β1It-τ21+α1It-τ2Λm-Yt-τ2-mY,where a=d+μ+γ.For model (3), Jackson and Chen-Charpentier [15] gave the basic reproduction number and found the equilibria of the model, and then they studied the stability of equilibria only for particular values of the parameters using numerical methods. Therefore, in this paper, we reconsider the plant disease model (3) in theoretical aspects, and we establish the stability of equilibria, the existence of Hopf bifurcation, and the stability, direction, and other properties of bifurcating periodic solution will also be discussed.This paper is organized as follows. In Section2, we discuss the stability of the equilibria and the existence of the Hopf bifurcations occurring at the endemic equilibrium. In Section 3, the formulae determining the direction of the Hopf bifurcations and the stability of bifurcating periodic solutions on the center manifold are obtained by using the normal form theory and the center manifold theorem by Hassard et al. [17]. In Section 4, we perform numerical simulations to illustrate the analytical results. We conclude with a brief discussion in Section 5.


## 2. Stability Analysis and Hopf Bifurcation

Letτ=max⁡{τ1,τ2}, and the initial conditions for (3) are (4)ϕ1θ,ϕ2θ,ϕ3θ∈C+=C-τ,0,R+3,ϕi0>0,i=1,2,3.By the fundamental theory of functional differential equations [18], it follows that, for any initial conditions (4), there is a unique solution (S(t),I(t),Y(t)) of (3) for all t≥0.LetΩ be the following subset of R+3: (5)Ω=S,I,Y∈R+3:S+I≤K,Y≤Λm.Using a proof process similar to that in [19, 20], we obtain the following lemma.Lemma 1.
The solutions of system (3) which satisfy the initial conditions (4) are positive. The set Ω is positively invariant.In [15], the basic reproduction number for (3) has been identified as (6)R0=ββ1ΛKm2a.Equation (3) always has a disease-free equilibrium E0(K,0,0). If R0>1, then (3) admits a unique endemic equilibrium E∗(S∗,I∗,Y∗), where (7)S∗=aαβ1KΛμ+mma+Kμβ1+α1m-dmβ1βΛμ+γ+aβ1Λμa+β1mμa+α1m2μa,I∗=β1βKΛμ-m2μaβ1βΛμ+γ+aβ1Λμa+β1mμa+α1m2μa,Y∗=β1βKΛμ-m2μamαmμa+βma+Kμβ1+α1m-dm.However, Jackson and Chen-Charpentier [15] did not give detailed dynamical analysis to this model. Theoretical analysis makes the model dynamics clear and enhances our understanding to the mathematical models. In this paper, we will give some analytic results of model (3).Linearizing system (3) at E0 gives characteristic equation (8)λ+μλ2+m+aλ+ma-ββ1KΛme-λτ=0,τ=τ1+τ2.It is clear that -μ is one root of (8). Let G(λ)=λ2+(m+a)λ+ma-(ββ1KΛ/m)e-λτ. If R0>1, we get G(0)=ma(1-R02)<0, and limλ→+∞G(λ)=+∞, then G(λ)=0 has one positive real root, and, hence, E0 is unstable.IfR0<1, it is easy to show that E0 is locally asymptotically stable when τ=0, and then, by Theorem 3.4.1 in Kuang [21], E0 is locally asymptotically stable for all τ≥0.Theorem 2.
IfR0≤1, then E0 is globally asymptotically stable for all τ1≥0 and τ2≥0.Proof.
Constructing the following Lyapunov functional:(9)Vt=It+βKmYt+∫t-τ1tβSρYρ1+αYρdρ+ββ1KΛm2∫t-τ2tIρ1+α1Iρdρ,then (10)V′t2=βSY1+αY-βKY-ββ1KIt-τ2Yt-τ2m1+α1It-τ2+ββ1KΛm2I1+α1I-aI≤βSY-βKY+ββ1KΛm2I1+α1I-aI≤ββ1KΛm2I-aI=R02-1aI≤0.If we set {(S,I,Y)∣V′(t)=0}, then the largest invariant set is the singleton {E0}. Therefore, by LaSalle’s invariance principle [22], E0 is globally asymptotically stable for all τ1≥0 and τ2≥0.
We now consider the local stability of the coexistence equilibriumE∗(S∗,I∗,Y∗) and the existence of Hopf bifurcation at E∗. The linearized system of (3) at E∗ is given by (11)dSdt=a11St+b11St-τ1+a12It+b13Yt-τ1,dIdt=b21St-τ1+b23Yt-τ1+a22It,dYdt=b32It-τ2+a33Yt+b33Yt-τ2,with (12)a11=-μ,b11=-βY∗1+αY∗,a12=d,b13=-βS∗1+αY∗2,b21=βY∗1+αY∗,b23=βS∗1+αY∗2,a22=-a,b32=β11+α1I∗2Λm-Y∗,a33=-m,b33=-β1I∗1+α1I∗.Therefore, we obtain the following characteristic equation: (13)λ3+m2λ2+m1λ+m0+e-λτ1n2λ2+n1λ+n0+e-λτ2p2λ2+p1λ+p0+e-λτ1+τ2q1λ+q0=0,where (14)m2=-a11+a22+a33,m1=a11a22+a11a33+a22a33,m0=-a11a22a33,n0=b21a33a12+a22,n1=-b21a12+a22+a33,n2=b21,p2=-b33,p1=b33a11+a22,p0=-b33a11a22,q1=-b21b33+b23b32,q0=b21b33a12+a22+b23b32a11.Case 1 (τ1=τ2=0).
Characteristic equation (13) becomes (15)λ3+m12λ2+m11λ+m10=0,where (16)m12=m2+n2+p2,m11=m1+n1+p1+q1,m10=m0+n0+p0+q0.
Note that(17)β1I∗1+α1I∗Λm-Y∗=mY∗,βY∗S∗1+αY∗=aI∗,then we get (βS∗/1+αY∗)(β1/(1+α1I∗))(Λ/m-Y∗)=ma, and thus (18)m12=μ+a+m+βY∗1+αY∗+β1I∗1+α1I∗>0,m11=μa+μm+ma+βY∗1+αY∗m+μ+γ+β1I∗1+α1I∗a+μ+ββ1Y∗I∗1+αY∗1+α1I∗-βS∗1+αY∗2β11+α1I∗2Λm-Y∗>0,m10=maμ+m+β1I∗1+α1I∗μ+γβY∗1+αY∗+μaβ1I∗1+α1I∗-βS∗1+αY∗2β11+α1I∗2Λm-Y∗μ>0.Since m11>μm+(βY∗/(1+αY∗))(μ+γ)+β1I∗a/(1+α1I∗) and m12>μ+a+m+β1I∗/(1+α1I∗), it is easy to show that m11m12-m10>0, and, thus, all roots of (15) have negative real parts. That is, E∗(S∗,I∗,Y∗) is locally asymptotically stable. Actually, by a similar proof as in [23], we can show that E∗(S∗,I∗,Y∗) is globally asymptotically stable for τ1=τ2=0.Case 2 (τ1>0,τ2=0).
Characteristic equation (13) becomes (19)λ3+m22λ2+m21λ+m20+n22λ2+n21λ+n20e-λτ1=0,where (20)m22=m2+p2,m21=m1+p1,m20=m0+p0,n22=n2,n21=n1+q1,n20=n0+q0.Suppose λ=iω(ω>0) is a root of (19), similar discussion as those in [24], and we have (21)gρ=ρ3+e22ρ2+e21ρ+e20=0,where ρ=ω2, e22=m222-n222-2m21, e21=m212-n212-2m20m22+2n20n22, and e20=m202-n202.Note thatg(0)=e20 and limρ→+∞g(ρ)=+∞, and then, by [25], we have the following lemma.Lemma 3.
For polynomial equation (21), we have the following results:
(1) If (H21)e20≥0 and Δ=e222-3e21≤0, then (21) has no positive root.
(2) If (H22)e20≥0, Δ=e222-3e21>0, ρ∗=(-e21+Δ)/3>0, g′(ρ∗)=3(ρ∗)2+2e22ρ∗+e21≤0, or (H23) e20<0, then (21) has positive root.Suppose that (21) has positive roots, and we assume that (21) has three positive roots: ρ1,ρ2, and ρ3; then ωk=ρk,k=1,2,3. The corresponding critical value of time delay τ1k(j) is (22)τ1kj=1ωkarccos⁡A24ωk4+A22ωk2+A20B24ωk4+B22ωk2+B20+2πjωk,k=1,2,3;j=0,1,2,…,where (23)A24=n21-m22n22,A22=n20m22+m20-n21m21,A20=-n20m20,B24=n222,B22=n212-2n20n22,B20=n202.±iωk is a pair of purely imaginary roots of (19) with τ1=τ1k(j). Let τ10∗=mink∈{1,2,3}{τ1k(0)}, when τ1=τ10∗, and (19) has a pair of purely imaginary roots ±iω10∗.We now verify the transversality condition, again by the analysis in [24], and we get (24)signdRe⁡λdτ1λ=iω10∗=signRe⁡dλdτ1-1λ=iω10∗=g′ω10∗2n21ω10∗2+n20-n22ω10∗22,assuming that (25)H24g′ω10∗2≠0.Therefore, we have the following result.Theorem 4.
For system (3), if τ2=0,
(1) if (H21) holds, then the endemic equilibriumE∗(S∗,I∗,Y∗) is locally asymptotically stable for all τ1≥0,
(2) if (H22) or (H23) and (H24) hold, then asτ1 increases from zero, there is a value τ10∗ such that the endemic equilibrium E∗(S∗,I∗,Y∗) is locally asymptotically stable when τ1∈[0,τ10∗) and unstable when τ1>τ10∗. Furthermore, system (3) undergoes a Hopf bifurcation at E∗ when τ1=τ10∗.Case 3 (τ1=0,τ2>0).
With similar analysis as to Case2, we get the following theorem.Theorem 5.
For system (3), if τ1=0,
(1) if (H31) holds, then the endemic equilibriumE∗(S∗,I∗,Y∗) is locally asymptotically stable for all τ2≥0,
(2) if (H32) or (H33) and (H34) hold, then asτ2 increases from zero, there is a value τ20∗ such that the endemic equilibrium E∗(S∗,I∗,Y∗) is locally asymptotically stable when τ2∈[0,τ20∗) and unstable when τ2>τ20∗. Furthermore, system (3) undergoes a Hopf bifurcation at E∗ when τ2=τ20∗.Assumptions (H31)–(H34) are very similar to (H21)–(H24), so we omit them.Case 4 (τ1>0,τ2∈[0,τ20∗)).
We consider (3) with τ2 in its stable interval and regard τ1 as a parameter. Let λ=iω1(ω1>0) be a root of (13), separating real and imaginary parts, and we have the following: (26)A41sin⁡ω1τ1+A42cos⁡ω1τ1=A43,A41cos⁡ω1τ1-A42sin⁡ω1τ1=A44,where (27)A41=-n1ω1+q0sin⁡ω1τ2-q1ω1cos⁡ω1τ2,A42=n2ω12-n0-q0cos⁡ω1τ2-q1ω1sin⁡ω1τ2,A43=-ω12m2+p2cos⁡ω1τ2+p1ω1sin⁡ω1τ2+p0cos⁡ω1τ2+m0,A44=-ω13+p2ω12sin⁡ω1τ2+ω1m1+p1cos⁡ω1τ2-p0sin⁡ω1τ2.
From (26), we have (28)ω16+e42ω14+e41ω12+e40+cos⁡ω1τ2c44ω14+c42ω12+c40+sin⁡ω1τ2c45ω15+c43ω13+c41ω1=0,where (29)e40=m02+p02-q02-n02,e41=m12-2m0m2-2p0p2+p12-n12+2n0n2-q12,e42=m22-2m1+p22-n22,c40=2m0p0-n0q0,c41=2m0p1-p0m1+n1q0-n0q1,c42=2m1p1-m0p2-m2p0-n1q1+q0n2,c43=2m1p2-m2p1+p0+n2q1,c44=2m2p2-p1,c45=-2p2.We make the following assumption.
( H 41 ) Equation (28) has finite positive roots ω1(1),ω1(2),…,ω1(k).
For every fixedω1(i),i=1,2,…,k, there exists a sequence τ1i(j) such that (28) holds, where (30)τ1ij=1ω1iarccos⁡A41A44+A42A43A412+A422+2πjω1i,i=1,2,…,k;j=0,1,2,….Let τ10=min⁡{τ1i(0)∣i=1,2,…,k}, when τ1=τ10, and (13) has a pair of purely imaginary roots ±iω10.
In addition to(H41), we further assume that (31)H42ddτ1Re⁡λλ=iω10≠0.Therefore, by the Hopf bifurcation theorem for functional differential equations in Hale [18], the following result holds.Theorem 6.
For system (3), suppose (H41) and (H42) are satisfied, and τ1>0 and τ2∈[0,τ20∗). Then the positive equilibrium E∗(S∗,I∗,Y∗) is asymptotically stable when τ1∈(0,τ10) and unstable when τ1>τ10. Furthermore, system (3) undergoes a Hopf bifurcation at E∗ when τ1=τ10.For the casesτ2>0, and τ1∈[0,τ10∗), we can get similar results as those in Theorem 6.


## 3. Direction and Stability of the Hopf Bifurcation

In this section, we shall study the direction of the Hopf bifurcation and stability of bifurcating periodic solutions by using the normal form theory and the center manifold theorem due to Hassard et al. [17]. In the previous section, we have shown that system (3) undergoes the Hopf bifurcation at τ1=τ10, without loss of generality, and we assume that τ2∗<τ10, where τ2∗∈(0,τ20∗).Letx1=S-S∗,x2=I-I∗,x3=Y-Y∗, τ1=τ10+v, and x¯i(t)=xi(τ1t) and, dropping the bars for simplification of notations, system (3) is transformed into a functional differential equation in C=C([-1,0],R3) as (32)x˙t=Lvxt+Fv,xt,where x(t)=(x1(t),x2(t),x3(t))T∈R3 and Lv:C→R3 and F:R×C→R3 are given, respectively, by (33)Lvϕ=τ10+vA~ϕ0+τ10+vB~ϕ-τ2∗τ1+τ10+vC~ϕ-1,(34)Fv,ϕ=τ10+vF1,F2,F3T,where (35)A~=a11a1200a22000a33,B~=0000000b32b33,C~=b110b13b210b23000,F1=K11ϕ1-1ϕ3-1+K12ϕ32-1+K13ϕ1-1ϕ32-1+K14ϕ33-1+⋯,F2=K21ϕ1-1ϕ3-1+K22ϕ32-1+K23ϕ1-1ϕ32-1+K24ϕ33-1+⋯,F3=K31ϕ22-τ2∗τ1+K32ϕ2-τ2∗τ1ϕ3-τ2∗τ1+K33ϕ23-τ2∗τ1+K34ϕ22-τ2∗τ1ϕ3-τ2∗τ1+⋯,K11=-β1+αY∗2,K12=αβS∗1+αY∗3,K13=αβ1+αY∗3,K14=-α2βS∗1+αY∗4,K21=β1+αY∗2,K22=-αβS∗1+αY∗3,K23=-αβ1+αY∗3,K24=α2βS∗1+αY∗4,K31=α1β1Y∗-Λ/m1+α1I∗3,K32=-β11+α1I∗2,K33=α12β1Λ/m-Y∗1+α1I∗4,K34=α1β11+α1I∗3.By the Riesz representation theorem, there exists a functionη(θ,v) of bounded variation for θ∈[-1,0], such that (36)Lvϕ=∫-10dηθ,vϕθforϕ∈C.In fact, we can choose(37)ηθ,v=τ10+vA~+B~+C~,θ=0,τ10+vB~+C~,θ∈-τ2∗τ1,0,τ10+vC~,θ∈-1,-τ2∗τ1,0,θ=-1.For ϕ∈C1([-1,0],R3), define (38)Avϕ=dϕθdθ,θ∈-1,0,∫-10dηs,vϕs,θ=0,Rvϕ=0,θ∈-1,0,Fv,ϕ,θ=0.Then system (32) is equivalent to (39)x˙t=Avxt+Rvxt,where xt(θ)=x(t+θ) for θ∈[-1,0].Forψ∈C1([0,1],(R3)∗), define (40)A∗ψs=-dψsds,s∈0,1,∫-10dηTt,0ψ-t,s=0,and a bilinear inner product (41)ψ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ω10τ10 are eigenvalues of A(0). Hence, they are also eigenvalues of A∗. We first need to compute the eigenvectors of A(0) and A∗ corresponding to iω10τ10 and -iω10τ10, respectively.Supposeq(θ)=(1,q2,q3)Teiω10τ10θ is the eigenvector of A(0) corresponding to iω10τ10, and then A(0)q(θ)=iω10τ10q(θ). Then, from the definition of A(0) and (33), we have (42)q2=a11-iω10iω10-a22-a12,q3=b32e-iω10τ2∗a11-iω10iω10-a33-b33e-iω10τ2∗iω10-a22-a12.Similarly, we can obtain the eigenvector q∗(s)=D(1,q2∗,q3∗)eiω10τ10s of A∗ corresponding to -iω10τ10, where (43)q2∗=-a11+iω10+b11eiω10τ10b21eiω10τ10,q3∗=b23a11+iω10b21a33+b33eiω10τ2∗+iω10.Choosing D¯ as D¯={1+q2q¯2∗+q3q¯3∗+τ10e-iω10τ10(b11+q¯2∗b21+q3b13+q3q¯2∗b23)+τ2∗q¯3∗e-iω10τ2∗(q2b32+q3b33)}-1, then, by (41), we see q∗(s),q(θ)=1.In the remainder of this section, we use the algorithms given in [17] and, using a computation process similar to that in [24–27], we get the coefficients used in determining the qualities of bifurcating periodic solutions:(44)g20=2τ10D¯q3e-2iω10τ101-q¯2∗K11+q3K12+q¯3∗q2e-2iω10τ2∗q2K31+q3K32,g11=τ10D¯1-q¯2∗K11q3+q¯3+2K12q3q¯3+q¯3∗K32q2q¯3+q¯2q3+2K31q2q¯2,g02=2τ10D¯e2iω10τ10q¯31-q¯2∗K11+q¯3K12+q¯3∗q¯2e2iω10τ2∗q¯2K31+q¯3K32,g21=2τ10D¯1-q¯2∗K11W201-12q¯3+W203-12eiω10τ10+W111-1q3+W113-1e-iω10τ10+K12W203-1q¯3eiω10τ10+2W113-1q3e-iω10τ10+K13q32+2q3q¯3e-iω10τ10+3K14q32q¯3e-iω10τ10+q¯3∗K31W202-τ2∗τ10q¯2eiω10τ2∗+2W112-τ2∗τ10q2e-iω10τ2∗+K32eiω10τ2∗W202-τ2∗/τ102q¯3+W203-τ2∗/τ102q¯2+e-iω10τ2∗W112-τ2∗τ10q3+W113-τ2∗τ10q2+3K33q22q¯2e-iω10τ2∗+K34q2q2q¯3+2q¯2q3e-iω10τ2∗,where (45)W20θ=ig20ω10τ10q0eiω10τ10θ+ig¯023ω10τ10q¯0e-iω10τ10θ+E1e2iω10τ10θ,W11θ=-ig11ω10τ10q0eiω10τ10θ+ig¯11ω10τ10q¯0e-iω10τ10θ+E2,E11=M11M1,E12=M12M1,E13=M13M1,E21=M21M2,E22=M22M2,E23=M23M2,M1=2iω10-a11-b11e-2iω10τ10-a12-b13e-2iω10τ10-b21e-2iω10τ102iω10-a22-b23e-2iω10τ100-b32e-2iω10τ2∗2iω10-a33-b33e-2iω10τ2∗,M11=2q3e-2iω10τ10K11+q3K12-a12-b13e-2iω10τ10-q3e-2iω10τ10K11+q3K122iω10-a22-b23e-2iω10τ10q2e-2iω10τ2∗q2K31+q3K32-b32e-2iω10τ2∗2iω10-a33-b33e-2iω10τ2∗,M12=22iω10-a11-b11e-2iω10τ10q3e-2iω10τ10K11+q3K12-b13e-2iω10τ10-b21e-2iω10τ10-q3e-2iω10τ10K11+q3K12-b23e-2iω10τ100q2e-2iω10τ2∗q2K31+q3K322iω10-a33-b33e-2iω10τ2∗,M13=22iω10-a11-b11e-2iω10τ10-a12q3e-2iω10τ10K11+q3K12-b21e-2iω10τ102iω10-a22-q3e-2iω10τ10K11+q3K120-b32e-2iω10τ2∗q2e-2iω10τ2∗q2K31+q3K32,M2=a11+b11a12b13b21a22b230b32a33+b33,M21=-K11q3+q¯3-2K12q3q¯3a12b13K11q3+q¯3+2K12q3q¯3a22b23-K32q2q¯3+q¯2q3-2K31q2q¯2b32a33+b33,M22=a11+b11-K11q3+q¯3-2K12q3q¯3b13b21K11q3+q¯3+2K12q3q¯3b230-K32q2q¯3+q¯2q3-2K31q2q¯2a33+b33,M23=a11+b11a12-K11q3+q¯3-2K12q3q¯3b21a22K11q3+q¯3+2K12q3q¯30b32-K32q2q¯3+q¯2q3-2K31q2q¯2.Thus, we can determineW20(θ) and W11(θ). Furthermore, we can compute g21 by (44). Then we can compute the following values: (46)c10=i2ω10τ10g20g11-2g112-g0223+g212,μ2=-Rec10Reλ′τ10,T2=-Imc10+μ2Imλ′τ10ω10τ10,β2=2Rec10.From [17], we know that μ2 determines the directions of the Hopf bifurcation: if μ2>0(μ2<0), then the Hopf bifurcation is supercritical (subcritical); β2 determines the stability of the bifurcating periodic solutions: the bifurcating periodic solutions are stable (unstable) if β2<0(β2>0); and T2 determines the period of the bifurcating periodic solutions: the period increases (decreases) if T2>0(T2<0).


## 4. Numerical Simulations

From Section3, we can determine the direction of a Hopf bifurcation and the stability of the bifurcating periodic solutions. In this section, we will give some numerical simulations of system (3) at different values of time delays.We choose the coefficients as follows:μ=0.1, K=50, β=0.1, α=0.2, d=0.2, β1=0.1, α1=0.1, Λ=15, γ=0.01, and m=0.3. α, d, α1, and γ are taken from [15]. Then system (3) has an endemic equilibrium E∗(19.5683,27.6652,35.5003). When τ2=0, we then have ω10∗=0.312 and τ10∗=4.4686. From Theorem 4, we know that E∗ is asymptotically stable when τ1<τ10∗, which is illustrated in Figure 1.Figure 1
E ∗ is asymptotically stable, where τ1=4.3<τ10∗=4.4686 and τ2=0.Whenτ1 passes through the critical value τ10∗, E∗ loses its stability and a Hopf bifurcation occurs; that is, a family of periodic solutions bifurcate from E∗ (see Figure 2). Similarly, we get ω20∗=0.6679 and τ20∗=2.4472.Figure 2
System (3) undergoes a Hopf bifurcation at the endemic equilibrium E∗, where τ1=7.3>τ10∗=4.4686 and τ2=0.Regardτ1 as a parameter, for τ2=2∈(0,2.4472), and then we have ω10=0.3116 and τ10=7.4719. Theorem 6 shows that E∗ is asymptotically stable when τ1<τ10 (see Figure 3) and unstable when τ1>τ10. From formulae (46) in Section 3, it follows that c1(0)=-5.3591×10-5+2.2474×10-5i, μ2=0.0086, β2=-1.0718×10-4, and T2=1.1427×10-4. Since μ2>0 and β2<0, the Hopf bifurcation is supercritical, and these bifurcating periodic solutions from E∗ at τ10 are stable, which are depicted in Figure 4.Figure 3
E ∗ is asymptotically stable, where τ1=7.1<τ10=7.4719 and τ2=2.Figure 4
System (3) undergoes a Hopf bifurcation at the endemic equilibrium E∗, where τ1=7.6>τ10=7.4719 and τ2=2.


## 5. Discussion

In this paper, we have studied the dynamics of a plant virus propagation model with two delays (3) proposed by Jackson and Chen-Charpentier [15]. The model describes the disease transmission dynamics between the insects and the plants.Jackson and Chen-Charpentier [15] studied model (3) using numerical methods. However, the problem of the theoretical analysis of this model remained unsolved and was an open problem.For this problem, first, by analyzing the characteristic equation, constructing a Lyapunov functional, and using LaSalle’s invariance principle, we prove that the disease-free equilibriumE0 is globally asymptotically stable if R0≤1 (Theorem 2), regardless of the length of the time delays, the sufficient conditions for the stability of the endemic equilibrium, and existence of Hopf bifurcation if R0>1 have been given, respectively. Then, by using normal form theory and center manifold theorem introduced by Hassard et al. [17], regarding τ1 as a parameter, we investigate the direction and stability of the Hopf bifurcation, and the explicit formulae which determine the direction and stability of the bifurcating periodic solutions are derived. Finally, the numerical simulation results in Figures 1–4 have verified the obtained analytic results.Our simulation results show that, for the parameter values considered, the disease will persist and exhibit oscillatory bahavior, and this manifests that the densities of the plants and insect-vectors will remain in an oscillatory case, and then agriculture workers must be alert to the virus even if they have noticed that fewer plants are becoming infected.



---

*Source: 7126135-2018-04-08.xml*