﻿ Stability Analysis of Chain, Mild and Passive Smoking Model

American Journal of Computational Mathematics
Vol.10 No.01(2020), Article ID:98185,12 pages
10.4236/ajcm.2020.101003

Stability Analysis of Chain, Mild and Passive Smoking Model

Abeer A. Alshareef, H. A. Batarfi*

Department of Mathematics, Faculty of Science, King Abdulaziz University, Jeddah, Saudi Arabia    Received: December 18, 2019; Accepted: February 7, 2020; Published: February 10, 2020

ABSTRACT

In this paper, the global stability of free smoking equilibrium point was evaluated and presented graphically. The linear stability of a developed mathematical model illustrates the effect on the population of chain, mild and passive smokers. MATLAB programming was used to simulate the solutions, the reproduction number ${R}_{0}$ and the nature of the equilibria.

Keywords:

Smoking Model, Stability Analysis, Global Stability and Lyapunov Method, Qualitative Behavior, Passive Smokers 1. Introduction

Castillo et al.  proposed a simple mathematical model for giving up smoking. They assumed that the total population consists of three populations: potential smokers, i.e. people who do not smoke yet and may become smokers in the future (P), smokers (S) and people who have already quit smoking permanently (Q). Since then, many articles were published further developing this model. The modifications were produced either by adding a factor or a new population (see e.g.    ). However, in 2008, Sharomi and Gumel  modified the model by including a class of smokers who temporarily quit smoking. They studied the global stability of smoker-free equilibrium (SFE). In 2011, Zaman  derived the giving-up smoking model taking into account the occasional smoking population and presented its qualitative behavior. In 2014, Z. Alkhudhari et al.  modified the model by assuming that the smoking classes are divided into occasional and chain smokers and they studied the existence and stability of equilibrium points. In 2017, Matintu  proposed a new modification by considering the moderate, chain smokers, potential smokers, and temporarily quit smokers. In 2019, Ana et al.  constructed a new smoking model where the total population was partitioned into three classes: passive smoker ${E}_{1}$ (those at risk from others’ smoking), infected I (people who are addicted to tobacco and can pass on the habit) and ${E}_{2}$ (people who have stopped smoking but are at risk of relapse) described as:

${\stackrel{˙}{E}}_{1}=\Delta -\frac{\beta }{N}{E}_{1}I-\left(\mu +{v}_{1}\right){E}_{1},$ (1a)

$\stackrel{˙}{I}=\frac{\beta }{N}{E}_{1}I-\left(\mu +ϵ+\theta \right)I,$ (1b)

${\stackrel{˙}{E}}_{2}=\theta I-\frac{\beta }{N}{E}_{2}I-\left(\mu +{v}_{2}+\omega \right)I,$ (1c)

where ${E}_{1}={E}_{1}\left(t\right)$,$I=I\left(t\right)$ and ${E}_{2}={E}_{2}\left(t\right)$. The parameters ( $\mu ,{v}_{1},{v}_{2},\omega ,\beta ,ϵ$ ) and $\theta$ represent respectively the rate of natural death, exit rate from ${E}_{1}$ to healthy population, exit rate form ${E}_{2}$ to healthy population, death rate of ex-smoker, infection rate of smoking, death rate of smokers and exit rate from I to ${E}_{2}$.$\Delta$ is the average number of healthy people who are at risk of becoming active smokers.

Our paper is divided into: section (2), evaluating the global stability of the model in . The modification of the model was stated in section (3) to include two classes of smokers: chain smokers ${S}_{1}$ and mild smokers ${S}_{2}$, with the full study of stability analysis and effect of each group on the population behavior (qualitative behavior).

2. The Global Stability of Equilibria

The equilibria of smoking-free equilibrium point in system (1a)-(1b), are given in  as

${P}_{1}=\left({E}_{1,0},{I}_{0},{E}_{2,0}\right)=\left(\frac{\Delta }{\mu +{v}_{2}},0,0\right),$ (2)

and the basic reproduction number ${\mathcal{R}}_{0}$ is defined as ${\mathcal{R}}_{0}=\frac{\beta }{\mu +ϵ+\theta }$.

Theorem 1.

If $\beta \le {\delta }_{4}$ and ${\mathcal{R}}_{0}\le 1$, then ${P}_{1}$ is globally asymptotically stable.

Proof.

The Lyapunov function is defined as

${V}_{1}={E}_{1,0}\left(\frac{{E}_{1}}{{E}_{1,0}}-1-\mathrm{ln}\frac{{E}_{1}}{{E}_{1,0}}\right)+{E}_{2}+I,$

hence,

${\stackrel{˙}{V}}_{1}=\left(1-\frac{{E}_{1,0}}{{E}_{1}}\right)+\stackrel{˙}{I}+{\stackrel{˙}{E}}_{2},$ (3)

By substituting (1a)-(1b) in (3),

$\begin{array}{c}{\stackrel{˙}{V}}_{1}=\left(1-\frac{{E}_{1,0}}{{E}_{1}}\right)\left(\Delta -\frac{\beta }{N}{E}_{1}I-\left(\mu +{v}_{1}\right){E}_{1}\right)+\frac{\beta }{N}I\left({E}_{1}+{E}_{2}\right)\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}-\left(\mu +ϵ+\theta \right)I+\theta I-\frac{\beta }{N}{E}_{2}I-\left(\mu +{v}_{2}+\omega \right){E}_{2},\end{array}$

