![MathJax Logo](/templates/jsp/_style2/_tandf/pb2/images/math-jax.gif)
Abstract
In this paper, we investigated a new heroin–cocaine epidemic model which incorporates spatial heterogeneity and nonlinear incidence rate. The main project of this paper is to explore the threshold dynamics in terms of the basic reproduction number , which was defined by applying the next-generation operator. The threshold type results shown that if
, then the drug-free steady state is globally asymptotically stable. If
, then heroin–cocaine spread is uniformly persistent. Furthermore, the globally asymptotic stability of the drug-free steady state has been established for the critical case of
by analysing the local asymptotic stability and global attractivity.
1. Introduction
Over the past several decades, mathematical modelling as an important tool has been successfully applied to investigate infectious diseases which have posed imminent health threats to human. However, the abuse of illicit drugs is also detrimental on society, such as heroin, cocaine, which need to obtain more attention from the public health agency. When modelling heroin addiction and spread, Mackintosh and Stewart [Citation24] have shown that classical epidemic models can also be applied to study drug addiction and informing control strategies. Based on the research, White and Comiskey [Citation41] proposed an ODE model to consider the heroin epidemics by separating the population into three parts: susceptible individuals (S), drug users not in treatment () and drug users in treatment (
). The model takes the following form:
(1)
(1) The details of parameters of model (Equation1
(1)
(1) ) are described in [Citation41]. The authors defined the basic reproduction number
of the model and the sensitivity analysis are carried out. Moreover, the existence of backward bifurcation at the endemic equilibrium investigated for
, and the authors concluded that effective prevention is indeed better than cure for blocking drugs abuse. Motivated by this work, Ma et al. [Citation21–23] proposed and studied some models which incorporating relapse, media impact and psychological addicts, and the global dynamics and existence of bifurcation have been analysed.
When taking time delay for interpreting the time needed for relapse from treatment individuals () turn to the untreated state (
) into consideration on the base of model (Equation1
(1)
(1) ), some DDE models have been proposed and studied (see [Citation6,Citation9,Citation17,Citation18,Citation20,Citation30,Citation46,Citation47]) and the authors mainly focused on the global dynamics and hopf bifurcation of these models. Moreover, to investigate the impact of treat age or susceptibility age on the heroin abuse, age-structured models have been proposed on the base of model (Equation1
(1)
(1) ) by researchers (see [Citation2,Citation3,Citation5,Citation10,Citation11,Citation19,Citation39]). However, the above-mentioned models are all involved in a homogeneous environment. Actually, the mobility of host population and spatial heterogeneity are also important and meaningful factors when modelling disease spreading. For example, Allen et al. [Citation1] have proposed a diffusive SIS epidemic model and studied the effect of spatial heterogeneity on the persistence and extinction of the disease. Song et al. [Citation33] have investigated an SEIR model in heterogeneous environment and pointed that the movement and spatial heterogeneous of the exposed and recovered individuals have a great effect on the disease dynamics. For more details about such models, one can refer literature [Citation8,Citation13,Citation29,Citation36,Citation37,Citation45] and references cited in. Hence, it is necessary to incorporate spatial heterogeneity for modelling. When considering the effect of spatial heterogeneity on heroin spreading, motivated by model (Equation1
(1)
(1) ), Duan et al. proposed and studied an age-structured heroin epidemic model with reaction–diffusion and the threshold dynamics of the model has been investigated in term of the corresponding basic reproduction number
[Citation8]. Wang and Sun in [Citation37] proposed and investigated a heroin epidemic model with reaction–diffusion and homogeneous Neumann boundary conditions in a bounded domain (Ω) based on model (Equation1
(1)
(1) ). By using the theory of the next-generation operator, the basic reproduction number
has been defined, which is regard as the threshold that determines whether the heroin extinct or not. The threshold dynamics in the light of
are achieved. Furthermore, the dynamics for a critical case of
equal to 1 has been improved by Xu and Geng in [Citation44].
The above-mentioned studies of drug abuse dynamics models mostly focus on one kind of drugs such as heroin. In fact, heroin and cocaine users are overlapping populations, i.e. many heroin users are also use cocaine at the same time [Citation15]. When considering the effect of interaction between the epidemic of heroin and cocaine spread, Chekroun et al. [Citation3] have proposed and studied an age-structured model with both the two kinds of drugs, the threshold dynamics have been studied in term of the basic reproduction number. However, spatial heterogeneity and movement of the host population have not been taken into consideration in [Citation3]. To the best of our knowledge, the heroin–cocaine co-infection models in which population density depends on both time and space variables are rarely studied. Motivated by [Citation3,Citation37], we propose an reaction–diffusion model which describing the dynamics of heroin and cocaine epidemics, the model takes the following form:
(2)
(2) with the homogeneous Neumann boundary conditions
(3)
(3) and initial conditions
(4)
(4) Here
is a bounded domain with smooth boundary denoted by
.
,
,
,
denote the density of susceptible, heroin users, cocaine users and recovered individuals at time t and location x, respectively.
represent the corresponding diffusion rates of
,
,
,
.
and
are the incidence rates of heroin users and cocaine users, respectively.
is the recruitment rate.
is the rate that an individual recovered from the consumption of heroin becomes a cocaine user,
is the rate that an individual recovered from the consumption of cocaine becomes a heroin user.
,
,
and
are the death rates of
,
,
and
, respectively.
and
are the recovery rates of heroin users and cocaine users, respectively.
and
are the rates at which a recovery individual becomes a heroin user and cocaine user, respectively. All the location-dependent parameters are continuous and strictly positive defined on
. The incidence rates of heroin users and cocaine users are given by general nonlinear functions
and
satisfy the following conditions:
(A1) |
| ||||
(A2) |
|
The rest of this paper is organized as follows. The well-posedness of model (Equation2(2)
(2) ) is established and the basic reproduction number
is defined in Section 2. In Section 3, we investigate the threshold dynamics of model (Equation2
(2)
(2) ) in term of
. A brief conclusion ends the paper.
2. Well-posedness
Let be the Banach space with the supremum norm
and
. Then
is a ordered Banach space. Denote
,
,
,
.
Denote ,
,
,
, and
are the Green function associated with
subject to the Neumann boundary condition. Define
are the
semigroup associated with
subject to the Neumann boundary condition. Then we have
It follows from [Citation31] that
is compact and strongly positive for all t>0, then there exists a constant M>0 such that
where
are the principle eigenvalue of
subject to the Neumann boundary condition.
Define operators by
where
. Let
, then we can rewrite model (Equation2
(2)
(2) ) as
(5)
(5) where
and
. It then follows from [Citation27] that the following lemma holds.
Lemma 2.1
For any , model (Equation2
(2)
(2) ) has a unique mild solution
on the interval of
where
. Moreover, this solution is a classical solution.
Based on Lemma 2.1, we are ready to present the following result.
Theorem 2.2
Model (Equation2(2)
(2) ) admits a unique solution
on
for any
with
. Furthermore, the solution semiflow
of model (Equation2
(2)
(2) ) is defined by
admits a global compact attractor.
Proof.
Based on Lemma 2.1, assume , then we have
as
by Theorem 2 in [Citation27]. It follows from the first equation of model (Equation2
(2)
(2) ) that
By the comparison principle and Lemma 2 in [Citation13], there exists a constant
such that
, for
. Then we have
Considering the following system:
(6)
(6) It follows from the standard Krein–Rutman theorem [Citation7] that model (Equation6
(6)
(6) ) admits one principle eigenvalue
which corresponding to a strongly positive eigenfunction
, and hence, model (Equation6
(6)
(6) ) has a solution
for
, where
is a constant satisfying
for
. Then it follows from the comparison principle that
Thus there exists a constant
such that
,
,
, for
, which leads to a contradiction.
Denote be the eigenvalue of
with eigenfunction of
which satisfying
. It follows from [Citation12] that
. Since
is uniformly bounded, then there exists a
such that
, for t>0.
In addition, assume that be the eigenvalue of
subject to the Neumann boundary condition, which satisfies
. It then follows from [Citation35] that
for
. Furthermore, for some
we have
It follows from the comparison principle [Citation13] that there exists a
and
such that
,
,
. Now, we define
then we have
where
. Thus there exist
and
such that
for
. Let
, for all
, utilizing the comparison principle and (Equation5
(5)
(5) ), we obtain
which implies that
Similarly, there exists two positive constants
such that
,
.
The above discussions imply that the system is point dissipative. In addition, it follows from [Citation42] that the solution semiflow is compact for t>0. Then, it follows from [Citation14] that
admits a global compact attractor.
According to [Citation13] that model (Equation2(2)
(2) ) admits a unique drug-free steady state
, linearizing model (Equation2
(2)
(2) ) yields that
(7)
(7) Let
, it then follows from model (Equation7
(7)
(7) ) that
(8)
(8) It follow from Theorem 7.6.1 in [Citation31] that eigenvalue problem (Equation8
(8)
(8) ) admits a principle eigenvalue
with a strictly positive eigenfunction. According to the idea of [Citation38], we will define the basic reproduction number for model (Equation2
(2)
(2) ). To this end, let
and
Assume
,
and
be the
-semigroup generated by the following system:
(9)
(9) Assume that the heroin and cocaine drug users are introduced at time t = 0, where the distribution of initial drug users is
, then
is the distribution of those drug users as time evolves to time t under the synthetical influences of mobility, mortality and transfer of individuals in infected compartments, and hence,
represents the distribution of new drug users at time t. Then, by [Citation38], the total number of new drug users is given by
it follows from [Citation38] that
is a continuous and positive operator which maps the initial distribution of drug users to the total drug users produced during the infection period. Furthermore, the basic reproduction number is defined by
where
represents the spectral radius of the operator
. Consequently, the following result holds followed by [Citation38].
Lemma 2.3
The principal eigenvalue and
have the same sign, and the drug-free steady state
is asymptotically stable when
.
3. Stability and persistence
Based on the above discussions, we first to investigate the global stability of disease-free steady state .
Theorem 3.1
If , then the drug-free steady state
of model (Equation2
(2)
(2) ) is globally asymptotically stable.
Proof.
Since , according to Lemma 2.3 that
. By continuity, there exists a small
such that
, where
is the principle eigenvalue
(10)
(10) It follows from the first equation (Equation2
(2)
(2) ) that
Then, by the comparison principle, there exists a
such that
for all
and
. Furthermore, it follows from model (Equation2
(2)
(2) ) that
(11)
(11) In addition, for any
, assume that there exists a
such that
where
is a strongly positive eigenfunction corresponding to
. By the comparison principle, we can obtain that
Since
, it then follows that
Thus the equation of S in model (Equation2
(2)
(2) ) is asymptotic to
It follows from Lemma 2 in [Citation13] and Corollary 4.3 in [Citation34] that
, and hence, it completes the proof.
Lemma 3.2
Let be the solution of model (Equation2
(2)
(2) ) with
. If there exists some
such that
,
and
, then
,
and
for
and
. Furthermore, there exists some
such that
for any
.
Proof.
By applying the similar method of [Citation13,Citation45], it follows from Lemma 2.1 and equation of and R that
and
and
It follows from the strong maximum principle and Hopf boundary lemma in [Citation28] that
,
,
. Furthermore, according to the proof of Theorem 2.2, there exists a positive constant
and
such that
,
,
for
. Then, from model (Equation2
(2)
(2) ) we have
Considering the following comparison system:
(12)
(12) it follows from the Lemma 2 in [Citation13] that system (Equation12
(12)
(12) ) admits a unique positive steady state
which is globally asymptotically stable. Thus, the standard parabolic comparison implies that
, and hence, the proof is complete.
Theorem 3.3
For every initial value function , assume model (Equation2
(2)
(2) ) admits a solution
on
. If
, then there exists a constant
such that for any
with
,
and
, then we have
Moreover, model (Equation2
(2)
(2) ) admits at least one positive steady state.
Proof.
Let
and
For
, it follows from Lemma 3.2 that
,
and
. That is,
. Moreover, assume that
and
is the omega limit set of orbit
.
Claim 1. for
.
Since , then it follows that
, and hence,
or
or
. For
, it then follows from model (Equation2
(2)
(2) ) that
, then we have
. Furthermore, it follows from model (Equation2
(2)
(2) ) that the equation of S is asymptotic to
which implies that
. If
for some
and
, then Lemma 3.2 implies that
for
and
, and hence,
or
. For the case of
, then model (Equation2
(2)
(2) ) implies that
, and thus
, which leads a contradiction. Furthermore, considering the case where
for some
then
for
, and thus
. It follows from (Equation2
(2)
(2) ) that
, which is also a contradiction. Thus the claim is valid.
Recall that , then Lemma 2.3 implies that
. By continuity, there exists a sufficient small
such that
, where
is the principle eigenvalue of
(13)
(13) Claim 2.
is uniformly weak repeller in the sense that
,
.
If it is not true, then there exists such that
, and hence, there exists a
such that for
,
,
,
,
. It follows from model (Equation2
(2)
(2) ) that
Assume
is the strongly positive eigenfunction associated with the principle eigenvalue
. Suppose
for some
, and hence,
where
is a solution of the following linear system:
Then we have
,
,
, this arrives a contradiction.
Define a function by
It follows that
and have the property that if either
or
and
, then
. This implies p is a generalized distance function for the semiflow
(see [Citation32]). Moreover, it follows from the above discussion that any forward orbit of
in
converges to
, thus
is isolated in
and
, where
is the stable set of
. Thus there is no cycle in
from
to
. It follows from Theorem 2.1 that the semiflow
has a global compact attractor in
. By [Citation32], there exists a
such that
for
, which implies that
,
,
. Accordingly, choosing
, and by Lemma 3.2, it then follows that
,
,
. Furthermore, it follows from Lemma 3.2 and [Citation25] that
admits at least a steady state in
, which is a positive steady state of model (Equation2
(2)
(2) ). This completes the proof.
Furthermore, motivated by the idea of [Citation4,Citation26,Citation43], we investigate the global stability of in the critical case of
.
Theorem 3.4
If , and
and
are Lipschitz continuous, then
is globally asymptotically stable.
Proof.
First, we prove the local stability of of model (Equation2
(2)
(2) ). Assume
and set
with
for a small
. Define
Then recall that
, and hence, we have
Assume the operator
admits a positive semigroup
, then it follows from [Citation16] that
for some
. Therefore,
can be rewritten by
(14)
(14) where
.
Then we have
where
. Moreover, let
be the semigroup for the following system:
It then follows from [Citation40] that there exists a positive constant
such that
for
, and hence, we have
where
and
are Lipschitz constants for
and
, respectively. Denote
Then, we have
where
, thus we have
according to Gronwall's inequality, it follows that
Thus, from the first equation of model (Equation2
(2)
(2) ), we can obtain
Assume
is the solution of the following system:
(15)
(15) Then, the comparison principle implies that
for
,
. Let
be the positive steady state of (Equation15
(15)
(15) ) and
, then
satisfies
It then follows that
which implies that
where
. By Gronwall's inequality, we have
By selecting
sufficiently small such that
, then we have
(16)
(16) It follows that
Recall that
, then we have
, and hence,
Since
, then we can selecting a sufficiently small σ such that
which implies the local asymptotic stability of
.
Next, we prove the is globally attractive. According to Theorem 2.2,
has a global attractor
. Define
Claim 1.
, the omega limit set
.
Similar to [Citation4,Citation26,Citation36,Citation43,Citation45], we define
and hence,
for t>0. It claims that
is strictly decreasing. To this end, fixed a
and let
,
and
for
, such that
for
, where
satisfies
(17)
(17) It follows from the comparison principle that
,
and
for
,
. Due to the arbitraryness of
, thus
is strictly decreasing.
Denote and
, then there exists a sequence
with
such that
. Since,
which implies
for
. If
or
or
, a similar argument procedures to the above leads to that
is strictly decreasing. This is a contradiction and hence
, then we have
. Consequently, it follows from the first equation of model (Equation2
(2)
(2) ) that
as
.
Claim 2. .
It follows from the discussions above that is globally attractive in
, and
forms the only compact invariant subset in
. Since the omega limits set
is compact invariant and
for
, thus
. Since the global attractor
is compact invariant, it follows from Lemma 3.11 in [Citation43] that
. Global attractivity and local asymptotic stability immediately lead to the global asymptotically stability of
. This completes the proof.
4. Conclusion
In this paper, a diffusion epidemic model for heroin and cocaine was formulated. The basic reproduction number was defined by using the method of next generation operators. The threshold type dynamics of the model in terms of have performed by applying the comparison arguments and persistence theory: if
, then the drug-free steady state is globally asymptotically stable, while the drug spread is persistent if
. Moreover, by means of Gronwall's inequality, comparison theorem and the properties of semigroup, the global stability of the drug-free steady state in the critical case
was investigated.
Though the spatial heterogeneity was introduced in model (Equation2(2)
(2) ), the age structure was not been considered in the model. As mentioned in [Citation3,Citation8,Citation10,Citation11,Citation19], treat age or susceptibility age on drug abuse is also should be taken into consideration. Thus a model with age structure and spatial heterogeneity is an interesting project, we leave this in future study.
Disclosure statement
No potential conflict of interest was reported by the author(s).
Additional information
Funding
References
- L.J.S. Allen, B.M. Bolker, Y. Lou, and A.L. Nevai, Asymptotic profiles of the steady states for an SIS epidemic reaction-diffusion model, Discrete Contin. Dyn. Syst. Ser. A. 21 (1) (2008), pp. 1–20.
- S. Bentout, S. Djilali, and B. Ghanbari, Backward, Hopf bifurcation in a heroin epidemic model with treat age, Int. J. Model. Simul. Sci. Comput. 12(02) (2021), pp. 1–22.
- A. Chekroun, M.N. Frioui, T. Kuniya, and T.M. Touaoula, Mathematical analysis of an age structured heroin-cocaine epidemic model, Discrete Contin. Dyn. Syst. B. 25(11) (2020), pp. 4449–4477.
- R. Cui, K. Lam, and Y. Lou, Dynamics and asymptotic profiles of steady states of an epidemic model in advective environments, J. Differ. Equ. 263 (4) (2017), pp. 2343–2373.
- A. Din and Y. Li, Controlling heroin addiction via age-structured modeling, Adv. Differ. Equ. 2020 (1) (2020), pp. 1–17.
- S. Djilali, A. Zeb, and T. Saeed, Effect of occasional heroin consumers on the spread of heroin addiction, Fractals 30 (05) (2022), pp. 1–12.
- Y. Du, Order Structure and Topological Methods in Nonlinear Partial Differential Equations, Vol. 1, in: Maximum Principles and Applications, World Scientific, New Jersey, 2006.
- X. Duan, X. Li, and M. Martcheva, Qualitative analysis on a diffusive age-structured heroin transmission model, Nonlinear Anal. RWA. 54 (2020), pp. 1–26.
- B. Fang, X. Li, M. Martcheva, and L.M. Cai, Global stability for a heroin model with two distributed delays, Discrete Contin. Dyn. Syst. Ser. B. 19 (2014), pp. 715–733.
- B. Fang, X. Li, M. Martcheva, and L.M. Cai, Global asymptotic properties of a heroin epidemic model with treat-age, Appl. Math. Comput. 263 (2015), pp. 315–331.
- B. Fang, X. Li, M. Martcheva, and L.M. Cai, Global stability for a heroin model with age-dependent susceptibility, J. Syst. Sci. Complex. 28 (6) (2015), pp. 1243–1257.
- R.B. Guenther and J.W. Lee, Partial Differential Equations of Mathematical Physics and Integral Equations, Dover Publica-tions, 1996.
- Z. Guo, F.B. Wang, and X. Zou, Threshold dynamics of an infective disease model with a fixed latent period and nonlocal infections, J. Math. Biol. 65 (6-7) (2012), pp. 1387–1410.
- J.K. Hale, Asymptotic Behavior of Dissipative Systems, American Mathematical Society, Providence, 1988.
- P.T. Harrell, B.E. Manchaa, H. Petras, R.C. Trenz, and W.W. Latimer, Latent classes of heroin and cocaine users predict unique HIV/HCV risk factors, Drug. Alcohol. Depend. 122 (3) (2012), pp. 220–227.
- P Hess, Pitman Research Notes in Mathematics Series; Vol. 247, Periodic-Parabolic Boundary Value Problems and Positivity. Longman Scientific and Technical; New York: 1991.
- G. Huang and A. Liu, A note on global stability for a heroin epidemic model with distributed delay, Appl. Math. Lett. 26 (7) (2013), pp. 687–691.
- S. Kundu, N. Kumari, and S. Kouachi, et al. Stability and bifurcation analysis of a heroin model with diffusion, delay and nonlinear incidence rate, Model. Earth Syst. Environ. 8 (1) (2022), pp. 1351–1362.
- L. Liu and X. Liu, Mathematical analysis for an age-structured heroin epidemic model, Acta. Appl. Math. 164 (1) (2019), pp. 193–217.
- L. Liu, X. Liu, and J. Wang, Threshold dynamics of a delayed multigroup heroin epidemic model in heterogeneous populations, Discrete Contin. Dyn. Syst. Ser. B. 21 (8) (2016), pp. 2615–2630.
- M. Ma, S. Liu, and J. Li, Bifurcation of a heroin model with nonlinear incidence rate, Nonlinear Dyn.88(1) (2017), pp. 555–565.
- M. Ma, S. Liu, and J. Li, Does media coverage influence the spread of drug addiction?, Commun. Nonlinear Sci. Numer. Simul. 50 (2017), pp. 169–179.
- M. Ma, S. Liu, and H. Xiang, et al. Dynamics of synthetic drugs transmission model with psychological addicts and general incidence rate, Physica A. 491 (2018), pp. 641–649.
- D.R. Mackintosh and G.T. Stewart, A mathematical model of a heroin epidemic: implications for control policies, J. Epidemiol. Community Health. 33(4) (1979), pp. 299–304.
- P. Magal and X.Q. Zhao, Global attractors and steady states for uniformly persistent dynamical systems, SIAM J. Math. Anal. 37 (1) (2005), pp. 251–275.
- P. Magal, G. Webb, and Y. Wu, On a vector-host epidemic model with spatial structure, Nonlinearity31 (12) (2018), pp. 5589–5614.
- R.H. Martin and H.L. Smith, Abstract functional differential equations and reaction-diffusion systems, Trans. Amer. Math. Soc. 321 (1990), pp. 1–44.
- M.H. Protter and H.F. Weinberger, Maximum Principles in Differential Equations, Springer-Verlag, New York, 1984.
- X. Ren, Y. Tian, L. Liu, and X. Liu, A reaction-diffusion within host HIV model with cell-to-cell transmission, J. Math. Biol. 76 (7) (2018), pp. 1831–1872.
- G.P. Samanta, Dynamic behaviour for a nonautonomous heroin epidemic model with time delay, J. Appl. Math. Comput. 35 (1-2) (2011), pp. 161–178.
- H.L. Smith, Monotone dynamic systems: an introduction to the theory of competitive and cooperative systems, in Math Surveys Monogr, Am Math Soc Providence RI, 1995.
- H.L. Smith and X.Q. Zhao, Robust persistence for semidynamical systems, Nonlinear Anal. 47 (9) (2001), pp. 6169–6179.
- P. Song, Y. Lou, and Y. Xiao, A spatial SEIRS reaction-diffusion model in heterogeneous environment, J. Differ. Equ. 267 (9) (2019), pp. 5084–5114.
- H.R. Thieme, Convergence results and a Poincar–Bendixson trichotomy for asymptotically autonomous differential equations, J. Math. Biol. 30 (7) (1992), pp. 755–763.
- M. Wang, Nonlinear Elliptic Equations, Science Press, Beijing, 2010.
- J. Wang and R. Cui, Analysis of a diffusive host-pathogen model with standard incidence and distinct dispersal rates, Adv. Nonlinear Anal. 10(1) (2021), pp. 922–951.
- J. Wang and H. Sun, Analysis of a diffusive heroin epidemic model in a heterogeneous environment, Complexity 2020 (2020), pp. 1–12.
- W. Wang and X.Q. Zhao, Basic reproduction numbers for reaction–diffusion epidemic models, SIAM J. Appl. Dyn. Syst. 11 (4) (2012), pp. 1652–1673.
- J. Wang, J. Wang, and T. Kuniya, Analysis of an age-structured multi-group heroin epidemic model, Appl. Math. Comput. 347 (2019), pp. 78–100.
- G Webb, Monographs and Textbooks in Pure and Applied Math Series; Vol. 89, Theory of Nonlinear Age-Dependent Population Dynamics. Dekker; New York: 1985.
- E. White and C. Comiskey, Heroin epidemics, treatment and ODE modelling, Math. Biosci. 208 (1) (2007), pp. 312–324.
- J. Wu, Theory and Applications of Partial Functional Differential Equations, Springer-Verlag, New York, 1996.
- Y. Wu and X. Zou, Dynamics and profiles of a diffusive host-pathogen system with distinct dispersal rates, J. Differ. Equ. 264 (8) (2018), pp. 4989–5024.
- J. Xu and Y. Geng, Global stability for a heroin epidemic model in a critical case, Appl. Math. Lett.121 (2021), pp. 1–6.
- Y. Yang, J. Zhou, and C.H. Hsu, Threshold dynamics of a diffusive SIRI model with nonlinear incidence rate, J. Math. Anal. Appl. 478(2) (2019), pp. 874–896.
- Z. Zhang and Y. Wang, Hopf bifurcation of a heroin model with time delay and saturated treatment function, Adv. Differ. Equ. 2019 (2019), pp. 1–16.
- Z. Zhang, F. Yang, and W. Xia, Hopf bifurcation analysis of a synthetic drug transmission model with time delays, Complexity 2019 (2019), pp. 1–17.