Global Dynamics of a Discretized Heroin Epidemic Model with Time Delay Academic research paper on "Mathematics"

Hindawi Publishing Corporation Abstract and Applied Analysis Volume 2014, Article ID 742385, 10 pages http://dx.doi.org/10.1155/2014/742385

Research Article

Global Dynamics of a Discretized Heroin Epidemic Model with Time Delay

Xamxinur Abdurahman, Ling Zhang, and Zhidong Teng

College of Mathematics and System Sciences, Xinjiang University, Urumqi, Xinjiang 830046, China Correspondence should be addressed to Xamxinur Abdurahman; xamxinur@sina.com Received 24 June 2014; Accepted 15 September 2014; Published 20 October 2014 Academic Editor: Wanbiao Ma

Copyright © 2014 Xamxinur Abdurahman et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

We derive a discretized heroin epidemic model with delay by applying a nonstandard finite difference scheme. We obtain positivity of the solution and existence of the unique endemic equilibrium. We show that heroin-using free equilibrium is globally asymptotically stable when the basic reproduction number R0 <1, and the heroin-using is permanent when the basic reproduction number R0 >1.

1. Introduction

As we all know, the use ofheroin and other drugs in Europe and more specifically in Ireland and the resulting prevalence are well documented [1-3]. It shows that the use of heroin is very popular and causes many preventable deaths. Heroin is so soluble in the fat cells that it crosses the blood-brain barrier within 15-20 seconds, rapidly achieving a high level syndrome in the brain and central nervous system which causes both the "rush" experience by users and the toxicity. Heroin-related deaths are associated with the use of alcohol or other drugs [4]. Treatment ofheroin users is a huge burden on the health system of any country.

We often study infectious diseases with mathematical and statistical techniques; see, for example, [5-11]; however, little has been done to apply this method to the heroin epidemics. In 1979, Mackintosh and Stewart [9] considered an exponential model which is simplified from infectious disease model of Kermack and McKendrick to illustrate how the heroin-using spreads in epidemic fashion. They arranged a numerical simulation to show how the dynamics of spread are influenced by parameters in the model. White and Comiskey [5] attempted to extend dynamic disease modeling to the drug-using career and formulated an ordinary differential equation. They divided the whole population into three classes, namely, susceptible, heroin users, and heroin users undergoing treatment. Their model allows a steady state

(constant) solution which represents an equilibrium between the number of susceptible, heroin users, and heroin users in treatment. Furthermore, this ODE model was revisited by Mulone and Straughan [12]; the authors proved that this equilibrium solution is stable both linearly and nonlinearly under the realistic condition in which relapse rate of those in treatment returning to untreated drug use is greater than the prevalence rate of susceptible becoming drug users. Recently, the study of the global properties and permanence of continuous heroin epidemic models attracted the researchers and have some very good results; see [13-16]. Specially, Samanta [15] considered a model with time-dependent coefficients and with different removal rates for three different classes, introduced some new threshold values R* and R*, and obtained the permanence of heroin-using career.

Motivated by Samanta [15] and Zhang and Teng [8], we alter a nonautonomous heroin epidemic model with time delay to an autonomous heroin epidemic model. For convenience, we replace U1 and U2 by U and V, respectively. Thus, we obtain the following continuous heroin epidemic model with a distributed time delay:

S(t) = X-ß1 (U)S(t) I U(t-s)dq(s)

U(t) = ß1 (U)S(t) \ U(t-s)dq(s) Jo

+ ß3U (t)V (t) - fa + P + i;1)U(t), V (t) - PU (t) - ß3U (t) V (t) - fa + Ç2) V (t),

where S(t), U(t), and V(t) represent the number of susceptible, heroin users not in treatment, and heroin users in treatment, respectively. We assume that the time taken to become heroin user is s. The function q(s) : [0,h] ^ [0, >x) is nondecreasing and has bounded variation such that

n(s)ds = n(h) - n(0) = 1.