$\begin{array}{c}{\stackrel{˙}{V}}_{1}=\left(1-\frac{{E}_{1,0}}{{E}_{1}}\right)\left(\Delta -\left(\mu +{v}_{1}\right){E}_{1}\right)-\frac{\beta }{N}{E}_{1}I+\frac{\beta }{N}{E}_{0,1}I+\frac{\beta }{N}I\left({E}_{1}+{E}_{2}\right)\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}-\left(\mu +ϵ+\theta \right)I+\theta I-\frac{\beta }{N}{E}_{2}I-\left(\mu +{v}_{2}+\omega \right){E}_{2},\end{array}$

Since ${E}_{1,0}=N$, then we have

${\stackrel{˙}{V}}_{1}=\left(1-\frac{{E}_{1,0}}{{E}_{1}}\right)\left(\Delta -{\delta }_{1}{E}_{1}\right)-{\delta }_{2}{E}_{2}+\left(\beta -\left(\mu +ϵ\right)\right)I.$

where ${\delta }_{1}=\mu +{v}_{1}$ and ${\delta }_{2}=\mu +{v}_{2}+\omega$. Using (2) we get

${\stackrel{˙}{V}}_{1}=-{\delta }_{1}\frac{{\left({E}_{1,0}-{E}_{1}\right)}^{2}}{{E}_{1}}-{\delta }_{2}{E}_{2}+\left(\beta -\left(\mu +ϵ\right)\right)I.$

By re-arranging the equation, we get

${\stackrel{˙}{V}}_{1}=-{\delta }_{1}\frac{{\left({E}_{1,0}-{E}_{1}\right)}^{2}}{{E}_{1}}-{\delta }_{2}{\mathcal{R}}_{0}{E}_{2}+\frac{{\delta }_{2}}{{\delta }_{3}}\left({\delta }_{3}\left({\mathcal{R}}_{0}-1\right)\right){E}_{2}+\left(\beta -{\delta }_{4}\right)I,$ (4)

where ${\delta }_{3}=\mu +ϵ+\theta$ and ${\delta }_{4}=\mu +ϵ$.

It is clear that ${\stackrel{˙}{V}}_{1}\le 0$ when $\beta \le {\delta }_{4}$ and ${\mathcal{R}}_{0}\le 1$ for all ${E}_{1},I,{E}_{2}>0$. Hence, the solutions of systems (1) is limited to $\Omega$, the largest invariant subset of ${\stackrel{˙}{V}}_{1}=0$, where $\Omega$ is the region of solutions given in . From (4), we see that ${\stackrel{˙}{V}}_{1}=0$ if and only if ${E}_{1}={E}_{1,0}$ and $I={E}_{2}=0$. By using LaSalle’s invariance principle, the point ${P}_{1}$ is globally asymptotically stable.

3. The Modified Model of Smoking

In this section, system (1) was modified to include two classes of smokers, chain smokers ${S}_{1}$ and mild smoker ${S}_{2}$ (see Figure 1). Accordingly, our model of ordinary differential equations (ODEs) are

Figure 1. Graphic presentation of system (6).

${\stackrel{˙}{E}}_{1}\left(t\right)=\text{Δ}-\frac{\beta }{N}{E}_{1}\left(t\right)\left(\eta {S}_{1}\left(t\right)+{S}_{2}\left(t\right)\right)-\left(\mu +{v}_{1}\right){E}_{1}\left(t\right),$ (5a)

${\stackrel{˙}{S}}_{1}\left(t\right)=a\frac{\beta }{N}\left(\eta {S}_{1}\left(t\right)+{S}_{2}\left(t\right)\right)\left({E}_{1}\left(t\right)+{E}_{2}\left(t\right)\right)-\left(\mu +{ϵ}_{1}+{\theta }_{1}\right){S}_{1}\left(t\right),$ (5b)

${\stackrel{˙}{S}}_{2}\left(t\right)=\left(1-a\right)\frac{\beta }{N}\left(\eta {S}_{1}\left(t\right)+{S}_{2}\left(t\right)\right)\left({E}_{1}\left(t\right)+{E}_{2}\left(t\right)\right)-\left(\mu +{ϵ}_{2}+{\theta }_{2}\right){S}_{2}\left(t\right),$ (5c)

${\stackrel{˙}{E}}_{2}\left(t\right)={\theta }_{1}{S}_{1}\left(t\right)+{\theta }_{2}{S}_{2}\left(t\right)-\frac{\beta }{N}{E}_{2}\left(t\right)\left(\eta {S}_{1}\left(t\right)+{S}_{2}\left(t\right)\right)-\left(\mu +{v}_{2}+\omega \right){E}_{2}\left(t\right),$ (5d)

