Influence of time delay on bifurcation of a synthetic drug transmission model with psychological addicts

A synthetic drug transmission model with psychological addicts and time delay is proposed in this paper. By analyzing the corresponding characteristic equation and choosing the time delay as the bifurcation parameter, a set of sufficient criteria guaranteeing local stability of the synthetic drug addiction equilibrium and the appearance of a Hopf bifurcation of the model is established. Further, the direction and stability of the Hopf bifurcation are investigated with the aid of normal form theory and center manifold theory. Finally, numerical simulations are performed to support the analytical results.


Introduction
It is well known that increased use of synthetic drugs such as crystal methamphetamine, ketamine, ecstasy and others, which are a new type of mental drugs based on chemical synthesis is an issue of concern in many parts of the world. Synthetic drugs are becoming increasingly popular in drug markets as they mainly appear in entertainment venues. Taking China for example, the number of reported synthetic drug abusers increased from 1.459 millon in the end of 2014 to 1.538 million in the end of 2017 by the report of China's Drug Situations Report in the past five years (2014-2018) [1]. The number of methamphetamine abusers was up to 1.35 millon which accounted for 56.1% among the 2.404 million drug-users by the end of 2018 according to China's Drug Situations Report (2018) [2]. What is more serious is that synthetic drugs are more addictive because they can directly affect the central nervous system and show the stronger psychological dependence compared with traditional drugs including heroin, morphine and marijuana, etc. In addition, any drug abuse and dependence constitute one of the most important modes of transmitting human immunodeficiency virus (HIV) and the Hepatitis C virus (HCV) [3][4][5]. Obviously, it is urgent to take effective measures to control the spread of synthetic drugs in order to eradicate the tremendous damages and panics brought by synthetic drugs abuse to social and public health system.
Since the heroin addiction was first defined as epidemic in [1981][1982][1983] in Ireland, the heroin-using can be modeled in a similar way to the epidemic models [6] and mathematical modeling has been used extensively to address issues of public health importance. For example, Ma et al. developed different forms of heroin epidemic models [3,[6][7][8][9][10] to study the transmission of heroin epidemics. Sharomi et al. formulated different smoking models [11][12][13][14][15][16][17] for giving up smoking. There were also some drinking models [18][19][20][21] proposed by Mushayabasa et al. to analyze the influence of binge drinking to public health. Similarly, mathematical modeling can be also used to describe the spread of synthetic drugs. In [22,23], Nyabadza et al. analyzed the methamphetamine transmission in South Africa by constructing a suitable mathematical model. In [24], Liu et al. formulated a synthetic drugs transmission model with treatment and studied global stability and backward bifurcation of the model. Considering the effect of synthetic drugs transmission caused by psychology, Ma et al. [25] proposed the following synthetic drugs transmission model with psychological addicts: in which the total population (N(t)) is divided into four classes at time t such as the susceptible (S(t)), the psychological addicts (P(t)), the physiological addicts (H(t)) and the drug-users in treatment (T(t)). A is the recruitment rate of the susceptible; β(N) is the contact rate; d is the natural death rate of the populations; π , γ , θ and σ are state transition rates. Ma et al. [25] studied stability of system (1).
As stated in [25], a susceptible one is more likely to initiated drug abuse when he contacts with a physiological addict compared to a psychological addict. On the other hand, there should be a period for a drug user in treatment before he relapses into the physiological addicts due to the effect of treatment and his self-control. With an eye to such considerations and motivated by the work of some other dynamical systems with time delay [26][27][28][29][30], we investigate a more realistic synthetic drugs transmission model as follows: where β 1 and β 2 are the contact rates of the psychological addicts and the physiological addicts, respectively. τ is the time delay due to the period that a drug user in treatment uses before he relapses into the physiological addicts.
We organize this paper as follows: in the following section, we analyze the local stability and appearance of Hopf bifurcation of system (2). Section 3 is concerned with the formulas determining the direction and stability of the Hopf bifurcation. Some numerical simulations are executed to illustrate the analytical results in Sect. 4. Finally, conclusions are presented in Sect. 5.

Local stability and existence of Hopf bifurcation
The stability and existence of unique positive equilibrium of system (2) is related to basic reproductive number 0 on free drug equilibrium point (FEP) D 0 , which is determined with the help of next generation matrix method [31]. The free drug equilibrium point is Consider the following matrices for finding the basic reproductive number 0 : .
Now the Jacobians of F and V at D 0 are Here For 0 > 1, system (2) has a unique synthetic drug addiction equilibrium D * (S * , P * , . At the unique synthetic drug addiction equilibrium D * (S * , P * , H * , T * ), system (2) can be linearized as whose characteristic equation is and When τ = 0, the characteristic equation becomes We denote and suppose that It follows from the Routh-Hurwitz theorem that if the condition (C 1 ) holds, all roots of Eq. (5) have negative real parts.

Figure 3
The phase plots of system (33) for different θ at τ = 5.7850. The remaining parameters are taken as given in the text
Thus, it is easy to see that Eq. (34) has one positive root χ 0 = 1.1885. Then we obtain ω 0 = 1.0902 and τ 0 = 9.7367, h (χ 0 ) = 7.7692 > 0. According to Theorem 1, we can con-  (33) is locally asymptotically stable when τ ∈ [0, 9.7367), and a Hopf bifurcation appears when the value of τ passes through the critical value 9.7367, which can be indicated by bifurcation diagrams in Fig. 1. Furthermore, based on Eq. (32), μ 2 = 1.5613 > 0, β 2 = -2.07 < 0 and T 2 = 0.3611 > 0 can be obtained. Accordingly, through Theorem 2, it can be concluded that the Hopf bifurcation at τ 0 is supercritical; the bifurcating periodic solutions are stable and increasing. In addition, Fig. 2 manifests that the numbers of susceptibles and drug-users in treatment decrease, whereas the numbers of the psychological and physiological addicts increase as we increase θ from 0.5 to 0.6. Also, due to the increase in θ from 0.5 to 0.6, system (33) loses its stability and shows limit cycle behavior, which is illustrated in Fig. 3. It is interesting to note that Fig. 4 reveals the reverse effect of increasing σ in a particular range on the four populations in system (33). Figure 5 shows that system (33) changes its behavior from limit cycle to stable focus as we increase the value of σ .

Conclusions
The prevention and control of synthetic drugs have received great attention due to their increasingly serious consequences on population. Quite a few mathematical models describing synthetic drugs transmission have been put up by scholars at home and abroad in