For understanding more realistic phenomenon of heroin, a little complicated epidemic model is helpful. By applying Micken's nonstandard discretization method [17] to continuous heroin epidemics model with time delay (1), we derive the following discretized heroin epidemic model with a distributed time delay:

The initial conditions of the system (2) are given by

c - J1) u -vp v

an - fn , Un - fn , Vn - fn ,

for n - -h, -h + 1,

where > 0 (n = -h, -h + 1,... ,0, i = 1,2, 3). Again, by biological meaning, we further assume that > 0 for all i = 1,2,3.

The paper is organized as follows. In Section 2, we prove the positivity and boundedness of the solution of system (2). In Section 3, we deal with the global asymptotic stability of the heroin-using free equilibrium. In Section 4, we consider the permanence of the discrete epidemic model applying Wang's technique. In the discretized epidemic model, sufficient condition for global asymptotic stability and permanence are the same as for the original continuous epidemic model. We give some numerical examples and conclusion in Sections 5 and 6.

Sn+1 -Sn -X-ßi (u„) Sn+1 lun_ktlk

-ViSn+i +ïlU„+1 +Ï2 Vn+V h

Un+1 -Un -ßi (un)sn+i Yun-knk

+ p3Un+1Vn+1 -fa +P + Z1)Un+1,

Vn+1 -Vn = PUn+1 - &Un+1Vn+1 - fa + & Vn+1,