where the parameters $\beta ,\omega ,{v}_{1},{v}_{2}$ and $\mu$ are defined as in the system (1a)-(1c), ${ϵ}_{1}$ and ${ϵ}_{2}$ are death rate of chain smokers ${S}_{1}$ and mild smokers ${S}_{2}$ respectively, while parameters ${\theta }_{1}$ and ${\theta }_{2}$ are exit rates from chain smokers and mild smokers to the healthy population (outside population N). We assume that the exposed people become either a chain or mild smoker at probabilities (1 − a) and a with $0 ). The chain smokers have a higher probability of generating more new smokers (by a factor $\eta \ge 1$ relative to the mild smokers). A simplified form of the model (5) can be written as

${\stackrel{˙}{E}}_{1}\left(t\right)=\text{Δ}-\frac{\beta }{N}{E}_{1}\left(t\right)\left(\eta {S}_{1}\left(t\right)+{S}_{2}\left(t\right)\right)-{d}_{1}{E}_{1}\left(t\right),$ (6a)

${\stackrel{˙}{S}}_{1}\left(t\right)=a\frac{\beta }{N}\left(\eta {S}_{1}\left(t\right)+{S}_{2}\left(t\right)\right)\left({E}_{1}\left(t\right)+{E}_{2}\left(t\right)\right)-{d}_{2}{S}_{1}\left(t\right),$ (6b)

${\stackrel{˙}{S}}_{2}\left(t\right)=\left(1-a\right)\frac{\beta }{N}\left(\eta {S}_{1}\left(t\right)+{S}_{2}\left(t\right)\right)\left({E}_{1}\left(t\right)+{E}_{2}\left(t\right)\right)-{d}_{3}{S}_{2}\left(t\right),$ (6c)

${\stackrel{˙}{E}}_{2}\left(t\right)={\theta }_{1}{S}_{1}\left(t\right)+{\theta }_{2}{S}_{2}\left(t\right)-\frac{\beta }{N}{E}_{2}\left(t\right)\left(\eta {S}_{1}\left(t\right)+{S}_{2}\left(t\right)\right)-{d}_{4}{E}_{2}\left(t\right),$ (6d)

where ${d}_{1}=\mu +{v}_{1}$,${d}_{2}=\left(\mu +{ϵ}_{1}+{\theta }_{1}\right)$,${d}_{3}=\left(\mu +{ϵ}_{2}+{\theta }_{2}\right)$ and ${d}_{4}=\left(\mu +{v}_{2}+\omega \right)$. The initial conditions of system (6), are given by

${E}_{1}\left(\alpha \right)={\phi }_{1}\left(\alpha \right)$,${S}_{1}\left(\alpha \right)={\phi }_{2}\left(\alpha \right)$,${S}_{2}\left(\alpha \right)={\phi }_{3}\left(\alpha \right)$,${E}_{2}\left(\alpha \right)={\phi }_{4}\left(\alpha \right)$

where $\alpha \in R$,$\left({\phi }_{1}\left(\alpha \right),{\phi }_{2}\left(\alpha \right),{\phi }_{3}\left(\alpha \right),{\phi }_{4}\left(\alpha \right)\right)\in C{R}_{+}^{4}$ and C is the Banach space of continuous functions mapping the interval $\left(-\infty ,0\right]$ into ${R}_{+}^{4}$. By the fundamental theory of functional differential equations , system (6) has a unique solution satisfying the initial conditions above.

3.1. Non-Negativity and Boundedness of Solutions

Theorem 2. For system (6), there exist a positive number K such that the compact set $\text{Λ}$,

$\text{Λ}=\left\{\left({E}_{1},{S}_{1},{S}_{2},{E}_{2}\right)\in {R}_{+}^{4}:0\le {E}_{1},{S}_{1},{S}_{2},{E}_{2}\le K\right\}.$

is positively invariant.

Proof. We have

${{\stackrel{˙}{E}}_{1}\left(t\right)|}_{{E}_{1}=0}=\text{Δ}>0$,

${{\stackrel{˙}{S}}_{1}\left(t\right)|}_{{S}_{1}=0}=a\frac{\beta }{N}{S}_{2}\left({E}_{1}+{E}_{2}\right)>0$ for  ${S}_{2},{E}_{1}$ and ${E}_{2}>0$,

${{\stackrel{˙}{S}}_{2}\left(t\right)|}_{{S}_{2}=0}=\left(1-a\right)\frac{\beta }{N}\eta {S}_{1}\left({E}_{1}+{E}_{2}\right)>0$ for  ${S}_{1},{E}_{1}$ and ${E}_{2}>0$,

${{\stackrel{˙}{E}}_{2}\left(t\right)|}_{{E}_{2}=0}={\theta }_{1}{S}_{1}+{\theta }_{2}{S}_{2}>0$ for  ${S}_{1}$ and   ${S}_{2}>0$,

so the solutions are non-negative. Now, to prove the boundedness of the solutions, we let $N\left(t\right)={E}_{1}\left(t\right)+{S}_{1}\left(t\right)+{S}_{2}\left(t\right)+{E}_{2}\left(t\right)$, hence

