# Hopf Bifurcation Analysis of a Synthetic Drug Transmission Model with Time Delays

Authors: Zizhen Zhang; Fangfang Yang; Wanjun Xia
Journal: Complexity
 (2019)

---

## Abstract

This paper is concerned with the Hopf bifurcation of a synthetic drug transmission model with two delays. Firstly, some sufficient conditions of delay-induced bifurcation for such a model are captured by using different combinations of the two delays as the bifurcation parameter. Secondly, based on the center manifold theorem and normal form theory, some sufficient conditions determining properties of the Hopf bifurcation such as the direction and the stability are established. Finally, to underline the effectiveness of the obtained results, some numerical simulations are ultimately addressed.


---

## Body

## 1. Introduction

Increased use of heroin and other addictive drugs is an issue of concern all over the world. Drug abuse affects not only life quality of the general public but also the overall situation of social stability and economic development [1–3]. The data from the World Drug Report published by the United Nations (U.N.) showed that thirty-five million people worldwide suffer from drug abuse disorders, and only one-seventh has received treatment [4]. It also showed that, in 2017, millions of people around the world injected drugs, including 1.4 million people living with HIV and 5.6 million people suffering from hepatitis C. From the statistical data, it is clear that the adverse health effect due to drug use is more serious and widespread than previous anticipation, and it is urgent to control the prevalence of addictive drugs.As stated in [5], the spread of heroin habituation and addiction has similar characteristics to an epidemic, including rapid diffusion and clear geographic boundaries. With the help of epidemic models, one could try to simulate and reveal the nature of epidemics and provide theoretical rules and results for preventing and controlling infectious diseases [6]. In recent decades, mathematical modelling technologies based on the infectious disease models have been developed to understand and combat drug-addiction problems. In [7], White and Comiskey formulated a heroin epidemic model with a standard incidence rate based on principles of mathematical epidemiology. In the succession, Mulone and Straughan [8] found stability conditions for steady states of the proposed model by White and Comiskey. In the subsequent papers, Wang et al. [9–12] studied global stability of a heroin epidemic model with bilinear incidence rate, respectively. Some other works related to the dynamical behaviour of heroin epidemic models with nonlinear incidence rates can be found in [3, 13–15], and models with age structure can be found in [13, 15–18]. For the analytical study of stochastic heroin epidemic models or some other heroin epidemic models, one can also see [1, 2, 19–22].All the aforementioned models consider only traditional drugs. Compared with traditional drugs, the relevant works are few. On the other hand, synthetic drugs are addictive more easily because they can directly affect the central nervous system as a new type of mental drug. According to China’s Drug Situation Report [23], among drug abuse, methamphetamine is the most common. Methamphetamine abusers accounted for 56.1% of the existing 2.444 million drug addicts, and methamphetamine has replaced heroin as the most abused drug in China. Besides, in 2018, the people who relapsed and abused synthetic drugs were accounted for 57.3% of the total number of drug abuse among the relapsed addicts. Based on the above facts, Ma et al. [24] formulated a synthetic drug transmission model with psychological addicts and general contact rates, and they investigated the local and global stability of the proposed model. However, they assumed that the contact rates of psychological addicts and physiological addicts are equal for the sake of simple calculation and analysis, which is not reasonable because the susceptible individuals who have never taken any drugs are more likely to initiate drug abuse once they contact with the physiological addicts compared to the psychological ones. Based on the work in [24], Saha and Samanta [25] proposed a synthetic drug transmission model with general contact rate and Holling type-II functional responses among susceptible and drug addicts. Recently, Liu et al. [26] proposed the following synthetic drug transmission model with different susceptible compartments:(1)dStdt=Λ−β1StIt−μSt,dQtdt=εIt+δRt−μQt−β2QtIt,dItdt=β1StIt+β2QtIt+σRt−ε+γ+μIt,dRtdt=γIt−δ+σ+μRt,where the meanings of St, Qt, It, and Rt are described concisely in Table 1. All the parameters Λ, β1, μ, ε, δ, γ, β2, and σ are positive constants, and their meanings are presented in Table 2. Liu et al. [26] studied the global exponential stability of the drug-free equilibrium and the global stability of the drug-addiction equilibrium, and they showed that special psychological treatment played a very positive role in drug abusers.Table 1
The state variables for system (1).
Parameter
Description
S
t
Number of susceptible individuals who have never taken any drugs at timet
Q
t
Number of susceptible individuals who have history of drug abuse
I
t
Number of drug users not in treatment
R
t
Number of individuals who are receiving treatmentTable 2
The description of parameters for system (1).
Parameter
Description
Λ
Immigration rate of the susceptible
β
1
Probability of transmission fromS to I
μ
Natural death rate of all populations
Ε
Self-cure rate fromI to Q
δ
Successful treatment rate fromI to Q
γ
Procession rate fromI to R
β
2
Probability of transmission fromQ to I
σ
Probability of treatment failureIt should be pointed out that system (1) assumed that drug abusers can give up drugs instantaneous. This assumption seems not to be realistic since the synthetic drugs are addictive easily and it usually needs a period to give up drugs. Time delays have been incorporated into dynamical models about some other fields by many scholars [27–33]. Generally speaking, delay differential equations exhibit much more complicated dynamics than ordinary differential equations since a time delay could cause the equilibrium of a dynamical model to lose its stability. Hence, it is important to know the critical point at which a synthetic drug transmission model changes its stability. Based on the discussion above, we consider the following synthetic drug transmission model with two time delays:(2)dStdt=Λ−β1StIt−μSt,dQtdt=εIt−τ1+δRt−τ2−μQt−β2QtIt,dItdt=β1StIt+β2QtIt+σRt−γ+μIt−εIt−τ1,dRtdt=γIt−σ+μRt−δRt−τ2,where τ1 is the time delay due to the period that the drug abusers use to give up drugs through self-control. τ2 is the time delay due to the period used to give up drugs through successful treatment. The transfer diagram of system (2) is depicted as in Figure 1. This paper mainly concerns the effect of the delays τ1 and τ2 on the stability of system (2).Figure 1
The transfer diagram of system (2).The framework of the current paper is arranged as follows. In Section2, the existence of Hopf bifurcation is discussed in detail by using the different combinations of the two delays as the bifurcation parameters and analyzing the distribution of roots of the associated characteristic equations. In Section 3, the direction of Hopf bifurcation and stability of the bifurcating periodic solutions are determined with the help of the center manifold theorem and normal form theory. In Section 4, the effectiveness of the obtained theoretical findings is certified through numerical simulations. Finally, conclusions are drawn in Section 5.


