Approximating a Giving Up Smoking Dynamic on Adolescent Nicotine Dependence in Fractional Order.

In this work, we consider giving up smoking dynamic on adolescent nicotine dependence. First, we use the Caputo derivative to develop the model in fractional order. Then we apply two different numerical methods to compute accurate approximate solutions of this new model in fractional order and compare their results. In order to do this, we consider the generalized Euler method (GEM) and multi-step generalized differential transform method (MSGDTM). We also show the unique positive solution for this model and present numerical results graphically.

Introduction

Mathematical modeling is a tool which have a lot of applications in applied sciences that require a deep study for the different verities of methods used in applied mathematics. Smoking is the prodigious cause of many diseases specially of different type of cancers. As cigarette smoke contains over forty thousands dangerous chemicals, which cause harmful infections to human health. The life of smokers is ten to thirteen years shorter than that of non smokers and at the rat of 1:2 of smokers die from diseases launched by the cigarette smoking. Besides, smokers have 70% more chances of heart attack than non smokers. The incidence rate of lung cancer is 10% high in smokers. To secure the life expectancy every scientist, doctor and mathematician try to control smoking. Mathematicians try to make different smoking models for the best representation of cigarette smoking phenomena. So several authors proposed different smoking models, for example see the first model presented by Castillo-Garsow et al., [1] in which they studied different classes of smokers (potential smokers (P), smokers (S), and quit smokers (Q)). Then Sharami et al., [2] modified the model presented in [1] and introduced a new class named chain smokers. In their work, they presented the development and public health impact of smoking related diseases. Zeb et al., [3] presented the square-root dynamics of a giving up smoking model. In this work they discussed local and global stability of the model and its general solutions in which the interaction of occasional and potential smokers occur. Nowadays, every researcher tries to establish different epidemic models in fractional order. Perhaps, the reality of nature could be better translated by fractional calculus. Fractional calculus is used in many fields of sciences [4-15]. In this paper, we present the fractional order derivative and find analytic numeric solution of model presented in [16], and the model is as follow: under the initial conditions: where P, S1, S2, Q and Q all functions of t denote the numbers of potential, occasional and quit smokers. Here Λ is number of incoming population per unit time, b is incoming rate per unit time i = 1, 2, 3, 4, 5 : Average rate of incoming for a year, d: Outgoing per unit time i = 1, 2, 3, 4, 5, μ is natural death rate or death related to other diseases not smoking and d represents death rate of all classes caused by smoking diseases.

By establish fractional order derivative into the model [16], we obtain the following fractional order model: Here we consider the Caputo sense fractional derivative and α is the order of the fractional time-derivative, subject to the initial conditions given in Eq (2). In the present paper, we take in account the model presented in [16]. First, we will use the Caputo derivative to develop the model in fractional order. Then we will apply two different numerical methods to compute accurate approximate solutions of this new model in fractional order and compare their results. In order to do this, we consider generalized Euler method (GEM) and the multi-step generalized differential transform method (MSGDTM). We also show the unique positive solution for this model and present numerical results graphically.

This paper is organized as follows. In Section 2, we present positivity of the new model. The idea of GEM and MSGDTM for solution of the proposed model are presented shortly in Sections 3. Numerical simulation results graphically of GEM and MSGDTM are presented in Section 4 with comparisons with the results of Runge-Kutta method (RKM). The conclusion is given in Section 5. Appendix is devoted to present some basic definitions and results which is needed in this work are given in the supporting file.

1 Non-negative Solutions

Let and X = (P, S1, S2, Q, Q).

We need lemma [Generalized Mean Value Theorem] [17] that help us in the proof of subsequent theorem.

Theorem 1 If P > 0, S1 > 0, S2 > 0, Q > 0, Q > 0, all at (0) the solutions of P, S1, S2, Q, Q at (t) of the system Eq (3), are positive for all t > 0.

Proof 1 If the above condition does not satisfied, then at least one of the individuals may be negative. Then the individuals will satisfy one of the following conditions.

There exists a first time t1 But , which is contradiction to the above supposition.

There exists a first time t2 But , which is contradiction to the above supposition.

There exists a first time t3 But , which is contradiction to the above supposition.

There exists a first time t4 But , which is contradiction to the above supposition.

There exists a first time t5 But , which is contradiction to the above supposition.

2 Generalized Euler Method (GEM) and Multistep Generalized Differential Transform Method (MSGDTM)

Generally, it is impossible to determine the analytic solution of nonlinear differential equations. However, different authors used different methods for the approximate analytic and numeric solutions of nonlinear differential equations [7-11]. But all these methods are applicable for short time interval. So we have in mind to find the approximate solutions of problems Eqs (2) and (3) via the GEM and MSGDTM. The GEM is derived to get numerical solution of initial value problems with Caputo derivatives by Odibat et al., [11]. This method is a generalization of the classical Euler’s method for detail see [11]. Also in this section, we use the MSGDTM that to find numerical solution of the system Eq (3) of fractional order differential equations and the detail analysis of this method is found in [12]. The basics steps of the GDTM can be found in [12-15]. In this paper, we are not going in detail of these methods but we use them only for the numerical solution of the system Eq (3). Applying the MSGDTM, we obtain the following system:

Here, P(k), S1(k), S2(k), Q(k) and Q(k) are the differential transformation of P(t), S1(t), S2(t), Q(t) and Q(t), respectively. Also, PS1 S2 N(k) is the differential transformation of function and defined as follows: In next section, we will present the numerical results.

3 Numerical and Simulation Results

We applied the GEM and MSGDTM to solve the system Eq (1) for α = 1. In order to demonstrate the effectiveness of these two methods as an approximate tool for solving the nonlinear system of fractional differential Eq (3) for larger time t, we apply these two methods on the interval [0, 30]. It is to be noted that GEM results are obtained for h = 0.00001 while the MSGDTM results are obtained when K = 10 and M = 30000 and RKM for h = 0.001. All the results are calculated by using computer algebra package Mathematica and in this paper, we show only the graphically obtained results.

We assume the parameters of the system Eq (1) shown in Table 1. Figs 1–5 show the approximate solutions for P(t), S1(t), S2(t), Q(t) and Q(t) obtained for different values of α using GEM. From the graphical results given in Figs 1–5, it can be seen that the results obtained using GEM match the results of Runge-Kutta for integer case(α = 1,) very well when α = 1, which implies that GEM can predict the behavior of these variables accurately for the region under consideration.

Table 1

The parameters of system Eq (1).

Parameter	Range	Value
Λ	2 per week	2
u1	u1 ∈ (0, ∞)	0.047
bi(i = 1, …, 5)	bi ∈ (0,1)	0.2
u2	u2 ∈ (0, ∞)	0.023
di(i = 1, …, 5)	di ∈ [0, 1)	0
γ1	γ1 ∈ (0, ∞)	0.025
γ2	γ1 ∈ (0, ∞)	0.012
γ3	γ1 ∈ (0, ∞)	0.074
μ	520 weeks	0.0019
β2	β2 ∈ (0, ∞)	0.012
∝	∝ ∈ (0, 1)	0.05
Ψ	Ψ ≥ 1	2

Fig 1

The population of P(t) versus t: α = 1.0 (solid line), α = 0.95 (dashed line) and α = 0.85 (dot-dashed line).

Fig 5

The population of Q(t) versus t: α = 1.0 (solid line), α = 0.95 (dashed line) and α = 0.85 (dot-dashed line).

Discussion

In this paper, we used Generalized Euler method (GEM) and multi-step generalized differential transform method (MSGDTM) as a reasonable basis for studying the dynamics of a new giving up smoking model and modified the integer order model Eq (1) into a fractional-order model Eq (3). The results obtained show that when α → 1 the solutions of fractional model, {P}(t), {S1}(t), {S2}(t), {Q}(t) and {Q}(t), reduce to the standard solutions P(t), S1(t), S2(t), Q(t) and Q(t). Finally, the recent appearance of fractional differential equations as models in some fields of Science and Engineering makes it necessary to investigate analytical and numerical methods of solution for such equations. The present success of the proposed methods for the considered model verifies that it is a useful tool for these kind of models in Science and Engineering.

Appendix

Basics of Fractional Calculus

In this appendix, we give some basic definitions and properties of the fractional calculus theory which are used further in [4–6, 17, 18].

Definition 1 A real function f(x)(x > 0) is said to be in the space C(α ∈ R) iff it can be written as f(x) = x f1(x) for some p > α where f1(x) is continuous in [0, ∞), and it is said to be in the space if f( ∈ C, m ∈ N.

Definition 2 The Riemann–Liouville integral operator of order α > 0 with a ≥ 0 is defined as

Properties of this operator can be found in [4]. We only need here the following:

For f ∈ C, α, β > 0, a ≥ 0, c ∈ R and γ > −1, we have where B(α, γ+1) is the incomplete beta function which is defined as The Riemann-Liouville derivative has certain disadvantages when trying to model real-world phenomena with fractional differential equations. For example, the Riemann-Liouville derivative of a constant is not zero. In addition, if an arbitrary function is a constant at the origin, its fractional derivation has a singularity at the origin for instant exponential and Mittag-Leffler functions. Theses disadvantages reduce the field of application of the Riemann-Liouville fractional derivative. One of the great advantages of the Caputo fractional derivative is that it allows traditional initial and boundary conditions to be included in the formulation of the problem [18]. Therefore, we shall introduce a modified fractional differential operator proposed by Caputo in his work on the theory of viscoelasticity.

Definition 3 The Caputo fractional derivative of f(x) of order α > 0 with a ≥ 0 is defined as

The Caputo fractional derivative was investigated by many authors,

Theorem 2 (Generalized Taylor’s formula) [17] Suppose that for k = 0, 1, …, n + 1, where 0 < α ≤ 1. Then we have with a ≤ ξ ≤ x, ∀x ∈ (a, b], where .

For mathematical properties of fractional derivatives and integrals one can consult the mentioned references in this paper.