$\begin{array}{c}\stackrel{˙}{N}=\text{Δ}-\frac{\beta }{N}{E}_{1}\left(\eta {S}_{1}+{S}_{2}\right)-{d}_{1}\text{\hspace{0.17em}}{E}_{1}+a\frac{\beta }{N}\left(\eta {S}_{1}+{S}_{2}\right)\left({E}_{1}+{E}_{2}\right)-{d}_{2}{S}_{1}-{d}_{3}{S}_{2}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}+\left(1-a\right)\frac{\beta }{N}\left(\eta {S}_{1}+{S}_{2}\right)\left({E}_{1}+{E}_{2}\right)+{\theta }_{1}{S}_{1}+{\theta }_{2}{S}_{2}-\frac{\beta }{N}{E}_{2}\left(\eta {S}_{1}{S}_{2}\right)-{d}_{4}{E}_{2}\\ \le \text{Δ}-{d}_{1}{E}_{1}-{d}_{5}{S}_{1}-{d}_{6}{S}_{2}-{d}_{4}{E}_{2}\le \text{Δ}-\stackrel{¯}{d}N\left(t\right),\end{array}$

where ${d}_{5}=\mu +{ϵ}_{1}$,${d}_{6}=\mu +{ϵ}_{2}$ and $\stackrel{¯}{d}=\mathrm{min}\left\{{d}_{i}\right\}$,$i=1,4,5,6$. That is, if $N\left(0\right)\le 0$ then $N\left(t\right) where $K=\text{Δ}/\stackrel{¯}{d}$. It follows that if ${E}_{1}\left(0\right)+{S}_{1}\left(0\right)+{S}_{2}\left(0\right)+{E}_{2}\left(0\right)\le K$ then ${E}_{1}\left(t\right),{S}_{1}\left(t\right),{S}_{2}\left(t\right),{E}_{2}\left(t\right)\le K$.

3.2. The Equilibria and Basic Reproduction Number ${\mathcal{R}}_{0}$

Theorem 3. 1) If ${\mathcal{R}}_{0}\le 1$, then there exists only one equilibrium point called smoking-free equilibrium ${P}_{0}$.

2) If ${\mathcal{R}}_{0}>1$, there exist two equilibria that are smoking-free equilibrium ${P}_{0}$ and smoking present equilibrium ${P}_{1}$.

Proof.

The system (6) has two equilibrium points which are: Smoking-free equilibrium point

${P}_{0}=\left({E}_{1,0},{S}_{1,0},{S}_{2,0},{E}_{2,0}\right)=\left(\frac{\text{Δ}}{{d}_{1}},0,0,0\right),$ (7)

and smoking present equilibrium point

${P}_{1}=\left({E}_{1}^{\text{*}},{S}_{1}^{\text{*}},{S}_{2}^{\text{*}},{E}_{2}^{\text{*}}\right)$

where

${E}_{1}^{*}=\frac{{\mathcal{R}}_{0}\text{Δ}}{\left({\mathcal{R}}_{0}-1\right)\beta +{\mathcal{R}}_{0}\left(\mu +{v}_{1}\right)},$

${S}_{1}^{*}=a\frac{\beta \left({\mathcal{R}}_{0}-1\right)}{{\mathcal{R}}_{0}{d}_{2}}\left({x}_{1}+{E}_{2}^{*}\right),$

${S}_{2}^{*}=\left(1-a\right)\frac{\beta \left({\mathcal{R}}_{0}-1\right)}{{\mathcal{R}}_{0}{d}_{3}}\left({x}_{1}+{E}_{2}^{*}\right),$

${E}_{2}^{*}=\frac{\text{Δ}\beta \left({\mathcal{R}}_{0}-1\right)\left[\left(1-a\right){\theta }_{2}{d}_{2}+a{\theta }_{1}{d}_{3}\right]}{{d}_{2}{d}_{3}{x}_{2}\left[\beta \left({\mathcal{R}}_{0}-1\right)+{\mathcal{R}}_{0}{d}_{1}\right]},$

and

${x}_{1}=\left(\frac{\text{Δ}{\mathcal{R}}_{0}}{{\mathcal{R}}_{0}{d}_{1}+\beta \left({\mathcal{R}}_{0}-1\right)}\right),$

${x}_{2}={d}_{4}+\left(\frac{\beta \left({\mathcal{R}}_{0}-1\right)}{{\mathcal{R}}_{0}{d}_{2}{d}_{3}}\right)\left({d}_{2}{d}_{6}\left(1-a\right)+a{d}_{3}{d}_{5}\right).$

Using the next generation method , the matrices F and V are

$F=\left(\begin{array}{cccc}0& -\eta \beta & -\beta & 0\\ 0& a\eta \beta & a\beta & 0\\ 0& \left(1-a\right)\eta \beta & \left(1-a\right)\beta & 0\\ 0& 0& 0& 0\end{array}\right),\text{\hspace{0.17em}}V=\left(\begin{array}{cccc}{d}_{1}& 0& 0& 0\\ 0& {d}_{2}& 0& 0\\ 0& 0& {d}_{3}& 0\\ 0& 0& 0& {d}_{4}\end{array}\right),$

and hence evaluate the reproduction number ${\mathcal{R}}_{0}$ which is the largest eigenvalue of $F{\left(V\right)}^{-1}$,

${\mathcal{R}}_{0}=\beta \frac{\left(1-a\right){d}_{2}+a\eta {d}_{3}}{{d}_{2}{d}_{3}}.$ (8)

3.3. Local Stability of Smoker Free Equilibrium

Theorem 4.