## 2. The Existence of Hopf Bifurcation

According to the analysis in the literature [26], we know that system (2) has a unique positive synthetic drug-addiction equilibrium E∗S∗,Q∗,I∗,R∗ when(3)R0=β2Λδ+σ+μμε+γ+μδ+σ+μ−γσ>1,S∗=Λβ1I∗+μ,Q∗=εδ+σ+μ+δγδ+σ+μμ+β2I∗I∗,R∗=γδ+σ+μI∗,where I∗ is the positive root of the following equation:(4)aI∗2+bI∗+c=0,a=μβ1β2γ+δ+σ+μ,b=μ2β1+β2γ+δ+μ+σ+μεβ1δ+σ+μ+μβ1γδ−β1β2Λδ+σ+μ,c=u2δ+σ+μμ+ε+μ2γδ+μ−μβ1Λδ+σ+μ.The linear part of system (2) is(5)dStdt=L11St+L13It,dQtdt=L22Qt+L23It+M23It−τ1+N24Rt−τ2,dItdt=L31St+L32Qt+L33It+L34Rt+M33It−τ1,dRtdt=L43It+L44Rt+N44Rt−τ2,whose characteristic equation is(6)λ4+a3λ3+a2λ2+a1λ+a0+b3λ3+b2λ2+b1λ+b0e−λτ1+c3λ3+c2λ2+c1λ+c0e−λτ2+d2λ2+d1λ+d0e−λτ1+τ2=0,where(7)a0=L11L22L33L44−L34L43+L13L31L22L44,a1=L34L43−L33L44L11+L22−L13L31L22+L44−L11L22L33+L44,a2=L11+L22L33+L44+L11L22+L33L44−L13L31−L34L43,a3=−L11+L22+L33+L44,b0=L11L44L22M33−L32M23,b1=L32M23L11+L44−M33L11L22+L11L44+L22L44,b2=M33L11+L22+L44−L32M23,b3=−M33,c0=L11L22L33N44+L13L31L22N44+L11L32L43N24,c1=−N44L11L22+L11L33+L22L33−L13L31N44−L32L43N24,c2=N44L11+L22+L33,c3=−N44,d0=L11N44L22M33−L32M23,d1=L32M23N44−M33N44L11+L22,d2=M33N44,L11=−β1I∗+μ,L13=−β1S∗,L22=−μ+β2I∗,L23=−β2Q∗,L31=β1I∗,L32=β2I∗,L33=β1S∗+β2Q∗,L34=σ,L43=γ,L44=−σ+μ,M23=ε,M33=−ε,N24=δ,N44=−δ.Case 1.
τ
1
=
τ
2
=
0. When τ1=τ2=0, equation (6) becomes(8)λ4+K3λ3+K2λ2+K1λ+K0=0,where(9)K0=a0+b0+c0+d0,K1=a1+b1+c1+d1,K2=a2+b2+c2+d2,K3=a3+b3+c3.
According to Routh–Hurwitz criterion, the drug-addiction equilibriumE∗S∗,Q∗,I∗,R∗ is locally asymptotically stable provided the following conditions are satisfied: H1K0>0, K3>0, K2K3>K1K4, and K1K2K3>K12K4+K0K32.Case 2.
τ
1
>
0 and τ2=0.
We analyze the effect ofτ1 on bifurcation for system (2). When τ1>0 and τ2=0, equation (6) becomes(10)λ4+K13λ3+K12λ2+K11λ+K10+L13λ3+L12λ2+L11λ+L10e−λτ1=0,where(11)K10=a0+c0,K11=a1+c1,K12=a2+c2,K13=a3+c3,L10=b0+d0,L11=b1+d1,L12=b2+d2,L13=b3.
Assume thatλ=iω1ω1>0 is a purely imaginary root of equation (10), then it follows that(12)L10−L12ω12cosω1τ1−L13ω13−L11ω1sinω1τ1=K12ω12−ω14,L10−L12ω12sinω1τ1+L13ω13−L11ω1cosω1τ1=K11ω1−K13ω13.
By the aid of equation (12), we have(13)ω18+e13ω16+e12ω14+e11ω12+e10=0,where(14)e10=K102−L102,e11=K112+2L10L12−2K10K12−L112,e12=K122+2K10+2L11L13−2K11K13−L122,e13=K132−2K12−L132.
Making the substitutionv1=ω12, equation (13) can be rewritten as(15)v14+e13v13+e12v12+e11v1+e10=0.
Define(16)g1v1=v14+e13v13+e12v12+e11v1+e10.
In order to establish the main results of this article, it is assumed thatH21 equation (15) has at least one positive real root. Without loss of generality, assume that equation (15) has four positive roots, denoted by v11, v12, v13, and v14. Accordingly, ω1i=v1i i=1,2,3,4 are the roots of equation (13). From equation (12), one has(17)τ1ij=1ω1i×arccosf11ω1if12ω1i+2nπ,with i=1,2,3,4; n=0,1,2,…; and(18)f11ω1i=L12−K13L13ω1i6+K11L13+K13L11−K12L12−L10ω1i4+K10L12+K12L10−K11L11ω1i2−K10L10,f12ω1i=L132ω1i6+L122−2L11L13ω1i4+L112−2L10L12ω1i2+L102.
Denote(19)τ10=τ1i00=minτ1i0∣i=1,2,3,4,ω10=τ1i0.
Next, we derive the condition of the occurrence for Hopf bifurcation. Differentiating both sides of equation (10) with regard to τ1 and substituting λ=iω10 into the obtained expression, it can be achieved that(20)Redλdτ1τ1=τ10−1=g1′v10L13ω103−L11ω102+L10−L12ω1022,where v10=ω102. To ensure the condition of the occurrence for Hopf bifurcation, we educe the following hypothesis: H22 g1′v10≠0, which suggests that transversality condition is matched. In conclusion, we have the following results.Theorem 1.
For system (2), if the conditions H1, H21, and H22 hold, then drug-addiction equilibrium E∗S∗,Q∗,I∗,R∗ is locally asymptotically stable when τ1∈0,τ10; system (2) undergoes a Hopf bifurcation at the drug-addiction equilibrium E∗S∗,Q∗,I∗,R∗ when τ1=τ10, and a family of periodic solutions bifurcate from the drug-addiction equilibrium E∗S∗,Q∗,I∗,R∗.Remark 1.
Although the assumptionsH21 and H22 seem to be tedious, however, one can verify the assumptions in numerical simulations.Case 3.
τ
1
=
0 and τ2>0.Remark 2.
Whenτ1=0 and τ2>0, the analysis of the effect of τ2 on bifurcation for system (2) is similar as that in Case 2. Therefore, we omit it here.Case 4.
τ
1
=
τ
2
=
τ
>
0.
We analyze the effect of τ on bifurcation for system (2). When τ1=τ2=τ>0, equation (6) becomes(21)λ4+K33λ3+K32λ2+K31λ+K30+L33λ3+L32λ2+L31λ+L30e−λτ+M32λ2+M31λ+M30e−2λτ=0,with(22)K30=a0,K31=a1,K32=a2,K33=a3,L30=b0+c0,L31=b1+c1,L32=b2+c2,L33=b3+c3,M30=d0,M31=d1,M32=d2.
By multiplyingeλτ on both sides of equation (21), it can be procured that(23)L33λ3+L32λ2+L31λ+L30+λ4+K33λ3+K32λ2+K31λ+K30eλτ+M32λ2+M31λ+M30e−λτ=0.
Supposeiωω>0 is the root of equation (23), then one has(24)J41ωcosωτ+J42ωsinωτ=J45ω,J43ωsinωτ−J44ωcosωτ=J46ω,where(25)J41ω=ω4−K32+M32ω2+K30+M30,J42ω=K33ω3−K31−M31ω,J43ω=ω4−K32−M32ω2+K30−M30,J44ω=K33ω3−K31+M31ω,J45ω=L32ω2−L30,J46ω=L33ω3−L31ω.
Using Cramer’s rule to solve the above equations, one obtains(26)cosωτ=l6ω6+l4ω4+l2ω2+l0ω8+j6ω6+j4ω4+j2ω2+j0,(27)sinωτ=l7ω7+l5ω5+l3ω3+l1ωω8+j6ω6+j4ω4+j2ω2+j0,with(28)j0=K302−M302,j2=K312−M312−2K30K32+2M30M32,j4=2K30−M322−2K31K33,j6=K332−2K32,l0=−K30−M30L30,l1=K31+M31L30−K30+M30L31,l2=K32−M32L30+K30−M30L32−K31−M31L31,l3=K32+M32L31+K30+M30L33−K31+M31L32−K33L30,l4=K31−M31L33−K32−M32L32+K33L31−L30,l5=K33L32−L31−K32−M32L33,l6=L32−K33L33,l7=L33.
Squaring both sides of equations (26) and (27), respectively, and adding them together, one has(29)ω16+k7ω14+k6ω12+k5ω10+k4ω8+k3ω6+k2ω4+k1ω2+k0=0,where(30)k0=j02−l02,k1=2j0j2−2l0l2−l12,k2=j22+2j0j4−l22−2l0l4−2l1l3,k3=2j0j6+2j2j4−2l0l6−2l2l4−l32−2l1l5,k4=j42+2j0+2j2j6−l42−2l2l6−2l1l7−2l3l5,k5=2j2+2j4j6−2l4l6−l52−2l3l7,k6=j62+2j4−l62−2l5l7,k7=2j6−l72.
Letω2=v, then equation (29) becomes(31)v8+k7v7+k6v6+k5v5+k4v4+k3v3+k2v2+k1v+k0=0.
Suppose thatH41 equation (29) has at least one positive root. Without loss of generality, assume that equation (31) has eight positive roots, denoted by v1,v2,…,v8. Accordingly, ωi=vi i=1,2,…,8 are the roots of equation (29). Thus, one has(32)τij=1ωi×arccosl6ωi6+l4ωi4+l2ωi2+l0ωi8+j6ωi6+j4ωi4+j2ωi2+j0+2nπ,with i=1,2,…,8; n=0,1,2,…; and(33)τ0=τi00=minτi0∣i=1,2,…,8,ω0=τi0.
Differentiating both sides of equation (23) with respect to τ yields(34)dλdτ−1=−f31λf32λ−τλ,where(35)f31λ=3L33λ2+2L32λ+L31+4λ3+3K33λ2+2K32λ+K31eλτ+2M32λ+M31e−λτ,f32λ=λλ4+K33λ3+K32λ2+K31λ+K30eλτ−λM32λ2+M31λ+M30e−λτ.
Hence, one obtains(36)Redλdττ=τ0−1=U3RV3R+U3IV3IV3R2+V3I2,with(37)U3R=L31−3L33ω02+K31+M31−3K33ω02cosω0τ0−2K32ω0−2M32ω0−4ω03sinω0τ0,U3I=2L32ω0+K31−M31−3K33ω02sinω0τ0+2K32ω0+2M32ω0−4ω03cosω0τ0,V3R=K33ω04−K31ω02+M31ω02cosω0τ0−ω05−K32ω03−M32ω03+K30ω0+M30ω0sinω0τ0,V3I=K33ω04−K31ω02−M31ω02sinω0τ0+ω05−K32ω03+M32ω03+K30ω0−M30ω0cosω0τ0.
In order to obtain the main results, it is necessary to make the following extra assumptionH42 U3RV3R+U3IV3I≠0, which ensures that transversality condition holds.Theorem 2.
For system (2), if the conditions H1, H41, and H42 hold, then drug-addiction equilibrium E∗S∗,Q∗,I∗,R∗ is locally asymptotically stable when τ∈0,τ0; system (2) undergoes a Hopf bifurcation at the drug-addiction equilibrium E∗S∗,Q∗,I∗,R∗ when τ=τ0, and a family of periodic solutions bifurcate from the drug-addiction equilibrium E∗S∗,Q∗,I∗,R∗.Remark 3.
Although the assumptionsH41 and H42 seem to be tedious, however, one can verify the assumptions in numerical simulations.Case 5.
τ
1
∈
0
,
τ
10 and τ2>0.
Motivated by the work in [34], in this case, we analyze the effect of τ2 on bifurcation for system (2) and fix τ1∈0,τ10. Suppose iω2′ω2′>0 is the root of equation (6). For convenience, we still denote ω2′ as ω2. Then, one has(38)J41ω2sinω2τ2+J42ω2cosω2τ2=J43ω2,J41ω2cosω2τ2−J42ω2sinω2τ2=J44ω2,where(39)J41ω2=c1ω2−c3ω23+d1ω2cosω2τ1−d0−d2ω22sinω2τ1,J42ω2=c0−c2ω22+d1ω2sinω2τ1+d0−d2ω22cosω2τ1,J43ω2=a2ω22−ω24−a0−b1ω2−b3ω23sinω2τ1−b0−b2ω22cosω2τ1,J44ω2=a3ω23−a1ω2−b1ω2−b3ω23cosω2τ1+b0−b2ω22sinω2τ1.
Squaring both sides of the above equations, respectively, and adding them together, one has(40)r40ω2+2r41ω2cosω2τ1+2r42ω2sinω2τ1=0,where(41)r40ω2=ω28+a32−2a2−b32−c32ω26+a22+2a0−2a1a3+b22−2b1b3−c22+2c1c3−d22ω24+a12−2a0a2−2b0b2+b12−c12+2c0c2−d12+2d0d2ω22+a02+b02−c02−d02,r41ω2=a3b3−b2ω26+a2b2−a1b3−a3b1+b0−c2d2+c3d1ω24+a1b1−a2b0−a0b2+c0d2−c1d1+c2d0ω2+a0b0−c0d0,r42ω2=−b3ω27+a2b3+b1−a3b2+c3d2ω25+a3b0+a1b2−a2b1−a0b3−c1d2+c2d1−c3d0ω23+a0b1−a1b0−c0d1+c1d0ω2.
Suppose thatH51 equation (40) has finite positive roots, denoted by ω21∗,ω22∗,…,ω2i∗. Thus,(42)τ2i∗j=1ω2i∗×arccosJ41ω2i∗×J44ω2i∗+J42ω2i∗×J43ω2i∗J412ω2i∗+J422ω2i∗+2nπ,with i=1,2,…,k; n=0,1,2,…; and(43)τ2∗=τ2i∗00=minτ2i∗0∣i=1,2,…,k.
Differentiating equation (6) with respect to τ2, one has(44)dλdτ−1=−f41λf42λ−τ2λ,with(45)f41λ=4λ3+3a3λ2+2a2λ+a1+3b3−τ1b2λ2−τ1λ3+2b2−τ1b1λ+b1−τ1b0e−λτ1+2d2−τ1d1λ−τ1d2λ2+d1−τ1d0e−λτ1+τ2+3c3λ2+2c2λ+c1e−λτ2,f42λ=λc3λ3+c2λ2+c1λ+c0e−λτ2+λd2λ2+d1λ+d0e−λτ1+τ2.
Thus, one obtains(46)Redλdτ2τ2=τ2∗−1=U4RV4R+U4IV4IV4R2+V4I2,with(47)U4R=2c2ω2∗sinω2∗τ2∗+c1−3c3ω2∗2cosω2∗τ2∗+a1−3a3ω2∗2+2d2−τ1d1ω2∗sinω2∗τ1cosω2∗τ2∗+cosω2∗τ1sinω2∗τ2∗+τ1d2ω2∗2+d1−τ1d0cosω2∗τ1cosω2∗τ2∗−sinω2∗τ1sinω2∗τ2∗+τ1ω2∗3+2b2−τ1b1ω2∗sinω2∗τ1+τ1b2−3b3ω2∗2+b1−τ1b0cosω2∗τ1,U4I=2c2ω2∗cosω2∗τ2∗−c1−3c3ω2∗2sinω2∗τ2∗+2a2ω2∗−4ω2∗3+2d2−τ1d1ω2∗cosω2∗τ1cosω2∗τ2∗−sinω2∗τ1sinω2∗τ2∗−τ1d2ω2∗2+d1−τ1d0sinω2∗τ1cosω2∗τ2∗+cosω2∗τ1sinω2∗τ2∗+τ1ω2∗3+2b2−τ1b1ω2∗cosω2∗τ1−τ1b2−3b3ω2∗2+b1−τ1b0sinω2∗τ1,V4R=d0ω2∗−d2ω2∗3sinω2∗τ1cosω2∗cosω2∗+cosω2∗τ1sinω2∗τ2∗−d1ω2∗2cosω2∗τ1cosω2∗τ2∗−sinω2∗τ1sinω2∗τ2∗,V4I=d0ω2∗−d2ω2∗3cosω2∗τ1cosω2∗cosω2∗−sinω2∗τ1sinω2∗τ2∗+d1ω2∗2sinω2∗τ1cosω2∗τ2∗+cosω2∗τ1sinω2∗τ2∗.
As in Case 4, we make the following extra assumptionH52U34V4R+U4IV4I≠0, which ensures that transversality condition holds.Theorem 3.
For system (2), if τ1∈0,τ10 and τ2>0 and the conditions H1, H51, and H52 hold, then drug-addiction equilibrium E∗S∗,Q∗,I∗,R∗ is locally asymptotically stable when τ2∈0,τ2∗; system (2) undergoes a Hopf bifurcation at the drug-addiction equilibrium E∗S∗,Q∗,I∗,R∗ when τ2=τ2∗, and a family of periodic solutions bifurcate from the drug-addiction equilibrium E∗S∗,Q∗,I∗,R∗.Remark 4.
Although the assumptionsH51 and H52 seem to be tedious, however, one can verify the assumptions in numerical simulations.