where Sn is the susceptible class, Un is the class of heroin users not in treatment, and Vn is the class of heroin users in treatment at nth step. Since the sufficient condition can be obtained, independently of the choice of a time step-size, we let the time step-size be one for the sake of simplicity. The nonnegative constants p2, and ^3 denote the death rate of the susceptible, heroin users not in treatment, and heroin users in treatment class, respectively. Throughout the paper, it is biologically natural to assume that ^1 < min{^2,^3|. The constant X > 0 denotes the recruitment rate of susceptible population from the general population. Constant P > 0 is the proportion of heroin users who enter the treatment class. The individuals in treatment who stop using heroin are susceptible at a constant rate > 0. Constant p3 represents the transmission rate from heroin users in treatment to untreated heroin users. p1(Un) is the probability per unit time and the transmission is used with

the form ^(Un)Sn+1 Ho Un- _kqk, which includes various delays. By a natural biological meaning, we assume that ^1(U) is a positive function and that there exists a constant Up > 0 such that p1 (U) is nondecreasing on the interval [0, Up]. The integer h> 0 is the time delay. The sequence % : < Vk <

(k = 0,1,... ,h) is nondecreasing and has bounded.

2. Basic Properties

For system (2), the heroin-using free equilibrium is given by

E0 -(S0,0,0), S0 --.

Define a positive constant A = Xk=o Ik. The stability of E0 is studied by using the next generation method in [7]. The associated matrix F (of the new heroin-using terms) and the M-matrix V (of the remaining transfer terms) are given as follows, respectively:

ßi (0)XA

p2 + P + ti

Clearly, F is nonnegative, V is a nonsingular M-matrix, and V - F has Z sign pattern. The associated basic reproduction number, denoted by ^0,isthengivenby^0 = p(FV-1),where p is the spectral radius of FV- . It follows that

R - ßi (0) AA 0 Hi {V2 +P+tiY

Now, we will consider the positivity and boundedness of solution to system (2). For most continuous epidemic models, positivity of the solution is clear, but, for system (2), the positivity of the sequences Sn, Un, and Vn holds in some condition.

Lemma 1. Let (Sn, Un, Vn) be any solution of system (2) with initial condition (3); then (Sn, Un, Vn) is positive for any n e N and VQ < P(1 + + £2)/p3.

Proof. Let (Sn, Un, Vn) be any solution of system (2) with initial condition (3). It is evident that system (2) is equivalent to the following iteration system:

Un+1 =

Vn+1 =

$n+1 =

P1 (Un)Sn+1 tk=Q Un-k1k +Un

1-psVn+l + + P + $1 ' ^n + PUn+1

1 + Ps Un+1 +H-3 2

X + Sn + ZiUn+1 + ty,n+1 1 + 01 (Un) lhk=Q Un-kVk + to1

In the following, we will use the induction to prove the positivity of solution. When n = 0,we have

V1 = S1 =

h (Uq)S1 tk=Q U-k*lk + Uq

l-ps V1 +^2 +P + Z2 '

Vq + PU1

1 + PsU1 + 3 +&

X + Sq + № +Z2V1 1 + 01 (Uq) IHk=Q U-kVk + to1

From (8)-(10), we see that, as long as U1 is obtained, V1 and S1 will be obtained too.

If U1 > 0, from (9), we directly obtain V1 > 0 and, from (10), we further obtain that S1 > O.Furthermore, wealsohave N1 =S1 +U1 +V1 > 0.

Let x = U1; then, from (8)-(10), we see that x satisfies the following equation:

y(x) = x-((l + l31 (Uq) JU-rfk + ^1)Uq

+P1X U-knkW1 (x)) (11)

x((l + p1 (UQ)Jy-ktik + ^1)^2 (x)) ,

W1 (x) = X + Sq + Z1X + Z2 v1,

W2 (X)=1+fr +P+$1

VQ +Px

1+PsX + tos +&

Substituting V1 in W2(x) gets

W2 (x)=1+^2 +$1 +

P(l+to +^2)-PsVq 1+tos 2 + PsX

Since VQ < P(1 + + ¡;2)/p3, we have

W2 (0)=1+fo2 +P+$1 -

1+tos 2

> 0, (14)

and because of W1(0) = X + SQ + ¡;2VQ/(1 + + ) > 0, thus Y(0) = -((l + p1 (Uq) XU-knk + to1) Uq

(Uq) tu-k1kW1 (0) )

x((l+01 (Uo)XU-knk + to1 )W2 (0)) <0.

Substituting W1 and W2 in f(x) yields

y(x) =

ax2 + bx- c

(1+p1 (Uq) tl=Q U-knk + to1) (mx + n) Here, constants a, b, c, m, and n are as follows:

a=(l + fo2)(1+p1 (UQ)X_P-knk + (l+fo1)>0, b = (l+fos + Q (1 + fo1)(l + fo2 + P + ^ + (1 + ^1) (p1 (Uq) XU-kVk + PsSq - PsNq

-P1 (UQ)XU-k1kP(l+fo3)

1 (Uo)XU-knk (I + Nq),

c = 01 (Uq) YU-k*lk [(1 + tos + S2) (I + Nq) -(1+ to) Vq]

- Uq (l + to)$1,

m = ps (l + fo2 + $1) > 0,

n=(l + to + P + t1)(l + to + ^2)-PsVq >0.

We take the limit on both sides of the above equation: lim f(x)

, ax2 +bx- с = lim т-т-г-

* - (l + ß^Uo) iL U-k1k + fai ) (тх + n)

= +œ>;

this means that f(x) = 0 has at least one positive solution x e (0, +от). So, we have =x > 0. Therefore, the positivity of Sx >0, > 0, and V1 > 0 is finally obtained. When n = 2, we have

ßis2 !k=o ui-кПк + Щ l-ß3V2 + ^2 + P + V

^ + PU 2

l+ßsU2 + fa3 +&

A + St +^U2 +%2V2 1 + ßi Ito U1 -кПк + fa

a similar argument as in the above for , V1, and Sj; we also can obtain that U2 > 0, V2 > 0, and S2 > 0. Lastly, by using the induction, we can finally obtain that Sn > 0, Un > 0, and Vn > 0, for all n > 0.

Now, we define the total population as Nn = Sn + Un + Vn. Then, from system (2), we know that

K+1 -K = X- - faUn+1 - faVn+1- (20)

Notice the assumption that fa < min(fa, fa3); we obtain

N„ <

l+fai X

l+Иг

+ ••• +

l +fai

l +fai

= ±{l-fai

l +fai

< max { —, N0

If X/fa > N^itiseasytoseethat N,, < X/fa = S0,forall large n. If X/fa < N0, from right hand side of system (2), we obtain

X + N0 l +Hi

Hence, we have N1 < N0 and there exists i e N such that Nj < X/fa = SQ. Then, we may use this Nj as a starting value instead of N0. This argument leads to the following result. □

Lemma 2. For any solution (Sn, Un, Vn) of system (2), the total population Nn = Sn + Un + Vn satisfies

lim sup Nn < S0 =

П—+Ж fai

thus (Sn, Un, Vn) is ultimately bounded.

Let a = {(Sn,Un,Vn) : Sn,Un,Vn > 0, Sn + Un + Vn < X/fa}; then a is the positive invariant set to the solution of system (2).

In the following, we will examine the existence of endemic equilibrium for a special case of system (2).

Lemma 3. Assume that p1(U) = p1 > 0 is a constant. If R0 > 1, system (2) admits a heroin-using equilibrium E* = (S*,U*,V*) when V < P/p3, where E* satisfies following equality:

X - S*U* - faS* + ^U* + $2V* = 0, A^1S*U* + p3U* V* -(fa + P + ^) U* = 0, (24) PU* - p3U*V* - (fa + $2) V* = 0.

Proof. Consider the following equation:

X - A^SU - faS + ^U + $2V = 0, A^SU + p3UV - (p2 +P + $1)U = 0, (25) PU-faUV-fa +$2)V = 0.

From the first equation and the second equation ofthe system (25), we have

X - faS + $2V + ¡33UV - (fa +P)U = 0. (26)

From the second equation and the third equation of the system (25), we obtain

A^SU-fa + S1 )U-(K +^2)V = 0; (27)

(Иъ +Ъ)У AßiS -fa +ti.

Since U = 0, from the second equation of the system, we have

fa + P + ti - ß3V Aßi

Substituting S in (28), we obtain

(to +Ï2)V

Substituting U and S in (26) yields a quadratic equation of V as follows:

F (V) = aV +bV + c = 0, where the coefficients are given by

a = toPl - Ap!p3to> c = top (to + P + ti) - XA^P,

b = XA0& - toPs (to + Si) + ^fiito (to + P) + ^ito^i-

Since R0 = XAp1/to1(to2 + P + > 1, then it is easy to see that c < 0 and b > 0. According to Descartes' rule of signs, if a > 0, then F(V) = 0 has a positive solution; if a < 0, then F(V) = 0 has two positive solutions. From the expression of S and U, we note that V < P/p3. Since

F (0) = c < 0, P(P) = A^iPto2 (to >0, (33) \P3 ' P3

This means that F(V) = 0 has a unique positive solution V* e (0,P/p3 ). Therefore, there exists a unique positive solution (S* ,U* ,V*) of system (2). □

For the local stability of the equilibria, we refer to Theorem 2 in [7] and have the following results.

Theorem 4. Assume that ft (U) = p1t ft is a positive constant. The heroin-using free equilibrium E0 = (S0, 0,0) of system (2) is locally asymptotically stable if R0 < 1 and unstable if R0 > 1.

3. Global Asymptotic Stability of the Heroin-Using Free Equilibrium

In this section, we still assume that p1(U) = ft > 0, and obtain a sufficient condition for global asymptotic stability of the heroin-using free equilibrium E of system (2).

Using a Lyapunov function similar to that in [11], we can easily prove the global asymptotic stability of the heroin-using free equilibrium E .

Theorem 5. If R0 < 1, the drug-using free equilibrium E0 of system (2) is globally asymptotically stable.

Proof. Let us take the following Lyapunov function:

Hn = un + Clvn + c2£( tuX + 2(Sn+1 - s0)2,

k=0\l=n-k

where c; >0 (i = 1,2, 3) are the constants to be defined later and S0 = Xjto. Using system (2), the difference of Hn satisfies

AH = Hn+i - Hn

= Un+i - Un + Ci (Vn+1 -Vn) h

+ C2 X (Un+ink - Un-kVk)

+ 2 {(s*! s0)2-{sn-s0)2}.

From Sn < Nn < S0, for all n > 0,we have h

kH < YUn-kVk [PiSn+i -C2 - c3piSn+i (sn+i - S°)}

+ (l-Ci)p3Un+iVn+i

+ [C2A -to -P-ïi +CiP + C3Ï1 (sn+i - S°)] Un

- [01 (to +Ï2)-C3 Ï2 (sn+i -S0 )]vn_

-c3to (Sn+1 -s°)2

^ YUn-k*1k {p1Sn+1 -C2 - C3p1Sn+1 {Sn+1 - ^

+ (1-C1)p3Un+1Vn+1

+ M -to2 -P-$1 +<1 P] Un+1 - C1 (to3 + Q Vn+1 -C3to1 (Sn+1 -S°)2-

Let us choose ct > 0 (i = 1,2,3) such that these constants satisfy the following inequalities:

PA+1 -C2 -C3P1 Sn+1 (Sn+1 -S0)<0, (37)

1-c1 <0, (38)

G1P + c2A<to2 +P + S1. (39)

From (37), we have c3p1^2+1 + (c3p1S0 -p1)Sn+1 +c2 > 0;since Sn+1 > 0, then the following inequality is true:

ßi(l + c3S0)2 <4C2C3;

that is,

ß1(s0) c32 + (2ßls0 - 4c2) c3 +ß1 < 0. (41)

Since R0 <1, which implies that p1AS0 < to2 + P + ^1,we can choosec2 = ^1S0 + e;here,e (0 < e < (to2 -A/31 S0)/A) isa sufficiently small positive number such that p1AS0+Ae < to2+ (1-c1)P + ^1.Since^1S0 -2c2 < 0 and (ftS0-2c2)2 > (p1S0)2, we can choose c3 > 0 to satisfy (41). We may further choose c1 >1 to satisfy (38). Therefore, AH is negative definite and is

equal to zero if and only if Sn+1 = S0, Un+1 = 0, and Vn+1 = 0.

(34) The proof is complete.

4. Permanence of System (2)

The system (2) is said to be permanent if there are positive constants m and M such that

m < lim inf Sn < lim sup Sn < M

hold for any sequence Sn of the system (2), and the same inequalities hold for Un and Vn. For each class Sn, Un, and Vn, m and M are independent of initial conditions.

Following the method used by Wang in [6], we will prove the permanence of system (2) for the general case; that is, assume that ¡31 (U) is related to U.

Theorem 6. If R0 > 1, then system (2) is permanent for any initial condition (3).

Proof. Firstly, from system (2) and Lemmas 1 and 2, for any e0 > 0, there exists sufficiently large n0 > 0 such that Un < X/fa1 + e0 as n>n0 + h. Then, we have

^n+1 =

X + Sn + ^Un+1 + ^K+1

1+ßl (Un)!to Un-klk + fa1

Since R0 = ß1 (0) AX/fa1(fa2+*P+£1) > 1,there exist 0 < a <Uß and p > 0 such that

Aßi (0) X

fa2 + P + ^ fa1 + aß1 (a) A

1 \ph 1 + fa1 + aß1 (a) A J

note that

fa1 + aß1 (a) A

1 + fa1 + ocß1 (oc) A

-. (49)

We claim that it is impossible that Un < a holds for all n>n1 > [ph]. The function [x] gives the smallest integer not less than x. Suppose the contrary, for n>n1 + h. Consider

c _ X + Sn + $1 Un+1 Vn+1 „ X + Sn

^n+1 - . „ N >

1 +fa +ß1 (Un) th=0 Un_knk 1 + fa1 + «ß1 № A

1 + fa + ocß1 (a) A

1 + fa1 + ocß1 (a) A

1+fa1 +ß1 (Un) Ik=0 Un-kIk

Let ß™(e0) - maxU£[o,\lfil+ea]ß1(U)- Thus we have

Sn+1 > j

1 + fa +№ Zto Un-knk _X_

> 1 + fa + pf (X/fa +eo)A' Notice that e0 can be arbitrarily small. Then, we have

lim inf S^ >ms -

1+fa1 +aß1 (o)A) ni+h+1'

1 + fa1 + ocß1 (a) A

i-n, -h

n-nl-h-1

From Lemma 1, Sn satisfies A

^n+1 >

fa1 + aß1 (a) A

1 + fa1 + ocß1 (oc) A

n-n:-h

1 + fa1 + ßf4 (AX/fa) (45) and we have that, for n > n1 + h + [ph], we have

ff - max ft (U). 1 U£[0,X/P! ] 1

Next, let us consider the positive sequences Sn and Un of (2). According to these sequences, we define

Hn -Un +fa^lA £ ¡1 uV (46)

A k=0\l=n-k )

Then, for n > 0,we obtain

^n+1 >

fa1 + aß1 (a) A X

fa1 + aß1 (a) A

1 + fa1 + ocß1 (a) A J 1 ph

1 + fa1 + ocß1 (a) A

Hence, for n>n1 + h+ [ph], we have

AH - Hn+1 - Hn

-Un+1 -Un +

fa2 + P + ^ Y (Un+1*lk -Un-knk)

AH > (ß1 (Un) Sn+1 - fa2 + PA+^ ) Y Un-kflk v A ' k=0

>(ß1 (0)SA -fafa^+^\YUn-kVk

- (ß1 (Un) Sn+1 - fa2 + PA+^ ) YUn-kVk + ßsUn+1Vn+1. V A ' k=0

k=0 Ah

fa2 +P + Ï1 (Aß, (0)SA \YT

Let e = ming{Uni+[ph]+h+g; d = —h, —h+1,...,0}. Now, we will show that Un > e for all n>n1 + [ph] + h. In fact, there is an integer n > 0 such that

Un > e, n1 + [ph] + h <n<n1 + [ph] + h + n, Un+1 < e, n = n1 + [ph] + h + n.

However, for n = n1 + [ph] + h + n,we have

If n2 - n1 <h + [ph], since

we have

ß1 (Un-1) Un-k-1 Ik + Un-1

1+to + P+Ï1 -ßsVn

1+to + P + Ç1'

Un+1 -e =

ß1 (Pn)Sn+1 tk=o Un-knk + Un

1 + to + P + Ï1 -ßsVn+1

ß1 (Un)Sn+1 tk=o Un-knk + Un 1+to + P + Ï1

1+to + P + Ï1 — -£

1+to + P + Ï1

ß1 (0)SAA — (to2 + P + Ï1)

1 + to +P + Ç1

to + P + Ï1

1+to + P + Ï1

ß1 (0)SA A to +P + Ï1

■e>0.

Which is a contradiction. Thus, Un > e for n>n1 + [ph] + h. Therefore, for n>n1 + [ph] + h,

AH > (to + P + Ï1)^

Aß1 (0)SA

to +P + Ï1

- e > 0, (56)

which implies that Hn ^ as n ^ But, from

Lemma 2 and (46), there exists a sufficiently large integer n[ > 0 such that, for n> n[,

to A k=o\l=n-kto .

<±{l + (to2 + P + ^)(h+1)}, to1

which is a contradiction. Hence, the claim is proved.

In the rest, we only need to consider the following two cases:

(i) Un > a for all large n.

(ii) Un oscillates about a for all large n.

We show that Un > mu for all large n, where 0 < mu < a, is a constant which will be given later. Clearly, we only need to consider case (ii). Let positive integers n1 and n2 be sufficiently large that UHi > a, U > a, and Un < a, for n1 < n < n2.

Un 1 + to + P + Ï1

> m„ =

1 + to + P + Ï1 1

1 + to + P + Ï1

h+[ph]

Hence, Un > mu for ne [n1,n2].

If n2 - n1 > h + [ph], we can easily obtain that Un >

for n e [n1,n1 + h + [ph] ]. Assume that there exists an integer n> 0 such that

Un > mu, n1 + h+ [ph] <n<n1 +h + [ph] + n,

Un+1 < mu,

■ = n1 + h+ [ph] + n.

However, for n = n1 + h+ [ph] + n,

Un+1 —mu >

to + P + Ï1 1+to +P + Ï1

ß1 (0)SaA

to + P + Ï1

■mu >0.

This is a contradiction to the proposition that Un+1 < mu. Therefore, Un > mu for n e [n1 ,n2]. Since these positive integers n1 and n2 are chosen in an arbitrary way, we conclude that Un > mu for all large n in case (ii). Hence,

liminfn^JJn > mu.

Note that, from that third equation of system (2), we have P P

lim inf V„ > mV =

to + & + xßi/to1 u to + & + ßss0

From Lemma 2 and the discussion above, we have ms < lim inf Sn < lim sup Sn < S°, mu < lim inf Un < lim sup Un < S0,

mV < lim inf Vn < lim sup Vn < S0. The proof is completed.

20 40 60 80 100 0 5 10 15 20

0.35 -,-,-,-

0.3 - ■ -

0.25 0.2 0.15

0.1 - • -0.05 - ■ , -

0I_■_' • •.......

0 5 10 15 20

Figure 1

5 4.5 4 3.5 3 2.5 2 1.5 1

0.18 0.16 0.14 0.12 0.1 0.08 0.06 0.04 0.02 0

5. Numerical Example

In order to confirm the validity of our results, we consider the following heroin epidemic model with a discrete time delay:

- Sn = X - piSn+1Un_h - faSn+i + %iUn+1 + un+1 -Un = Sn+iUn_h + p3un+iVn+i -(fa +P + Zi )un+i,

Vn+i -Vn = PUn+i - l33Un+iVn+i - (fa + Q vn+i-

Now, we present a numerical example. For the sake of simplicity, we choose the parameters as pi = 0.9, p3 = 0.8, X = 2, fai = 0.1, fa2 = 0.2, fa3 = 0.1, P = 0.4, ^ = 0.1, and = 0.2; we get R0 = 25.7143 < 1. Figure 1 shows that the disease free equilibrium E of the system (64) is globally asymptotically stable when R0 <1. Figure 2 shows that the system (64) is permanent when R0 >1.

6. Conclusions

In this paper, we have modified the Samanta heroin epidemic model into an autonomous heroin epidemic model with distributed time delay. Further, we established a discretized heroin epidemic model with time delay, sufficient conditions have been obtained to ensure the global asymptotic stability of heroin-using free equilibrium when R0 <1 and pi(U) is replaced by a positive constant. We also carried out some discussion about the heroin-using equilibrium, but our results are only restricted to the existence of this equilibrium for ^i(U) = > 0,a special case of system (2). The stability of heroin-using equilibrium is yet to be studied. As a main result of this paper, we obtained the permanence of the system (2). From the expression of R0 = ^i(0)XA/fai(fa2 + P + we see that a decrease in pi (transmission coefficient from susceptible population) will cause a decrease of the same proportion in Rq. If the rate of migration or recruitment is restricted into susceptible community, the spread of the disease can also be kept under control by reducing X and thereby decreasing R0. The spread of the heroin users can also be controlled by educators, epidemiologists, and treatment

9.5 9 8.5 8 7.5 7 6.5 6 5.5 5

Figure 2

providers to increase the values of £ (removal rate of heroin users not in treatment who stop using heroin but are susceptible) and P (proportion of heroin users who enter treatment) and thereby to decrease R0. This analysis tells us that prevention is better than cure; efforts to increase prevention are more effective in controlling the spread of habitual drug use than efforts to increase the numbers of individuals accessing treatment.

Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.

Acknowledgment

This work was supported by the National Natural Science Foundation of China (Grants nos. 11261056, 11261058, and 11271312).

References

[1] European Monitoring Centre for Dru gs and Drug Addiction (CMCDDA), "Annual Report, 2005," http://www.emcdda .europa.eu/publications/annual-report/2005.

[2] C. Comiskey, "National prevalence of problematic opiate use in Ireland," Tech. Rep., EM-CDDA, 1999.

[3] A. Kelly, M. Carvalho, and C. Teljeur, Prevalence of Opiate Use in Ireland 2002-2001. A 3-Source Capture Recapture Study. A Report to the National Advisory Committee on Drugs, Subcommittee o n Prevalence, Small Area Health Research Unit, Department of Public, 2003.

[4] K. A. Sporer, "Acute heroin overdose," Annals of Internal Medicine, vol. 130, no. 7, pp. 584-590,1999.

[5] E. White and C. Comiskey, "Heroin epidemics, treatment and ODE modelling," Mathematical Biosciences, vol. 208, no. 1, pp. 312-324, 2007.

[6] W. Wang, "Global behavior of an SEIRS epidemic model with time delays," Applied Mathematics Letters, vol. 15, no. 4, pp. 423428, 2002.

[7] P. van den Driessche and J. Watmough, "Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission," Mathematical Biosciences, vol. 180, pp. 29-48, 2002.

[8] T. Zhang and Z. Teng, "Global behavior and permanence of SIRS epidemic model with time delay," Nonlinear Analysis: Real World Applications, vol. 9, no. 4, pp. 1409-1424, 2008.

[9] D. R. Mackintosh and G. T. Stewart, "A mathematical model of a heroin epidemic: Implications for control policies," Journal of Epidemiology and Community Health, vol. 33, no. 4, pp. 299304, 1979.

[10] H. Su and X. Ding, "Dynamics of a nonstandard finite-difference scheme for Mackey-Glass system," Journal of Mathematical Analysis and Applications, vol. 344, no. 2, pp. 932-941, 2008.

[11] M. Sekiguchi and E. Ishiwata, "Global dynamics of a discretized SIRS epidemic model with time delay," Journal of Mathematical Analysis and Applications, vol. 371, no. 1, pp. 195-202, 2010.

[12] G. Mulone and B. Straughan, "A note on heroin epidemics," Mathematical Biosciences, vol. 218, no. 2, pp. 138-141, 2009.

[13] J. Liu and T. Zhang, "Global behaviour of a heroin epidemic model with distributed delays," Applied Mathematics Letters, vol. 24, no. 10, pp. 1685-1692, 2011.

[14] G. Huang and A. Liu, "A note on global stability for a heroin epidemic model with distributed delay," Applied Mathematics Letters, vol. 26, no. 7, pp. 687-691, 2013.

[15] G. P. Samanta, "Dynamic behaviour for a nonautonomous heroin epidemic model with time delay," Journal of Applied Mathematics and Computing, vol. 35, no. 1-2, pp. 161-178, 2011.

[16] X. Wang, J. Yang, and X. Li, "Dynamics of a Heroin epidemic model with very population," Applied Mathematics, vol. 2, no. 6, pp. 732-738, 2011.

[17] R. E. Mickens, "Discretizations of nonlinear differential equations using explicit nonstandard methods," Journal of Computational and Applied Mathematics, vol. 110, no. 1, pp. 181-185,1999.

Copyright of Abstract & Applied Analysis is the property of Hindawi Publishing Corporation and its content may not be copied or emailed to multiple sites or posted to a listserv without the copyright holder's express written permission. However, users may print, download, or email articles for individual use.