If ${\mathcal{R}}_{0}<1$, then ${P}_{0}$ is locally asymptotically stable.

Proof. The Jacobian matrix of system (6) at ${P}_{0}$ is

$J\left({P}_{0}\right)=\left(\begin{array}{cccc}-{d}_{1}& -\eta \beta & \beta & 0\\ 0& \eta a\beta -{d}_{2}& a\beta & 0\\ 0& \left(1-a\right)\beta \eta & \left(1-a\right)\beta -{d}_{3}& 0\\ 0& {\theta }_{1}& {\theta }_{2}& -{d}_{4}\end{array}\right),$

with the eigenvalues ${\lambda }_{1}=-{d}_{1}$,${\lambda }_{2}=-{d}_{4}$, and ${\lambda }_{3,4}=c±b$, where

$c=\left(\frac{1}{2}\left[\beta \left(1-a\right)-{d}_{3}\right]+\left[\beta a\eta -{d}_{2}\right]\right),$

and

$b=\frac{1}{2}\sqrt{2\left({\mathcal{R}}_{0}-1\right)+{\left({d}_{2}-a\beta \eta \right)}^{2}+{\left({d}_{3}-\left(1-a\right)\beta \right)}^{2}+2a\eta \beta \left(\beta \left(1-a\right)\right)}.$

It is easy to show that $c<0$ when ${\mathcal{R}}_{0}<1$. To study the nature of b, we start with the fact that,

$\beta \left(\left(1-a\right){d}_{2}+a\eta {d}_{3}\right)<{d}_{2}{d}_{3},$

hence,

$\beta \left(1-a\right){d}_{2}<{d}_{2}{d}_{3}$, and $\beta a\eta {d}_{3}<{d}_{2}{d}_{3}$,

$\beta \left(1-a\right)<{d}_{3}$, and $\beta a\eta <{d}_{2}$.

Next, we use the numerical approach by assuming two small positive real numbers ${\xi }_{1}$ and ${\xi }_{2}$ such that $\beta \eta a+{\xi }_{1}={d}_{2}$ and $\beta \left(1-a\right)+{\xi }_{2}={d}_{3}$ hence

$\begin{array}{c}b=\frac{1}{2}\sqrt{2\left({\mathcal{R}}_{0}-1\right)+{\xi }_{1}^{2}+{\xi }_{2}^{2}+2\text{\hspace{0.17em}}\left({d}_{2}-{\xi }_{1}\right)\left({d}_{3}-{\xi }_{2}\right)}\\ =\frac{\sqrt{2}}{2}\sqrt{\left({\mathcal{R}}_{0}-1\right)+\frac{1}{2}{\xi }_{1}^{2}+\frac{1}{2}{\xi }_{2}^{2}+{d}_{2}{d}_{3}-{\xi }_{2}{d}_{2}-{\xi }_{1}{d}_{3}+{\xi }_{1}{\xi }_{2}}\\ =\frac{\sqrt{2}}{2}\sqrt{\left({\mathcal{R}}_{0}-1\right)+\frac{1}{2}{\left({\xi }_{1}+{\xi }_{2}\right)}^{2}+{d}_{2}{d}_{3}-{\xi }_{2}{d}_{2}-{\xi }_{1}{d}_{3}}.\end{array}$

Since $\frac{1}{2}{\left({\xi }_{1}+{\xi }_{2}\right)}^{2}\cong 0$ and ${d}_{2}{d}_{3}<{\xi }_{1}{d}_{2}+{\xi }_{2}{d}_{3}$ then ${d}_{2}{d}_{3}<{d}_{2}+{d}_{3}$ and

$\left({\mathcal{R}}_{0}-1\right)+\frac{1}{2}{\left({\xi }_{1}+{\xi }_{2}\right)}^{2}+{d}_{2}{d}_{3}-{\xi }_{2}{d}_{2}-{\xi }_{1}{d}_{3}<0.$

Thus ${\lambda }_{3,4}$ are complex numbers. Finally, we see that all the real parts of ${\lambda }_{i},i=1,2,3,4$ are negative so the point ${P}_{0}$ is locally asymptotically stable.

3.4. The Global Stability of Smoker Free Equilibrium

Theorem 5. If $\beta \eta \le {d}_{5},\beta \le {d}_{6}$ and ${\mathcal{R}}_{0}\le 1$, then ${P}_{0}$ is globally asymptotically stable.

Proof. As done in Theorem (1), we defined the Lyapunov function as

${V}_{2}={E}_{1,0}\left(\frac{{E}_{1}}{{E}_{1,0}}-1-\mathrm{ln}\frac{{E}_{1}}{{E}_{1,0}}\right)+{S}_{1}+{S}_{2}+{E}_{2},$

hence

${\stackrel{˙}{V}}_{2}=\left(1-\frac{{E}_{1,0}}{{E}_{1}}\right){\stackrel{˙}{E}}_{1}+{\stackrel{˙}{S}}_{1}+{\stackrel{˙}{S}}_{2}+{\stackrel{˙}{E}}_{2}.$ (9)

By substituting (6) in (9),