## 3. Properties of Hopf Bifurcation

In this section, by employing the center manifold theorem and normal form theory, the properties of the Hopf bifurcation at the critical valueτ2∗ are determined. Let t=sτ, St=u1sτ, Qt=u2sτ, It=u3sτ, Rt=u4sτ, and τ=τ2∗+ϱ where ϱ∈R. Throughout this section, we assume that τ1∗<τ2∗, where τ1∗∈0,τ10. Then, system (2) becomes the following form:(48)u˙t=Lϱut+Fϱ,ut,where Lϱ: C⟶R4 and F: R×C⟶R4 are defined, respectively, by(49)Lϱϕ=τ2∗+ϱM1maxϕ0+M2maxϕ−τ1∗τ2∗+M3maxϕ−1,Fϱ,ϕ=−β1ϕ10ϕ30−β2ϕ20ϕ30β1ϕ10ϕ30+β2ϕ20ϕ300,with(50)M1max=L110L1300L22L230L31L32L33L3400L43L44,M2max=000000M23000M3300000,M3max=0000000N240000000N44.By the Riesz representation theorem in [35], there exists a function ηθ,ϱ of bounded variation for θ∈−1,0 such that(51)Lϱϕ=∫−10dηθ,ϱϕθ.In fact, one can choose(52)ηθ,ϱ=τ2∗+ϱM1max+M2max+M2max,θ=0,τ2∗+ϱM2max+M3max,θ∈τ1∗τ2∗,0,τ2∗+ϱM3max,θ∈−1,−τ1∗τ2∗,0,θ=−1,where δ is defined by(53)δθ=0,θ≠0,1,θ=0.Forϕ∈C−1,0,R4, define(54)Aϱϕ=dϕθdθ,−1≤θ<0,∫−10dηθ,ϱϕθ,θ=0,Rϱϕ=0,−1≤θ<0,Fϱ,ϕ,θ=0.Then, system (48) is equivalent to(55)u˙t=Aϱut+Rϱut,where ut=ut+θ.The adjoint operatorA∗ϱ of Aϱ is defined by(56)A∗φ=dϕsds,0<s≤1,∫−10dηTs,0ϕ−s,s=0.Forϕ∈C−1,0,R4 and φ∈C10,1,R4∗, define(57)φs,ϕθ=φ¯0ϕ0−∫θ=−10∫ξ=0θφ¯ξ−θdηθϕξdξ,where ηθ=ηθ,0.From Section2, it can be seen that ±iτ2∗ω2∗ are the eigenvalues of A0, so ±iτ2∗ω2∗ are also the eigenvalues of A∗0. Suppose that qθ=1,q2,q3,q4Teiτ2∗ω2∗θ and q∗s=D1,q2∗,q3∗,q4∗Teiτ2∗ω2∗s be the eigenvectors for A0 and A∗0 corresponding to +iτ2∗ω2∗ and −iτ2∗ω2∗, respectively. Then, one can obtain(58)q2=L33iω2∗−L11L13iω2∗−L22−M23e−iτ1∗ω2∗−N24e−iτ2∗ω2∗,q3=iω2∗−L11L13,q4=L43iω2∗−L11L13iω2∗−L44−N24e−iτ2∗ω2∗,q2∗=L11L32+M33eiτ1∗ω2∗iω2∗+L31iω2∗+L22+M23eiτ1∗ω2∗,q3∗=−L11iω2∗+L31,q4∗=−N24eiτ2∗ω2∗q2∗+L34q3∗iω2∗+L44+N44eiτ2∗ω2∗.In order to guaranteeq∗s,qθ=1, the value of D needs to be determined. In view of equation (57), we have(59)D¯=1+∑i=14qi∗q¯i+τ1∗eiτ1∗ω2∗q¯2M23q2∗+M33q3∗+τ2∗eiτ2∗ω2∗q¯4N24q2∗+N44Q4∗−1.In what follows, based on the algorithm in [36] and the similar computation process as that in [37–40], we can obtain the expressions of g20, g11, g02, and g21 as follows:(60)g20=2D¯τ2∗q¯3∗−1β1q3+q¯3∗−q¯2∗β2q2q3,g11=D¯τ2∗q¯3∗−1β1q3+q¯3+q¯3∗−q¯2∗β2q2q¯3+q¯2q3,g02=2D¯τ2∗q¯3∗−1β1q¯3+q¯3∗−q¯2∗β2q¯2q¯3,g21=2D¯τ2∗q¯3∗−1β1W1110q3+12W2010q¯3+W1130+12W2030+q¯3∗−q¯2∗β2W1120q3+12W2020q¯3+W1130q2+12W2030q¯2,with(61)W20θ=ig20q0τ2∗ω2∗eiτ2∗ω2∗θ+ig¯02q¯03τ2∗ω2∗e−iτ2∗ω2∗θ+E1e2iτ2∗ω2∗θ,W11θ=−ig11q0τ2∗ω2∗eiτ2∗ω2∗θ+ig¯11q¯0τ2∗ω2∗e−iτ2∗ω2∗θ+E2.E
1 and E2 can be solved by(62)E1=2L11∗0−L1300L22∗L23∗−N24e−2iτ2∗ω2∗L31−L32L33∗L3400−L43L44∗−1×−β1q3−β2q2q3β1q3+β2q2q30,E2=2L110L1300L22L23+M23N24L31L32L33+M33L3400L43L44+N44−1×−β1q3+q¯3−β2q2q¯3+q¯2q3β1q3+q¯3+β2q2q¯3+q¯2q30,where(63)L11∗=2iω2∗−L11,L22∗=2iω2∗−L22,L23∗=−L23−M23e−2iτ1∗ω2∗,L33∗=2iω2∗−L33−M33e−2iτ1∗ω2∗,L44∗=2iω0−L44−N44e−2iτ2∗ω2∗.Thus, we have(64)C10=i2τ2∗ω2∗g11g20−2g112−g0223+g212,μ2=−ReC10λ′τ2∗,β2=2ReC10,T2=ImC10+μ2Imλ′τ2∗τ2∗ω2∗.In conclusion, we can obtain the following results based on the fundamental results about Hopf bifurcation in the literature [32].Theorem 4.
For system (2), if μ2>0 μ2<0, then the Hopf bifurcation is supercritical (subcritical); if β2<0 β2>0, then the bifurcating periodic solutions are stable (unstable); if T2>0 T2<0, then the period of the bifurcating periodic solutions increases (decreases).