$\begin{array}{c}{\stackrel{˙}{V}}_{2}=\left(1-\frac{{E}_{1,0}}{{E}_{1}}\right)\left(\text{Δ}-\frac{\beta }{N}{E}_{1}\left(\eta {S}_{1}+{S}_{2}\right)-{d}_{1}{E}_{1}\right)+\frac{a\beta }{N}\left(\eta {S}_{1}+{S}_{2}\right)\left({E}_{1}+{E}_{2}\right)-{d}_{2}{S}_{1}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}+\frac{\left(1-a\right)\beta }{N}\left(\eta {S}_{1}+{S}_{2}\right)\left({E}_{1}+{E}_{2}\right)-{d}_{3}{S}_{2}+{\theta }_{1}{S}_{1}+{\theta }_{2}{S}_{2}-\frac{\beta }{N}{E}_{2}\left(\eta {S}_{1}+{S}_{2}\right)-{d}_{4}{E}_{2}\\ =\left(1-\frac{{E}_{1,0}}{{E}_{1}}\right)\left(\text{Δ}-{d}_{1}{E}_{1}\right)+\frac{\beta }{N}{E}_{1,0}\left(\eta {S}_{1}+{S}_{2}\right)-{d}_{2}{S}_{1}-{d}_{3}{S}_{2}+{\theta }_{1}{S}_{1}+{\theta }_{2}{S}_{2}-{d}_{4}{E}_{2}\end{array}$

By using the equilibrium conditions $\text{Δ}={d}_{1}{E}_{1,0}$ and ${E}_{1,0}=N$,

$\begin{array}{c}{\stackrel{˙}{V}}_{2}=-{d}_{1}\left(\frac{{\left({E}_{1}-{E}_{1,0}\right)}^{2}}{{E}_{1}}\right)+\left(\beta \eta -{d}_{5}\right){S}_{1}+\left(\beta -{d}_{6}\right){S}_{2}-{d}_{4}{E}_{2}\\ =-{d}_{1}\left(\frac{{\left({E}_{1}-{E}_{1,0}\right)}^{2}}{{E}_{1}}\right)+\left(\beta \eta -{d}_{5}\right){S}_{1}+\left(\beta -{d}_{6}\right){S}_{2}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}+\frac{{d}_{4}}{{c}_{1}}\left({c}_{1}\left[{\mathcal{R}}_{0}-1\right]\right){E}_{2}-{\mathcal{R}}_{0}{d}_{4}{E}_{2},\end{array}$ (10)

where ${c}_{1}={d}_{2}{d}_{3}$ and ${c}_{2}=\beta \left(1-a\right){d}_{2}+a\eta {d}_{3}$. Clearly ${\stackrel{˙}{V}}_{2}\le 0$, if $\beta \eta \le {d}_{5}$,$\beta \le {d}_{6}$ and ${\mathcal{R}}_{0}\le 1$, for all ${E}_{1},{S}_{1},{S}_{2},{E}_{2}>0$. Hence, the solutions of system (6) are limited to $\omega$, the largest invariant subset of ${\stackrel{˙}{V}}_{2}=0$. From (10) we see that ${\stackrel{˙}{V}}_{2}=0$ if and only if ${E}_{1}={E}_{1,0}$ and ${S}_{1}={S}_{2}={E}_{2}=0$. Using LaSalle’s invariance principle, we show that ${P}_{0}$ is globally asymptotically stable.

3.5. Numerical Simulations and Discussion

Using MATLAB programming, the solution of system (6) is evaluated by assuming values for the parameters as $\text{Δ}=50$,${ϵ}_{1}=0.03$,${ϵ}_{2}=0.02$,${\theta }_{1}=0.4$,${\theta }_{2}=0.5$,$\mu =0.05$,${v}_{1}=0.2$,${v}_{2}=0.3$,$N=200$, and $\omega =0.001$ with two sets of values for $\beta ,\eta$ and a.

Set (1): ${\mathcal{R}}_{0}\le 1$,$\beta =0.05$,$\eta =1.5$, and $a=0.4$.

Set (2): ${\mathcal{R}}_{0}>1$,$\beta =0.5$,$\eta =2.5$, and $a=0.6$.

A three sets of values are applied and tested for the initial conditions:

IC1: ${E}_{1}\left(0\right)=890,{S}_{1}\left(0\right)=120,{S}_{2}\left(0\right)=90,{E}_{2}\left(0\right)=40$.

IC2: ${E}_{1}\left(0\right)=480,{S}_{1}\left(0\right)=60,{S}_{2}\left(0\right)=40,{E}_{2}\left(0\right)=30$.

IC3: ${E}_{1}\left(0\right)=290,{S}_{1}\left(0\right)=30,{S}_{2}\left(0\right)=20,{E}_{2}\left(0\right)=10$.

In Figures 2-5, the behavior of each smoking population is plotted as time increases. Figure 2 shows the fast decrease of passive smokers for both ${\mathcal{R}}_{0}\le 1$ and ${\mathcal{R}}_{0}>1$ within a short period of time and reaches or approaches the steady

Figure 2. The plot of the passive smokers ( ${E}_{1}$ ), for set (1), ${\mathcal{R}}_{0}<1$ and set (2), ${\mathcal{R}}_{0}>1$ together with IC1 (dash), IC2 (solid) and IC3 (dash dot).