## 4. Numerical Simulation

In this section, numerical simulations are employed to confirm the efficiency of the theoretical analysis in the present paper. Considering the biological significance of the parameters in system (2), we refer to the range of the parameter values provided in the literature [26]. By extracting some values from the literature [26] and considering the biological significance of the parameters and conditions for the occurrence of Hopf bifurcation, we choose Λ=1, β1=0.01, μ=0.02, ε=0.1, δ=0.095, β2=0.008, σ=0.0011, and γ=0.095. Then, system (2) turns into(65)dStdt=1−0.01StIt−0.02St,dQtdt=0.1It−τ1+0.095Rt−τ2−0.02Qt−0.008QtIt,dItdt=0.01StIt+0.008QtIt+0.0011Rt−0.115It−0.1It−τ1,dRtdt=0.095It−0.0211Rt−0.095Rt−τ2.Then, one can obtainR0=2.3530>1, and equation (4) becomes the following form:(66)3.3776e−007I∗2−3.6335e−006I∗−1.3277e−005=0,from which one gets the unique positive root I∗=13.6396. Further, we obtain that system (65) has a unique drug-addiction equilibrium E∗6.3940,18.7755,13.6396,11.1607. By delicate calculation, it is obtained that ω10=0.1872 and τ10=16.0145; ω20=0.9243 and τ20=21.9566; ω0=2.9207 and τ0=8.6947; ω2∗=1.4648 and τ2∗=10.7568 when τ1∗=7.75∈0,τ10.By Theorem1, E∗6.3940,18.7755,13.6396,11.1607 is asymptotically stable for system (65) when τ1=15.2608<τ10, which is depicted in Figure 2. E∗6.3940,18.7755,13.6396,11.1607 is unstable for system (65) and Hopf bifurcation occurs when τ1=16.3642>τ10, which are simulated in Figure 3. Similar simulations can be shown as in Figures 4 and 5 for Theorem 2, Figures 6 and 7 for Theorem 3, and Figures 8 and 9 for Theorem 4, respectively.Figure 2
Waveform plots of system (65) with τ1=15.2608<τ10.
(a)
(b)
(c)
(d)Figure 3
Waveform plots of system (65) with τ1=16.3642>τ10.
(a)
(b)
(c)
(d)Figure 4
Waveform plots of system (65) with τ2=20.3665<τ20.
(a)
(b)
(c)
(d)Figure 5
Waveform plots of system (65) with τ2=23.9608>τ20.
(a)
(b)
(c)
(d)Figure 6
Waveform plots of system (65) with τ=8.4565<τ0.
(a)
(b)
(c)
(d)Figure 7
Waveform plots of system (65) with τ=8.7807>τ0.
(a)
(b)
(c)
(d)Figure 8
Waveform plots of system (65) with τ2=10.2508<τ2∗ and τ1=7.75∈0,τ10.
(a)
(b)
(c)
(d)Figure 9
Waveform plots of system (65) with τ=11.3691>τ2∗ and τ1=7.75∈0,τ10.
(a)
(b)
(c)
(d)Whereafter, by some complex calculations, we obtainC10=−13.0664+9.6207i and λ′τ2∗=1.0081−0.6309i. By the results in equation (64), it can be derived that μ2=12.9614>0, β2=−26.1328<0, and T2=−0.0916<0. It follows from Theorem 4 that the Hopf bifurcation is supercritical since μ2>0, the bifurcating periodic solutions are stable since β2<0, and the period of the periodic solutions decreases as τ2 increases since T2<0.