Figure 3. The plot of the chain smokers ( ${S}_{1}$ ), for set (1), ${\mathcal{R}}_{0}<1$ and set (2), ${\mathcal{R}}_{0}>1$ together with IC1 (dash), IC2 (solid) and IC3 (dash dot).

Figure 4. The plot of the mild smokers ( ${S}_{2}$ ), for set (1), ${\mathcal{R}}_{0}<1$ and set (2), ${\mathcal{R}}_{0}>1$ together with IC1 (dash), IC2 (solid) and IC3 (dash dot).

Figure 5. The plot of the ex-smokers ( ${E}_{2}$ ), for set (1), ${\mathcal{R}}_{0}<1$ and set (2), ${\mathcal{R}}_{0}>1$ together with IC1 (dash), IC2 (solid) and IC3 (dash dot).

state. Similar behavior was observed on the chain and mild populations for ${\mathcal{R}}_{0}\le 1$, Figure 3 and Figure 4 which reached zero as time increased. For ${\mathcal{R}}_{0}>1$, a peak appears at the beginning of the time and decreases to the steady state as time increases. This behavior was seen on the ex-smokers’ population for both ${\mathcal{R}}_{0}\le 1$ and ${\mathcal{R}}_{0}>1$. The choice of initial values of the three sets has no effect on the behaviour of the different population types. It tends to zero for the chain, mild and ex-smokers for ${\mathcal{R}}_{0}<1$ as time increases. ${E}_{2}$ behaves similarly for any initial conditions, where we observe that ${E}_{2}$ tends to zero ( ${\mathcal{R}}_{0}\le 1$ ) or steady state value ( ${\mathcal{R}}_{0}>1$ ). Similarly for other types of population which means that ${P}_{0}$ and ${P}_{1}$ are globally asymptomatically stable. Figure 6 shows the phase plane space plots that explain the relationships between the different populations. These relationships are the same for different initial conditions. Next, we consider the relationships of ex-smokers’ population versus the chain/mild smokers

Figure 6. The phase plane plot of the populations of system (6) using the values of the parameters of set (1) with the initial conditions ${E}_{1}\left(0\right)=30$,${S}_{1}\left(0\right)=290$,${S}_{2}\left(0\right)=20$ and ${E}_{2}\left(0\right)=10$.

${S}_{1}$ and ${S}_{2}$. As the chain smokers’ number decreases, we notice an increase on the population of ex-smokers ( ${E}_{2}$ ) up to a maximum value at around one third of chain smokers’ initial values and around two-thirds of the mild smokers’ initial value, Figure 5(c) and Figure 5(d). Finally, Figure 5(e) illustrates the impact of decreasing the number of chain smokers which leads to a decrease in the number of mild smokers.

As we observe from the Figures, we stress the importance of decreasing the number of chain smokers (either by using the media or educational direct seminars) on increasing the number of passive smokers (as the ex-smokers and mild are both joined) to create a healthy action.

Finally, we observe that when the contact between the two classes of smokers with the other populations is weak (small β) then the effect of the smoker’s classes (chain/mild) was not being enough to generate the new individuals to the smokers’ classes. So, the number of smoker’s groups decreased and approached to zero (because the death or stop smoking). On the other hand, if β increase, we found that the number of smoker’s groups increased also, since the rate of interaction between them and other groups (or the effect by smokers) was growing. However, we can see that, one of the solutions which is introduced to reduce the spread of smoking is prevent it in the larger places (such as educations) in order to minimize the rate of contact to be very small or zero.

4. Conclusion

The modified system/model shows interesting behavior between the different types of smoker’s populations. The main advantage was the more decrease in the size of heavy/mild smokers in comparison with earlier studies. This should encourage others (ex- or passive) smokers to interact better with chain/mild smokers to increase the number of quitters. By comparing our results with the results in , we see that they are approximately similar to giving the behavior of each group, and the difference is that, our model explains and shows the effect of both of the chain smokers and mild smokers separately on the exposed population. Hence, we get that the chain smokers affect more than the mild smokers on the exposed groups. So, our model is more real.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

Cite this paper

Alshareef, A.A. and Batarfi, H.A. (2020) Stability Analysis of Chain, Mild and Passive Smoking Model. American Journal of Computational Mathematics, 10, 31-42. https://doi.org/10.4236/ajcm.2020.101003

References

1. 1. Smith, A. (2008) Stroke of Genius for Gasoline Downsizing. Ricardo Q Rev, Q3.

2. 2. Sroka, Z.J. (2009) Micro Machines—The Concept of Engine Downsizing. Engine Technology International, 1, 58-62.

3. 3. Flierl, R., et al. (2011) Univalve—A Fully Variable Mechanical Valve Lift System for Future Internal Combustion Engines. Magazine ATZ, 4, 34-39.

4. 4. Merker, G.P., Schwartz, Ch. and Teichmann, R. (2009) Combustion Engines Development. Springer Edition.

5. 5. Bielaczyc, R. and Woodbrun, J. (2010) Global Trends in Emissions Regulation and Reduction. Combustion Engines Scientific Magazine, 142, 3-27.