## 5. Conclusions

In this paper, a delayed synthetic drug transmission model with relapse and treatment is investigated by incorporating two delays into the model proposed in the literature [26]. We consider not only the time delay due to the period that the drug abusers use to give up drugs through self-control but also the time delay due to the period used to give up drugs through successful treatment. Compared with the model in [26], the delayed synthetic drug transmission model in the present paper is more general because it usually needs a period for the drug abusers to give up drugs through either self-control or successful treatment.It has been shown that, under certain conditions, the drug-addiction equilibrium is locally asymptotically stable when the value of the time delay is below the critical value. In this case, system (2) is in ideal stable state and the synthetic drug transmission can be controlled easily. However, once the value of the time delay is above the critical value, system (2) will lose its stability and undergo a Hopf bifurcation at the corresponding critical value of the time delay, which is not welcomed in reality. The occurrence of Hopf bifurcation means that the existence of populations in system (2) changes from the drug-addiction equilibrium to a limit cycle, and in this case, the synthetic drug transmission is out of control. Therefore, it is vital to take some necessary measures to postpone and eliminate the occurrence of the Hopf bifurcation for system (2).Specially, the direction of the Hopf bifurcation and stability of the bifurcating periodic solutions are determined by employing the center manifold theorem and normal form theory. In addition, according to the numerical simulations, it is easily observed that the time delay due to the period that the drug abusers use to give up drugs through self-control is marked because the critical value ofτ1 is much smaller when we only consider it. From this point of view, it is strongly recommended that drug abusers should have strong will in the process of giving up drugs.



---

*Source: 3492589-2019-11-22.xml*