6. 6. Shakun, J.D., et al. (2012) Global Warming Preceded by Increasing Carbon Dioxide Concentrations during the Last Deglaciation. Nature, 484, 49-54. http://dx.doi.org/10.1038/nature10915

7. 7. Regulation, E.C. (2009) nr 2009/33/WE of the European Parliament and of the Council—On the Promotion of Clean and Energy-Efficient Road Transport Vehicles. Brussels.

8. 8. Merkisz, J. and Radzimirski, S. (2007) Pollutant Emissions from Road Transportation in the Year 1970-2005 in Poland. Combustion Engines Scientific Magazine, 4-13.

9. 9. Tong, R.E. and Sanyol, A. (2012) Lower Energy Penalty CO2 Capture System. Carbon Capture Journal, 4.

10. 10. Kulazynski, M. and Sroka, Z.J. (2011) Developing Engine Technology. Print PAP, Lodz.

11. 11. Oppenheim, A.K., Kuhl, A.L., Packard, A.K., Hedrick, J.K. and Johnson, W.P. (1996) Model and Control of Heat Release in Engines. American Journal of Engine Combustion Flow Diagnosis, 1157, 15-23.

12. 12. Sroka, Z.J. (2012) Some Aspects of Thermal Load and Operating Indexes after Downsizing for Internal Combustion Engine. Journal of Thermal Analysis and Calorimetry, 110, 51-58. http://dx.doi.org/10.1007/s10973-011-2064-x

13. 13. Guzzella, L. and Onder, Ch.H. (2010) Introduction to Modelling and Control of Internal Combustion Engine Systems. 2nd Edition, Springer Editor. http://dx.doi.org/10.1007/978-3-642-10775-7

14. 14. Dudek, K. (1998) Thermography as a Diagnosing System. Journal of Transdisciplinary Systems Science, 3, 68-79.

15. 15. Heywood, J.B. (1989) Internal Combustion Engine Fundamentals. McGraw-Hill International Edition, Singapore.

16. 16. Sroka, Z.J. (2012) Thermal Load of Tuned Piston. Archives of Civil and Mechanical Engineering, 12, 342-347. http://dx.doi.org/10.1016/j.acme.2012.05.004

17. 17. Blair, G.P. (1996) Design and Simulation of Four-Stroke Engines. Edition SAE, Warrendale. http://dx.doi.org/10.4271/R-161

18. 18. Ramos, J.I. (1989) Internal Combustion Engine Modelling. American Publication Corporation, Hemisphere.

19. 19. Castillo, C., Jordan, G. and Rodriguez, A. (1997) Mathematical Models for the Dynamics of Tobacco Use, Recovery, and Relapse. Technical Report Series.

20. 20. Awan, A., Sharif, A., Hussain, T. and Ozair, M. (2017) Smoking Model with Cravings to Smoke. Advanced Studies in Biology, 9, 31-41. https://doi.org/10.12988/asb.2017.61245

21. 21. Alkhudhari, Z., Al-Sheikh, S. and Al-Tuwairqi, S. (2014) Global Dynamics of a Mathematical Model on Smoking. Applied Mathematics, 2014, Article ID: 847075. https://doi.org/10.1155/2014/847075

22. 22. Sallew, G. (2016) Mathematical Analysis on the Spread and Control of Global Dynamics of Tobacco Smoking with Induced Death. Journal of Mathematics and Computer, 3, 15-26.

23. 23. Lahrouz, A., Omari, L., Kiouach, D. and Belmaati, A. (2011) Deterministic and Stochastic Stability of a Mathematical Model of Smoking. Statistics and Probability Letters, 81, 1276-1284. https://doi.org/10.1016/j.spl.2011.03.029

24. 24. Sharomi, O. and Gumel, A. (2008) Curtailing Smoking Dynamics: A Mathematical Modeling Approach. Applied Mathematics and Computation, 195, 475-499. https://doi.org/10.1016/j.amc.2007.05.012

25. 25. Zaman, G. (2011) Qualitative Behavior of Giving up Smoking Models. Bulletin of the Malaysian Mathematical Sciences Society, 34, 403-415.

26. 26. Alkhudhari, Z., Al-Sheikh, S. and Al-Tuwairqi, S. (2014) The Effect of Occasional Smokers on the Dynamics of a Smoking Model. International Mathematical Forum, 9, 1207-1222. https://doi.org/10.12988/imf.2014.46120

27. 27. Matintu, S.A. (2017) Smoking as Epidemic: Modeling and Simulation Study. American Journal of Applied Mathematics, 5, 31-38. https://doi.org/10.11648/j.ajam.20170501.14

28. 28. Pulecio-Montoya, A., Lopez-Monteneqro, L. and Benavides, L. (2019) Analysis of a Mathematical Model of Smoking. Contemporary Engineering Sciences, 12, 117-129. https://doi.org/10.12988/ces.2019.9517

29. 29. Hale, J. and Lunel, S. (1993) Introduction to Functional Differential Equations. Science and Business Media. https://doi.org/10.1007/978-1-4612-4342-7

30. 30. Heffernan, J., Smith, R. and Wahl, L. (2005) Perspectives on the Basic Reproduction Ratio. Journal of the Royal Society Interface, 2, 281-293. https://doi.org/10.1098/rsif.2005.0042