Processing math: 67%
Research article Special Issues

A mathematical model of tuberculosis (TB) transmission with children and adults groups: A fractional model

  • Received: 01 December 2019 Accepted: 19 February 2020 Published: 17 March 2020
  • MSC : 34A34

  • We consider a novel fractional model to investigates the (tuberculosis) TB model dynamics with two age groups of human, that is, the children and the adults. First, we formulate the model and present the basic results associated to the model. Then, using the fractional operators, Caputo and the Atangana-Baleanu and obtain a generalized model. Further, we give a novel numerical approach for the solution of the fractional model and obtain their approximate solution. We show graphical results with various values of the fractional order. A comparison of the two operators are shown graphically. The results obtained through Atangana-Baleanu operator is flexible than that of Caputo derivative. The infection in tuberculosis (TB) infected people decreases fast when decreasing the fractional order.

    Citation: Fatmawati, Muhammad Altaf Khan, Ebenezer Bonyah, Zakia Hammouch, Endrik Mifta Shaiful. A mathematical model of tuberculosis (TB) transmission with children and adults groups: A fractional model[J]. AIMS Mathematics, 2020, 5(4): 2813-2842. doi: 10.3934/math.2020181

    Related Papers:

    [1] Jing Li, Linlin Dai, Kamran, Waqas Nazeer . Numerical solution of multi-term time fractional wave diffusion equation using transform based local meshless method and quadrature. AIMS Mathematics, 2020, 5(6): 5813-5838. doi: 10.3934/math.2020373
    [2] Farman Ali Shah, Kamran, Zareen A Khan, Fatima Azmi, Nabil Mlaiki . A hybrid collocation method for the approximation of 2D time fractional diffusion-wave equation. AIMS Mathematics, 2024, 9(10): 27122-27149. doi: 10.3934/math.20241319
    [3] Fouad Mohammad Salama, Nur Nadiah Abd Hamid, Norhashidah Hj. Mohd Ali, Umair Ali . An efficient modified hybrid explicit group iterative method for the time-fractional diffusion equation in two space dimensions. AIMS Mathematics, 2022, 7(2): 2370-2392. doi: 10.3934/math.2022134
    [4] Emad Salah, Ahmad Qazza, Rania Saadeh, Ahmad El-Ajou . A hybrid analytical technique for solving multi-dimensional time-fractional Navier-Stokes system. AIMS Mathematics, 2023, 8(1): 1713-1736. doi: 10.3934/math.2023088
    [5] Mubashir Qayyum, Efaza Ahmad, Hijaz Ahmad, Bandar Almohsen . New solutions of time-space fractional coupled Schrödinger systems. AIMS Mathematics, 2023, 8(11): 27033-27051. doi: 10.3934/math.20231383
    [6] Amjad Ali, Iyad Suwan, Thabet Abdeljawad, Abdullah . Numerical simulation of time partial fractional diffusion model by Laplace transform. AIMS Mathematics, 2022, 7(2): 2878-2890. doi: 10.3934/math.2022159
    [7] Xiangmei Li, Kamran, Absar Ul Haq, Xiujun Zhang . Numerical solution of the linear time fractional Klein-Gordon equation using transform based localized RBF method and quadrature. AIMS Mathematics, 2020, 5(5): 5287-5308. doi: 10.3934/math.2020339
    [8] Yudhveer Singh, Devendra Kumar, Kanak Modi, Vinod Gill . A new approach to solve Cattaneo-Hristov diffusion model and fractional diffusion equations with Hilfer-Prabhakar derivative. AIMS Mathematics, 2020, 5(2): 843-855. doi: 10.3934/math.2020057
    [9] Tareq Eriqat, Rania Saadeh, Ahmad El-Ajou, Ahmad Qazza, Moa'ath N. Oqielat, Ahmad Ghazal . A new analytical algorithm for uncertain fractional differential equations in the fuzzy conformable sense. AIMS Mathematics, 2024, 9(4): 9641-9681. doi: 10.3934/math.2024472
    [10] Lynda Taleb, Rabah Gherdaoui . Approximation by the heat kernel of the solution to the transport-diffusion equation with the time-dependent diffusion coefficient. AIMS Mathematics, 2025, 10(2): 2392-2412. doi: 10.3934/math.2025111
  • We consider a novel fractional model to investigates the (tuberculosis) TB model dynamics with two age groups of human, that is, the children and the adults. First, we formulate the model and present the basic results associated to the model. Then, using the fractional operators, Caputo and the Atangana-Baleanu and obtain a generalized model. Further, we give a novel numerical approach for the solution of the fractional model and obtain their approximate solution. We show graphical results with various values of the fractional order. A comparison of the two operators are shown graphically. The results obtained through Atangana-Baleanu operator is flexible than that of Caputo derivative. The infection in tuberculosis (TB) infected people decreases fast when decreasing the fractional order.


    Many natural phenomena can be modeled mathematically to obtain approximate models[1]. Compared to the classical diffusion equation, the fractional diffusion equation may be more suitable for modelling anomalously slow transport processes with memory and inheritance. In recent years, fractional calculus has found widespread applications in many fields including turbulence, wave propagation, signal processing, porous media, and anomalous diffusion[2,3,4].

    Considering that the single-term time fractional derivative cannot adequately describe many complex physical or biological processes, recently, a multi-term time and time distributed order fractional equations have been developed. The time distribution order equation is also a generalisation of the multi-term time equation. Therefore, the study of multi-term time fractional differential equations becomes very important and useful in various applications[5,6]. For example, multi-term fractional diffusion equation has simplified the modelling of phenomena such as diffusion processes, viscoelastic damping materials, oxygen delivery through capillaries and anomalous relaxation of magnetic resonance imaging signal magnitude[7,8,9]. Because the numerical solution is the most important in practice, a great deal of research has been done in the study of numerical solutions of multi-term fractional diffusion equations[10,11,12]. Qiu has analyzed numerical solutions for the Volterra integrodifferential equations with tempered multi-term kernels[13]. Hu et al. have formulated a backward Euler difference scheme for the integro-differential equations with the multi-term kernels[14]. Guo et al. have proposed the alternating direction implicit numerical approaches for computing the solution of multi-dimensional distributed order fractional integrodifferential problems[15]. {Guo et al. have developed an efficient finite difference/generalized Hermite spectral method for the distributed-order time-fractional reaction-diffusion equation on one-, two-, and three-dimensional unbounded domains[16].

    In general, exact solutions of the fractional diffusion equation are rarely obtained in practical applications. Therefore, it is necessary to develop some effective numerical methods to solve the multi-term time-space fractional equations. A large number of numerical methods have been developed for two-dimensional time-space fractional order diffusion equations with a single time fractional order derivative as a special case of multi-term time-space fractional equations. Abd-Elhameed et al. have introduced a new set of orthogonal polynomials to effectively obtain numerical solutions of the nonlinear time-fractional generalized Kawahara equation by using the collocation algorithm[17]. Moustafa et al. have created Chebyshev polynomials for the time fractional fourth-order Euler-Bernoulli pinned-pinned beam based on the Petrov-Galerkin[18]. Peng et al. have developed a novel fourth-order compact difference scheme for the mixed-type time-fractional Burgers equation, using the L1 discretization formula and a nonlinear compact difference operator[19]. Marasi and Derakhshan have proposed a hybrid numerical method based on the weighted finite difference and the quintic Hermite collocation methods for the for solving the variable-order time fractional mobile-immobile advection-dispersion model[20].

    It is well known that meshless methods are a type of point set based numerical method that considers a set of scattered and uniform data points. Due to this property, the meshless method can be applied to high-dimensional models with irregular and complex domains[21,22,23]. However, due to the singularity of the spatial fractional operators, we only deal with problems on convex domains. A meshless method with Hermite splines of order five is used to discretize the Riesz fractional operator in the spatial direction, which gives higher accuracy with fewer points.

    From the last few decades, there are many methods to solve the single time fractional diffusion equations, for instance the finite difference method, interpolation, implicit stepping methods, etc. The Laplace transform is one of the powerful tools for solving differential equations in engineering and other scientific disciplines. However, solving differential equations with the Laplace transform sometimes results in solutions in the Laplace domain that are not easily invertible to the real domain by analytical methods. Therefore, we use numerical inversion methods to transform the obtained solutions from the Laplace domain to the real domain[24,25]. The Laplace transform overcomes the memory effect arising from the convolution integral expressions of the time fractional derivative term, and better results can be obtained in the case of more general smoothness.

    The remaining sections of our paper are organised as follows. Some important preliminary, definitions and lemmas are given in Section 2. The introduction of the model and the time discretization based on the Laplace transform are given in Section 3. The simplification and approximation theory of the equation is given in Section 4. The basis construction, the meshless method for solving the simplified space fractional equation and convergence analysis of the quintic Hermite spline are presented in Section 5. Meanwhile, numerical examples are given in Section 6. Section 7 explains the analysis and the results of the research. Finally, a brief conclusion is given in Section 8.

    In this section we will introduce some concepts and properties. Let Ω satisfy segment conditions of the form [26], let Υ be a rectangular domain containing Ω, let the symbol |Ω stand for restriction to Ω, and let Υ=[a,b]×[c,d]Ω.

    Definition 2.1. The left and right Caputo fractional derivational of order α on [a,b] is defined by

    CxDαLf(x)=1Γ(nα)xa(xη)nα1nf(η)ηndη,CxDαRf(x)=1Γ(nα)bx(tη)nα1nf(η)ηndη,

    where n1<αn, and n=α.

    Definition 2.2. The left and right Riemann-Liouvile fractional derivatives operator with respect to order α on [a,b] is defined by

    RLxDαLf(x)=1Γ(nα)nxnxa(xη)nα1f(η)dη,RLxDαRf(x)=(1)nΓ(nα)nxnbx(ηx)nα1f(η)dη,

    where n1<αn, and n=α.

    Definition 2.3. Let α>0,m=α, the connection between Riemann-Liouville derivatives and Caputo derivatives is

    CxDαLf(t)=RLxDαLf(t)m1k=0f(k)(a)Γ(kα+1)(ta)kα,CxDαRf(t)=RLxDαRf(t)m1k=0f(k)(b)Γ(kα+1)(bt)kα.

    Definition 2.4. The Riesz fractional derivative with order α>0 on a finite interval [a,b] is defined by

    α|x|αf(x)=cα(RLxDαLf(x)+RLxDαRf(x)),

    where

    cα=12cos(απ/2),

    α2k+1,k=0,1,, and for n1αn,nN.

    Definition 2.5. For given v: [0,]R, the definition of Laplace transform is

    L{v(t)}=V(s)=0estv(t)dt.

    Lemma 2.1. ([27,Lemma 1.2.]) Suppose that v(t)Cp[0,), the Laplace transform of Caputo fractional derivative about v(t) is

    L{αtv(t)}=sαV(s)p1i=0sαi1v(i)(0),p1<α<pZ+.

    We introduce some definitions of spaces. Let Ω be a domain in Rn,

    C[a,b]={u(x)|u(x)is a continuous function on[a,b]},Cm[a,b]={u(x)|u(m)(x)is a continuous function on[a,b]},uC2(Ω)=max{u(k,l)C(Ω),k,lN,k+l2},S5,2(π)={ϕC2[a,b]:ϕ|kjP5,j=1,2,3,,n},

    where P5 is the set of polynomial functions with order not greater than 5 over kj.

    Definition 2.6. For any nonnegative integer m let Cm(Ω) denote the vector space consisting of all functions f which, together with all their partial derivatives Dαf of orders αm, are continuous on Ω. We abbreviate C0(Ω)Ω. Let

    C=m=0(Ω).

    The subspaces C0(Ω) and C0(Ω) consist of all those functions in C(Ω) and C(Ω), respectively, that have compact support in Ω.

    Definition 2.7. Give a positive integer τ and a real number r(1r<). The Sobolev space Wτ,r is defined by

    Wτ,r(Ω)={uLr(Ω),

    the weak derivative DθuLr(Ω) for 0|θ|τ}, with norm

    uWτ,r=(0|θ|τDθurLr(Ω))1r,

    where θ=(θ1,θ2), |θ|=θ1+θ2 and θ1, θ2 are non-negative integers.

    Definition 2.8. ([26]) For yΩ, there exists a nonzero vector oy and a neighborhood Uy such that if zˉΩUy, then z+toyΩ for 0<t<1, and call that Ω satisfies the segment condition.

    Lemma 2.2. ([26]) If Ω satisfies the condition of Definition 2.8, then the set of restrictions to Ω of functions in C0(R2) is dense in Wτ,r(Ω).

    We consider the multi-term time-space fractional diffusion equations of the following form

    rq=0aq(C0Dαqt)u(x,y,t)=Δβ,γu(x,y,t)+f(x,y,t),(x,y)Ω,0<tT, (3.1)

    subject to

    u(x,y,0)=ϕ(x,y),(x,y)Ω,u(x,y,t)|Ω=0,t(0,T], (3.2)

    where 0<αq<<α1<α0<1 is the time fractional order, C0Dαqt is the Caputo fractional derivative with αq order given by Definition 2.1,

    rq=0aq=1,q=0,1,,r

    are the coefficients, and ΩR2 is bounded convex domain.

    The spatial operator Δβ,γ, 1<β,γ<2 is a Riesz fractional order operator given by Definition 2.4,

    Δβ,γu(x,y,t):=Kxβu(x,y,t)|x|β+Kyγu(x,y,t)|y|γ =Kxcβ(RLxDβLu(x,y,t)+RLxDβRu(x,y,t))+Kycγ(RLyDγLu(x,y,t)+RLyDγRu(x,y,t)),

    where the constants Kx>0,Ky>0 are diffusion coefficients. And the left side and right side Riemann-Liouville derivatives on x,y direction, respectively, are defined by Definition 2.2,

    RLxDβLu(x,y,t)=1Γ(2β)2x2xa(xv)1βu(v,y,t)dv,RLxDβRu(x,y,t)=1Γ(2β)2x2bx(vx)1βu(v,y,t)dv,RLyDγLu(x,y,t)=1Γ(2γ)2y2yc(yv)1γu(x,v,t)dv,RLyDγRu(x,y,t)=1Γ(2γ)2y2dy(vy)1γu(x,v,t)dv.

    Suppose that u(x,y,t)C1(Ω), using the Laplace transform on Eq (3.1) and owing to the property of Lemma 2.1, we have

    Lu(x,y,t)=U(x,y,s),Lf(x,y,t)=F(x,y,s),LC0Dαqtu(x,y,t)=sαqU(x,y,s)sαq1u(x,y,0)=sαqU(x,y,s)sαq1ϕ(x,y).

    So this equation could be

    rq=0aq(sαqU(x,y,s)sαq1ϕ(x,y))=Δβ,γU(x,y,s)+F(x,y,s),(x,y)Ω. (3.3)

    Equation (3.2) becomes with the boundary conditions

    U(x,y,s)|Ω=0,sC. (3.4)

    Then the methods for the inverse Laplace transform methods are based on numerical integration of the Bromwich complex contour integral. From [27,29], using the strategy of Talbot, the Bromwich line can be transformed into a contour that starts and ends in the left half plane,

    u(x,y,t)=L1U(x,y,s)=12πiσ+iσiestU(x,y,s)ds,σ>σ0,

    where σ0 is the convergence abscissa. Two simpler types of contours have mainly been proposed mainly:

    ● Parabolic path: s=μ(1+iz)2, z=γ+ic, where c>0,<γ<, then,

    s(γ)=μ((1c)2γ2)+2iμγ(1c).

    ● Hyperbolic path:

    s(γ)=ω+λ(1sin(δiγ))

    for <γ<.

    On either of the above contours, the Bromwich integral becomes

    u(x,y,t)=12πies(z)tU(x,y,s(z))s(z)dsι2πiLl=LezltU(x,y,s(zl))s(zl),zl=lι. (3.5)

    Let operator Δβ,γu(x,y,t) on Υ=[a,b]×[c,d], using the Laplace transform, it becomes as follows:

    Δβ,γU(x,y,s)=Kxcβ{CxDβLU(x,y,s)+U(a,y,s)Γ(2β)(xa)1β+CxDβRU(x,y,s)+U(b,y,s)Γ(2β)(bx)1β}+Kycγ{CyDγLU(x,y,s)+U(x,c,s)Γ(2γ)(yc)1γ+CyDγRU(x,y,s)+U(x,d,s)Γ(2γ)(dy)1γ}.

    To avoid the singularity of the operator Δβ,γ, let

    A(x,y)=(xa)β1(bx)β1(yc)γ1(dy)γ1,

    for sαqlC, the Eq (3.3) becomes

    A(x,y)(rq=0aqsαqlU(x,y,sl)Δβ,γU(x,y,sl))=A(x,y)(rq=0aq(sαq1lϕ(x,y))+F(x,y,sl)),

    then expand the sl and denote

    G=(G1,G2)T,W=(W1,W2)T,
    {G1U(x,y,sl)A(x,y)(rq=0aq(Re(sαql)Re(U(x,y,sl))Im(sαql)Im(U(x,y,sl)))Δβ,γRe(U(x,y,sl)))=A(x,y)(rq=0aqRe(sαq1l)ϕ(x,y)+Re(F(x,y,sl)))W1(x,y,sl),G2U(x,y,sl)A(x,y)(rq=0aq(Re(sαql)Im(U(x,y,sl))+Im(sαql)Re(U(x,y,sl)))Δβ,γIm(U(x,y,sl)))=A(x,y)(rq=0aqIm(sαq1l)ϕ(x,y)+Im(F(x,y,sl)))W2(x,y,sl),

    where Re(U(x,y,sl)) stands the real part of U, Im(U(x,y,sl)) stands the imaginary part of U. Then the Eq (3.3) becomes

    GU(x,y,sl)=W(x,y,sl). (4.1)

    Meanwhile, the Eq (3.4) becomes

    {Re(U(x,y,sl))|Ω=0,Im(U(x,y,sl))|Ω=0. (4.2)

    Let Ωr be rectangular domains containing Ω. Denote by S the set of 2-dimension polynomial functions.

    Lemma 4.1. ([28,Theorem 2.2]) C0(R2)|Ω is dense in C2(Ω).

    Lemma 4.2. ([6,Lemma 2.1]) S is dense in C(Ωr) with the norm C2(Ωr).

    If Ω is bounded and closed, then Ω contains the segment condition, so according the Lemmas 2.2, 4.1 and 4.2, we can obtain the theorem.

    Lemma 4.3. Assume that the closed domain Ω is bounded, then S is dense in C2(Ω).

    According to Lemma 4.3, we can obtain the polynomial dense theory.

    Remark 4.1. Let ΩΥ be an arbitrary domain. Then the set of restrictions to Ω of functions in S5,2(Π1)×S5,2(Π2) is dense in C2(Ω), which leads to the set of restrictions to Ω of functions in (S5,2(Π1)×S5,2(Π2))2 is dense in (C2(Ω))2.

    Let Π1=[a,b], then the division is

    Π1:a=x0x1xN=b,

    h is the max length of the division. Si(x),Vi(x) and Wi(x) denote the Hermite splines of

    Si(x)={[xi+1xxi+1xi]3(6(xi+1xxi+1xi)215(xi+1xxi+1xi)5+10),x[xi,xi+1],[xxi1xixi1]3(6(xxi1xixi1)215(xxi1xixi1)5+10),x[xi1,xi],0,elsewhere,
    Vi(x)={3(xi+1x)5(xi+1xi)47(xi+1x)4(xi+1xi)4+4(xi+1x)3(xi+1xi)2,x[xi,xi+1],3(xxi1)5(xixi1)4+7(xxi1)4(xixi1)34(xxi1)3(xixi1)2,x[xi1,xi],0,elsewhere,
    Wi(x)={0.5(xi+1x)5(xi+1xi)3(xi+1x)4(xi+1xi)2+0.5(xi+1x)3(xi+1xi),x[xi,xi+1],0.5(xxi1)5(xixi1)3(xxi1)4(xixi1)2+0.5(xxi1)3(xixi1),x[xi1,xi],0,elsewhere.

    For Si(x),Vi(x) and Wi(x) from above, we have the following properties:

    Si(xk)=δik,Si(xk)=0,Si(xk)=0,Vi(xk)=0,Vi(xk)=δik,Vi(xk)=0,Wi(xk)=0,Wi(xk)=0,Wi(xi)=δik.

    Remark 5.1. Hermite bases are derived from segmented Hermite interpolating basis functions by the division Π1. On the kth divisions [xk,xk+1],k=0,,N1, it satisfies

    P(i)(xk)=f(i)(xk),i=0,1,2;k=0,,N,

    where P(x) is the interpolation polynomial and f(x) is the interpolated function, then the number of interpolating basis functions 6 is obtained. Thus, the total number of basis functions 6N on Π1 is obtained. The number of Hermite spline functions 3(N+1) is obtained from the definition of Si(x),Vi(x),Wi(x) above, and at the interior points Si(x),Vi(x) and Wi(x) is a function with two segments and at the endpoints is only a function with one segment. The following theorem will prove that Si(x),Vi(x) and Wi(x) are the bases.

    Theorem 5.1.

    {Hi(x)}3N+2i=0={Si(x)}Ni=0{Vi(x)}Ni=0{Wi(x)}Ni=0

    is linearly independent and is the base of S5,2(Π1).

    Proof. First, we will show that Si(x),Vi(x) and Wi(x) is linearly independent. Assume that,

    Ni=0ciSi(x)+Ni=0diVi(x)+Ni=0eiWi(x)=0.

    Due to the properties of the Hermite splines, when x=xk, ck=0, k=0,1,,N, then take the derivative of the above

    Ni=0diVi(x)+Ni=0eiWi(x)=0,

    when x=xk, dk=0, k=0,1,,N, then take the derivative of the above

    Ni=0eiWi(x)=0,

    when x=xk, ek=0, k=0,1,,N, so Sk(xi),Vk(xi) and Wk(xi) are linearly independent.

    Next, we will verify that it is a base of S5,2(Π1). Due to the definition of the S5,2(Π1), so Si(x), Vi(x), Wi(x)C2[a,b]. On the other hand, Si(x), Vi(x) and Wi(x) are a piecewise quintic polynomial. Thus, Si(x), Vi(x), Wi(x)S5,2(Π1).

    Since

    dimS5,2(Π1)=6N3(N1)=3N+3

    and

    dim{Si(x),Vi(x),Wi(x)}=3(N+1),

    so, {Si(x),Vi(x),Wi(x)} is a base of S5,2(Π1).

    Then according to Theorem 5.1 and Remark 4.1, it yields a new base

    SxyH(x)×H(y)

    on Υ is dense on Ω. So,

    U(x,y,sl)3N+2i=03N+2j=0dijlHi(x)×Hj(y)UN(x,y,sl),(x,y)Ω, (5.1)

    then, using the inverse Laplace transform based the Talbots strategy from Eq (3.5), we could obtain the numerical solution uN(x,y,t).

    Definition 5.1. For any ε>0, if

    GU(x,y,sl)W(x,y,sl)(C(Ω))2=max(x,y)Ω|GU(x,y,sl)W(x,y,sl)|<ε,

    then, U(x,y,sl) is an ε-approximate solution for Eq (4.1).

    We will provide the method of obtaining the ε-approximate solution. First, the minimum bounding rectangle

    Υ=[a,b]×[c,d]

    containing Ω is given. Subsequently, we will calculate residuals of two parts:

    (1) The residual inside Ω is defined as

    L1U(x,y,sl)GU(x,y,sl)W(x,y,sl)(C(Ω))2=2j=1GjU(x,y,sl)Wj(x,y,sl)C(Ω).

    (2) The residual on the boundary Ω is defined as

    L2U(x,y,sl)(ReU(x,y,sl)C(ΥΩ)+ImU(x,y,sl)C(ΥΩ)).

    For any ε>0, if there exists UN(x,y,sl) such that

    LUN(x,y,sl)=(L1+L2)(UN(x,y,sl))ε,

    so, UN(x,y,sl) is residual approximate solution of Eq (4.1) on Ω. If

    L(UN(x,y,sl))=minUN(x,y,sl)(L1+L2)(UN(x,y,sl))ε, (5.2)

    then UN(x,y,sl) is called the best ε-approximate solution.

    Lemma 5.1. G: (C2(Υ))2(C(Υ))2 is a bounded operator.

    Proof. For sl=(κl+iωl), and denoted that U(x,y,sl)Ul(x,y),

    CxDβLU(x,y,sl)(C)2=1Γ(2β)xa(xv)1β2U(v,y,sl)x2dv(C)21Γ(2β)xa(xv)1βUl(C2)2dv(C)2θ1Ul(C2)2,

    where θ1 is constants, and it could be similarly obtained that

    CxDβRU(x,y,sl)(C)2θ2Ul(C2)2,
    CyDγLU(x,y,sl)(C)2θ3Ul(C2)2

    and

    CyDγRU(x,y,sl)(C)2θ4Ul(C2)2.
    A(x,y)Δβ,γUl(x,y)=Kxcβ{A(x,y)(CxDβLUl(x,y)+CxDβRUl(x,y))+Ul(a,y)Γ(2β)(bx)β1(yc)γ1(dy)γ1+Ul(b,y)Γ(2β)(xa)β1(yc)γ1(dy)γ1}+Kycγ{Ul(x,c)Γ(2γ)(xa)β1(bx)β1(dy)γ1+Ul(x,d)Γ(2γ)(xa)β1(bx)β1(yc)γ1+A(x,y)(CyDγLUl(x,y)+CyDγRUl(x,y))}.

    Since A(x,y), (xa)β1,(bx)β1,(yc)γ1,(dy)γ1 is continuous, it has

    A(x,y)Δβ,γU(x,y,sl)(C)2KxcβU(a,y,sl)Γ(2β)(bx)β1(yc)γ1(dy)γ1+A(x,y)CxDβLU(x,y,sl)(C)2+KxcβA(x,y)CxDβRU(x,y,sl)+U(b,y,sl)Γ(2β)(xa)β1(yc)γ1(dy)γ1(C)2+KycγA(x,y)CyDγLU(x,y,sl)+U(x,c,sl)Γ(2γ)(xa)β1(bx)β1(dy)γ1(C)2+KycγA(x,y)CyDγRU(x,y,sl)+U(x,d,sl)Γ(2γ)(xa)β1(bx)β1(yc)γ1(C)2θ5Ul(C2)2,
    Re(sαqlU(x,y,sl))CRe(sαql)CRe(U(x,y,sl))Cθ6UlC2,Im(sαqlU(x,y,sl))Cθ7UlC2

    and

    rq=0aq=1,
    G1Ul(x,y)CA(x,y)Crq=0aq(θ6UlC2+θ7)+θ8ReUlC2θ9UlC2(Υ),

    similarly,

    G2Ul(x,y)Cθ10UlC2,

    hence,

    GU(x,y,sl)(C)2θUl(C2)2,

    so G is bounded.

    Theorem 5.2. Let U(x,y,sl) be the exact solution of Eq (4.1) on Ω, UN(x,y,sl) be the ε-approximate solution. For every ε>0, there exists N1, when NN1, coefficients dijl of UN1(x,y,sl) from Eq (5.1) satisfy Eq (5.2).

    Proof. U(x,y,sl) could be approximated by UN1(x,y,sl) on ΥΩ. For each fixed ε>0, there exists N1 such that the residual L(UN1(x,y,sl)) satisfies Eq (5.2).

    Let UN1(x,y,sl) be residual approximate solutions, taking min {ε4G,ε4}, in which G is defined by

    G=sup{Gu:u(C2)2,u(C2)21},

    there exists N1 such that the following two parts hold. Inside Ω, we suppose that

    UN1(x,y,sl)U(x,y,sl)C2(Ω)ε4G,

    when (x,y)Ω,

    L1UN1(x,y,sl)=2j=1Gj(UN1(x,y,sl)U(x,y,sl))C(Ω)2j=1GjUN1(x,y,sl)U(x,y,sl)C2(Ω)ε2.

    On the Ω, from the boundary condition U(x,y,sl)=0, we suppose that

    UN1(x,y,sl)U(x,y,sl)C(ΥΩ)UN1(x,y,sl)U(x,y,sl)C2(ΥΩ)ε4,

    hence, when (x,y) on the Ω,

    L2UN1(x,y,sl)=(Re(UN1(x,y,sl)U(x,y,sl))C(ΥΩ)+Im(UN1(x,y,sl)U(x,y,sl))C(ΥΩ))ε2.

    so,

    LUN1(x,y,sl)=(L1+L2)UN1(x,y,sl)ε

    and

    LUN1(x,y,sl)=minUN1(x,y,sl)L(UN1(x,y,sl))ε,

    so the theorem holds.

    Theorem 5.3. If Eq (4.1) is well-posed, then UN1(x,y,sl) obtained from Theorem 5.2 is the approximate solution of Eq (4.1) on Ω.

    Proof. Since UN1(x,y,sl) is the ε-approximate solution, for every ε>0, it yields,

    UN1(x,y,sl)U(x,y,sl)(C(Ω))2G1C(Ω)GUN1(x,y,sl)GU(x,y,sl)(C(Ω))2G1C(Ω)GUN1(x,y,sl)W(x,y,sl)(C(Ω))2G1C(Ω)ε,

    where G is bounded. It implies that UN1(x,y,sl) is the approximate solution of Eq (4.1) on Ω.

    Let S5,2(Π1) and S5,2(Π2) be two quintic spline space with partition

    Π1:a=x0<x1<<xn=b,
    Π2:c=y0<y1<<ym=d

    and

    Π1×Π2=[a,b]×[c,d].

    The quintic spline bases have the following properties.

    Theorem 5.4. Let u(x)Cm[a,b], 1m5, then there exists z(x)S5,2(Π1), such that

    (z(x)u(x))(k)C[a,b]KuCm[a,b]hm+1k,k=0,1,2,

    which h is the partition of the spline space, and K is the constant.

    Proof. The division of [a,b] is

    :a=x0<x1<<xj<xj+1<<xN=b,

    h is the max length of the division, and set subinterval

    πj=[xj,xj+1],j=0,1,,N1.

    ljk(x),ˉljk(x),ˉˉljk(x) be Hermite interpolation polynomials, j=0,1,,N,k=0,1, and satisfy

    1k=0ljk(x)=1,ljk(xi)=δik,ljk(xi)=0,lj(xi)=0,ˉljk(xi)=0,
    ˉljk(xi)=δik,ˉljk(xi)=0,ˉˉljk(xi)=0,ˉˉljk(xi)=0,ˉˉljk(xi)=δik.

    First, on [xj,xj+1], we prove

    1k=0ljk(x)(xjkx)p+1k=0ˉljk(x)p(xjkx)p1+1k=0ˉˉljk(x)p(p1)(xjkx)p2=0,

    1pd5, xj0=xj,xj1=xj+1, when p=1, it has

    1k=0ˉˉljk(x)p(p1)(xjkx)p2=0.

    Consider that, u(y)=(yx)p, it could be interpolated as follows:

    (yx)p=1k=0ljk(y)(xjkx)p+1k=0ˉljk(y)p(xjkx)p1+1k=0ˉˉljk(y)p(p1)(xjkx)p2.

    Setting y=x, we obtain that

    1k=0ljk(x)(xjkx)p+1k=0ˉljk(x)p(xjkx)p1+1k=0ˉˉljk(x)p(p1)(xjkx)p2=0.

    Next, due to property of the Hermite interpolation polynomial,

    1k=0ˉljk(x)C(πj)=1k=0ˉljk(x)ˉljk(xi)C(πj)=1k=0ˉljk(ξ)(xix)C(πj)K0h.

    Similarly,

    1k=0ˉˉljk(x)C(πj)=1k=0ˉˉljk(x)00C(πj)=1k=0ˉˉljk(x)ˉˉljk(xi)ˉˉljk(xi)(xxi)C(πj)1k=012ˉˉljk(ξ)(xxi)2C(πj)K1h2.

    For any uCm[a,b], there has z(x)S5,2(π), suppose that

    z(xj)=u(xj),z(xj)=u(xj),z(xj)=u(xj),

    so,

    z(x)u(x)C[a,b]=N1j=01k=0ljk(x)u(xjk)+1k=0ˉljk(x)u(xjk)+1k=0ˉˉljk(x)u(xjk)(mp=01p!pu(x)xp(1k=0ljk(x)(xjkx)p+1k=0ˉljk(x)p(xjkx)p1+1k=0ˉˉljk(x)p(p1)(xjkx)p2))C(πj)N1j=01k=0ljk(x)C(πj)u(xjk)mp=01p!pu(x)xp(xjkx)pC(πj)+1k=0ˉljk(x)C(πj)u(xjk)mp=01(p1)!pu(x)xp(xjkx)p1C(πj)+1k=0ˉˉljk(x)C(πj)u(xjk)mp=01(p2)!pu(x)xp(xjkx)p2C(πj)1k=0ljk(x)C(πj)1(m+1)!m+1u(x)xm+1(xjkx)m+1C(πj)+1k=0ˉljk(x)C(πj)1m!m+1u(x)xm+1(xjkx)mC(πj)+1k=0ˉˉljk(x)C(πj)(m1)!m+1u(x)xm+1(xjkx)m1C(πj)M1(m+1)!u(m+1)C[a,b]hm+1+M2(m)!u(m+1)C[a,b]hm+1+M3(m1)!u(m+1)C[a,b]hm+1K2u(m+1)C[a,b]hm+1.

    Then, on πj set

    z(i)(xjk)=u(i)(xjk),k=0,1; i=0,1,2,

    let

    w(x)=u(x)z(x),

    so

    w(xjk)=0,k=0,1;w(xjk)=w(xjk)=0,k=0,1,

    then let

    g(x)=w(x),

    so

    ξ[xj,xj+1],

    such that g(ξ)=0, and

    g(xjk)=g(xjk)=0,k=0,1.

    Due to z(x)P4, so z(x) is the polynomial interpolation of u(x) at the point (ξ,u(ξ)),(xjk,u(xjk)),(xjk,u(xjk)),k=0,1, so

    z(x)u(x)C[a,b]N1j=0z(x)u(x)C(πj)K3hmu(m+1)C[a,b].

    Then, let

    h(x)=w(x),

    so

    h(xjk)=0,k=0,1;

    η1(xj,ξ),η2(ξ,xj+1), such that

    h(η1)=g(η1)=0,h(η2)=g(η2)=0,

    due to zP3, so z(x) is the cubic polynomial interpolation of z(x), so

    z(x)u(x)C[a,b]N1j=0z(x)u(x)C(πj)K4hm1u(m+1)C[a,b].

    Finally,

    D(k)(z(x)u(x))C[a,b]Ku(m+1)C[a,b]hm+1k,k=0,1,2.

    According to [30], the following lemma is given.

    Lemma 5.2. Let u(x,y)Cm(Ω), 2m6, then there exists

    z(x,y)S5,2(Π1)×S5,2(Π2),

    such that

    (zu)(k,l)(x,y)C(Ω)λu(m+1,m+1)Cm(Ω)hm(k+l),k,l=0,1,2,

    which h is the partition of the space, and λ is the constant.

    According to Theorem 5.4 and Lemma 5.2, we can infer that:

    Remark 5.2. Let u(x,y)C4(Ω), then there exists

    z(x,y)S5,2(Π1)×S5,2(Π2),

    such that

    z(x,y)u(x,y)C2(Ω)λuC4(Ω)h2,

    which h is the partition of the space, and λ is the constant.

    Theorem 5.5. The numerical solution ˜UN(x,y,sl) obtained from the proposed meshless method converges to the exact solution U(x,y,sl).

    Proof. Owing to Theorem 5.1, Hi(x) is the base of S5,2(π), so numerical solution UN,l(x,y) obtained from Eq (5.1) belongs to S5,2(π)×S5,2(π). From the Remark 5.2 and Theorem 5.3, we have

    G˜UN(x,y,sl)W(x,y,sl)(C)2λ1UlC2h2,˜UN(x,y,sl)(C(Ω))2λ2UlC2h2.

    Assume that

    G˜UN(x,y,sl)=W(x,y,sl),˜UN(x,y,sl)|Ω=w(x,y,sl),

    so

    \| \mathbb{W}^*(x, y, s_l)-\mathbb{W}(x, y, s_l) \|_{(C)^2} \le \lambda_1 \| \mathcal{U}_{l} \|_{C_2}h^2, \quad \| w^*(x, y, s_l)\|_{(C(\partial \Omega))^2} \le \lambda_2 \| \mathcal{U}_{l} \|_{C_2}h^2.

    Then \exists N , such that

    \begin{equation*} \begin{split} &\|\tilde{\mathcal{U}}_N(x, y, s_l)-\mathcal{U}(x, y, s_l)\|_{C_2}\le \lambda_3 \| \mathcal{U}_{l} \|_{C_2}h^2, \end{split} \end{equation*}

    where \lambda_1, \lambda_2, \lambda_3 are constants.

    In this section we give two examples to demonstrate the effectiveness of our theoretical analysis. The examples will discuss a single time fractional term and a multiple time fractional term on different domains, respectively. Calculate the

    L_{\infty}(t) = \max\limits_{N}\left|u(x, y, t)-u_N(x, y, t) \right|

    and

    E(t) = \left\|u(x, y, t)-u_N(x, y, t) \right\|_{L_2} = \left(\int_{\Omega}\left(u(x, y, t)-u_N(x, y, t)\right)^2 {\rm d}\Omega\right)^{1/2},

    where u(x, y, t) is the exact solution, u_N(x, y, t) is the approximate solution by our method. If t = 1, L_{\infty} = L_{\infty}(1) . Meanwhile, let the L = 10 of Bromwich be integrated by the inverse Laplace transform. The node

    n\triangleq N+1

    from Eq (5.1).

    Example 6.1. Consider the single term form Eq (3.1), where r = 1, K_x = K_y = 1 ,

    \begin{eqnarray*} {^C_0\!D}^{\alpha}_t u (x, y, t) = \frac{\partial^{\beta}u(x, y, t)}{\partial |x|^{\beta}}+\frac{\partial^{\gamma}u(x, y, t)}{\partial |y|^{\gamma}}+f(x, y, t), \; \; \; \; (x, y)\in\Omega, \; 0 < t\le 1, \end{eqnarray*}

    with

    u(x, y, t)|_{\partial \Omega} = 0, \; u(x, y, 0) = 0.

    Let \alpha = 2/3 , \beta = 3/2 , \gamma = 5/4 . L_{\infty} of Example 6.1 on rectangular domain and circular domain are shown in Table 1.

    Table 1.  L_{\infty} on difference domains for Example 6.1.
    Node n L_{\infty} in (1) of Example 6.1 L_{\infty} in (2) of Example 6.1
    2\times2 3.95185\times10^{-9} 1.54171\times10^{-8}
    3\times3 3.95185\times10^{-9} 1.56203\times10^{-8}
    4\times4 3.95184\times10^{-9} 1.58076\times10^{-8}

     | Show Table
    DownLoad: CSV

    It can be concluded that our method is valid in a verifiable way and that it gives better results in the general case of smoother time solutions.

    (1) When (x, y) on rectangular domains, \Omega = [0, 1]\times [0, 1], the true solution is

    u(x, y, t) = x^2(1-x)^2y^2(1-y)^2t^{\frac{4}{5}}.

    The error figure is shown in the Figure 1a at n = 3\times 3 . And the error L_{\infty} are shown in Table 1.

    Figure 1.  Error for Example 6.1, when \alpha = 2/3 , \beta = 3/2 , \gamma = 5/4 , t = 1 .

    (2) When (x, y) on circular domains \Omega ,

    \Omega = \{(x, y)|(x, y)\in (x-1/2)^2+(y-1/2)^2\le 1/4 \},

    the true solution is

    u(x, y, t) = (x-1/2)^2(y-1/2)^2t^{\frac{4}{5}}.

    The figure of error u(x, y, t)-u_N(x, y, t) when t = 1 at n = 3\times 3 is shown in the Figure 1b. And the error L_{\infty} are shown in Table 1.

    With the above two numerical examples we find that our method gets high accuracy on different regions, showing that our method can handle arbitrary convex regions. Our error convergence is second-order, and since the solution of u with respect to the space is x^2(1-x)^2y^2(1-y)^2 or (x-1/2)^2(y-1/2)^2 , we have fewer points to get a high accuracy error, which is in accordance with the theory. At the same time, the solution of u with respect to the time is t^{\frac{4}{5}} , Laplace transform can be used to deal with lower order smooth solutions.

    Example 6.2. [31] Consider the multi-term from Eq (3.1), {where r = 4, K_x = K_y = 1

    \begin{eqnarray*} \sum\limits_{q = 0}^{4}a_q \left({^C_0\!D}^{\alpha_q}_t \right)u (x, y, t) = \frac{\partial^{\beta}u(x, y, t)}{\partial |x|^{\beta}}+\frac{\partial^{\gamma}u(x, y, t)}{\partial |y|^{\gamma}}+f(x, y, t), \; \; \; \; (x, y)\in \Omega, \; 0 < t\le 1, \end{eqnarray*}

    with

    u(x, y, t)|_{\partial \Omega} = 0, \; u(x, y, 0) = 0,
    \begin{align*} &f(x, y, t) = \sum\limits_{i = 0}^{4}a_i t^{\frac{\alpha_0+1}{2}-\alpha_i}E_{1, \frac{\alpha_0+1}{2}-\alpha_i+1}(t)x^2(1-x)^2y^2(1-y)^2\\ &+\frac{t^{\frac{\alpha_0+1}{2}} E_{1, \frac{\alpha_0+1}{2}+1}(t)y^2(1-y)^2 }{\cos(\beta \pi/2)}\left\{ 2\frac{x^{2-\beta}+(1-x)^{2-\beta}}{\Gamma(3-\beta)} - 12 \frac{x^{3-\beta}+(1-x)^{3-\beta}}{\Gamma(4-\beta)} +24 \frac{x^{4-\beta}+(1-x)^{4-\beta}}{\Gamma(5-\beta)} \right\}\\ &+\frac{t^{\frac{\alpha_0+1}{2}} E_{1, \frac{\alpha_0+1}{2}+1}(t)x^2(1-x)^2 }{\cos(\gamma \pi/2)}\left\{ 2\frac{y^{2-\gamma}+(1-y)^{2-\gamma}}{\Gamma(3-\gamma)} -12\frac{y^{3-\gamma}+(1-y)^{3-\gamma}}{\Gamma(4-\gamma)} +24 \frac{y^{42-\gamma}+(1-y)^{4-\gamma}}{\Gamma(5-\gamma)} \right\} \end{align*}

    where

    E_{a, b}(t): = \sum\limits_{i = 0}^{\infty}\frac{t^i}{\Gamma(ai+b)}.

    Then the exact solution is

    u(x, y, t) = t^{\frac{\alpha_0+1}{2}}E_{1, \frac{\alpha_0+1}{2}+1}(t)x ^2(1-x)^2y^2(1-y)^2.

    (1) When (x, y) on rectangular domains,

    \Omega = [0, 1]\times [0, 1].

    When

    {\bf \alpha} = (0.05, 0.08, 0.1, 0.15, 0.2), \; \; \; {\bf a} = (3/10, 1/10, 3/20, 1/5, 1/4), \; \; \; \beta = 1.6, \gamma = 1.6.

    We calculate the error E(T) and compare it with [31] in Table 2 at T = 1 .

    Table 2.  Error E(T) when T = 1 for (1) of Example 6.2.
    Mesh length h [31] Node n E(T)
    1/8 1.3862\times10^{-4} 2\times2 1.44824\times10^{-6}
    1/16 3.1353\times10^{-5} 3\times3 1.44824\times10^{-6}
    1/24 1.3203\times10^{-5} 4\times4 1.44824\times10^{-6}

     | Show Table
    DownLoad: CSV

    (2) When (x, y) on circular domains, the

    \Omega = \{(x, y)|(x-1/2)^2+(y-1/2)^2\le 1/4\}.

    When

    {\bf \alpha} = (0.35, 0.45, 0.6, 0.7, 0.8), \; \; \; {\bf a} = (3/10, 1/5, 4/30, 1/6, 1/5),

    where \beta = 1.02 , \gamma = 1.02.

    The numerical solution and the absolute errors when t = 1 at n = 3\times 3 are shown in Figure 2.

    Figure 2.  Numerical solution(left) and absolute errors(right), when t = 1 , n = 3 .

    We calculated the L_2 error E(T) on the rectangular domain and compared it with [31]. It can be seen that we obtain higher accuracy with fewer points, which proves the high efficiency of our method. We also carry out experiments with different parameters \alpha, {\bf a} and \beta, \gamma on the circular domain and calculate the absolute errors at the n = 3 .

    From Figure 2, it can be seen that our method also achieves high error accuracy, indicating the applicability of our method. The high error accuracies obtained by our method in different regions and also with different parameters show the stability and efficiency of our method. Because of the high smoothness of u with respect to x, y , we get high error accuracy with fewer points, which is consistent with our theoretical analysis.

    In this paper, we proposed a meshless method of solving the minimum residual approximate solution for Eq (3.1). Different from previous methods, we use the Laplace transform method to deal with the multi-term time fractional operator, we transform the time into complex frequency domain by Laplace transform, Eq (3.1) is transformed into complex equation Eq (3.3). Then, on the spatial direction, we proposed a quintic Hermite meshless method to deal with space fractional operators on arbitrary convex region based on the theory of polynomial functions dense theorem. The approximate accuracies become higher by increasing number of Quintic Hermite spline functions. The minimum residual approximate solution of Eq (4.1) is obtained by Theorems 5.3 under the condition of well-posed equations. Meanwhile, using Theorem 5.4 and Lemma 5.2, it infers Remark 5.2, which is the convergence of the biquintic spline function. Then by using Remark 5.2 and Theorem 5.3, we can obtain Theorem 5.5 to show the convergence of the method in the spatial direction. We use numerical inversion methods to transform the obtained the minimum residual approximate solution from the Laplace domain to the real domain by using the strategy of Talbot through parabolic path.

    In Numerical examples, we fix the L = 10 in Eq (3.5) by parabolic path to get the numerical solution. First, we handle the single term time-space fractional diffusion equations, we can deduce that the method can deal with time fractions that are not sufficiently smooth, and we can get higher precision with fewer nodes in arbitrary convex region from Table 1 and Figure 1. This also proves that Laplace transform is effective for dealing with insufficiently smooth time-fractional operators. Then, we solve the multi-term time-space fractional diffusion equations with 4 terms. These results are compared with [31], and it is found that our method achieves better accuracy with fewer points. At the same time, we found that the accuracy of the single term is better than that of the multi-term. In addition, the accuracy is higher on rectangular areas than on circular areas. These experimental results are consistent with theoretical expectations and demonstrate the effectiveness and efficiency of our method.

    In this paper, the use of the extension theorem allows the meshless method to be applied to arbitrary convex regions in two dimensions, and the use of the Laplace transform allows to deal with multi-term low-order time solutions. In the future, through the study of spatial Riesz operators, we will investigate meshless methods for solving equations in arbitrary regions of higher dimensions. In addition, this method can also be used to study equations of time-distributed order.

    The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

    The authors would like to thank the referees for their helpful comments and suggestions, which lead to a much improved version of this paper.

    The authors declare that they have no conflicts of interest.



    [1] World Health Organization, Anti-tuberculosis Drug Resistance in the World: Third Global Report, No. WHO / HTM / TB / 2004. 343, World Health Organization, 2004.
    [2] J. Cohen, Extensively drug-resistant TB gets foothold in South Africa, Science 313 (2006), 1554.
    [3] L. B. Reichman, J. H. Tanne, Timebomb: The Global Epidemic of Multi-Drug Resistant Tuberculosis, 2002.
    [4] Y. Zhou, K. Khan, Z. Feng, et al. Projection of tuberculosis incidence with increasing immigration trends, J. Theor. Biol., 254 (2008), 215–228.
    [5] P. Rodrigues, M. G. M. Gomes, C. Rebelo, Drug resistance in tuberculosis: a reinfection model, Theor. Popul. Biol., 71 (2007), 196–212.
    [6] N. Blaser, C. Zahnd, S. Hermans, et al. Tuberculosis in Cape Town: an age-structured transmission model, Epidemics, 14 (2016), 54–61.
    [7] C. P. Bhunu, W. Garira, Z. Mukandavire, et al. Tuberculosis transmission model with chemoprophylaxis and treatment, B. Math. Biol., 70 (2008), 1163–1191.
    [8] Centers for Disease Control and Prevention, TB in Children in the United States, CDC, 2014. Available from: https://www.cdc.gov/tb/topic/populations/tbinchildren/default.htm.
    [9] S. M. Blower, P. M. Small, P. C. Hopewell, Control strategies for tuberculosis epidemics: new models for old problems, Science, 273 (1996), 497–500.
    [10] Fatmawati, H. Tasman, An optimal treatment control of TB-HIV coinfection, International Journal of Mathematics and Mathematical Sciences, 2016 (2016).
    [11] R. I. Hickson, G. N. Mercer, K. M. Lokuge, A metapopulation model of tuberculosis transmission with a case study from high to low burden areas, PLoS One, 7 (2012).
    [12] R. M. G. J. Houben, T. Sumner, A. D. Grant, et al. Ability of preventive therapy to cure latent Mycobacterium tuberculosis infection in HIV-infected individuals in high-burden settings, P. Natl. A. Sci., 111 (2014), 5325–5330.
    [13] R. Kaplan, J. Caldwell, K. Middelkoop, et al. Impact of ART on TB case fatality stratified by CD4 count for HIV-positive TB patients in Cape Town, South Africa (2009-2011), J. Acq. Imm. Def., 66 (2014), 487–494.
    [14] Fatmawati, U. D. Purwati, F. Riyudha, et al. Optimal control of a discrete age-structured model for tuberculosis transmission, Heliyon, 6 (2020).
    [15] S. G, Samko, A. A, Kilbas, O. I. Marichev, Fractional Integrals and Derivatives: Theory and Applications, Gordon and Breach Science Publishers, 1993.
    [16] I. Podlubny, Fractional Differential Equations, Academic Press, San Diego, California, USA, 1999.
    [17] T. Sardar, S. Rana, J. Chattopadhyay, A mathematical model of dengue transmission with memory, Commun. Nonlinear Sci., 22 (2015), 511–525.
    [18] J. Huo, H. Zhao, L. Zhu, The effect on backward bifurcation in a fractional orde HIV model, Nonlinear Analysis: Real World Applications, 26 (2015), 289–305.
    [19] M. Saeedian, M. Khalighi, N. Azimi-Tafreshi, et al. Memory effects on epidemic evolution: the susceptible-infected-recovered epidemic model, Phys. Rev. E., 95 (2017).
    [20] C. M. A. Pinto, A. R. M. Carvalho, The HIV/TB coinfection severity in the presence of TB multidrug resistant strains, Ecol. Complex., 32 (2017), 1–20.
    [21] Fatmawati, E. M. Shaiful, M. I. Utoyo, A fractional order model for HIV dynamics in a two-sex population, International Journal of Mathematics and Mathematical Sciences, 2018 (2018).
    [22] G. C. Wu, Z. G. Deng, D. Baleanu, et al. New variable-order fractional chaotic systems for fast image encryption, Chaos: An Interdisciplinary Journal of Nonlinear Science, 29 (2019).
    [23] G. C. Wu, T. Abdeljawad, J. Liu, et al. Mittag-Leffler stability analysis of fractional discrete-time neural networks via fixed point technique, Nonlinear Analysis: Modelling and Control, 24 (2019), 919–936.
    [24] G. C. Wu, D. Q. Zeng, D. Baleanu, Fractional impulisve differential equations: Exact solutions, integral equations and short memory case, Frac. Calc. Appl. Anal., 22 (2019), 180–192.
    [25] M. Itik, S. P. Banks, Chaos in a three-dimensional cancer model, Int. J. Bifurcat. Chaos, 20 (2010), 71–79.
    [26] M. Caputo, M. Fabrizio, A new definition of fractional derivative without singular kernel, Progr. Fract. Differ. Appl., 1 (2015), 1–13.
    [27] M. Caputo, M. Fabrizio, On the notion of fractional derivative and applications to the hysteresis phenomena, Meccanica, 52 (2017), 3043–3052.
    [28] T. Zhang, L. Xiong, Periodic motion for impulsive fractional functional differential equations with piecewise Caputo derivative, Appl. Math. Lett., 101 (2020), 106072.
    [29] A. Atangana, D. Baleanu, New fractional derivatives with nonlocal and non-singular kernel: theory and application to heat transfer model, Therm. Sci., 20 (2016).
    [30] A. Atangana, I. Koca, Chaos in a simple nonlinear system with Atangana-Baleanu derivatives with fractional order, Chaos, Solitons & Fractals, 89 (2016), 447–454.
    [31] A. Atangana, K. M. Owolabi, New numerical approach for fractional differential equations, Math. Model. Nat. Pheno., 13 (2018).
    [32] A. Atangana, J. F. Gomez-Aguilar, Decolonisation of fractional calculus rules: Breaking commutativity and associativity to capture more natural phenomena, The European Physical Journal Plus, 133 (2018).
    [33] M. A. Khan, S. Ullah, M. Farooq, A new fractional model for tuberculosis with relapse via Atangana-Baleanu derivative, Chaos, Solitons & Fractals, 116 (2018), 227–238.
    [34] E. Bonyah, Chaos in a 5-D hyperchaotic system with four wings in the light of non-local and non-singular fractional derivatives, Chaos, Solitons & Fractals, 116 (2018), 316–331.
    [35] K. Muhammad Altaf, A. Atangana, Dynamics of Ebola disease in the framework of different fractional derivatives, Entropy, 21 (2019).
    [36] R. Jan, M. A. Khan, P. Kumam, et al, Modeling the transmission of dengue infection through fractional derivatives, Chaos, Solitons & Fractals, 127 (2019), 189–216.
    [37] W. Wang, M. A. Khan, P. Kumam, et al. A comparison study of bank data in fractional calculus, Chaos, Solitons & Fractals, 126 (2019), 369–384.
    [38] M. A. Khan, The dynamics of a new chaotic system through the Caputo-Fabrizio and AtanaganBaleanu fractional operators, Adv. Mech. Eng., 11 (2019).
    [39] Fatmawati, M. A. Khan, M. Azizah, et al. A fractional model for the dynamics of competition between commercial and rural banks in Indonesia, Chaos, Solitons & Fractals, 122 (2019), 32–46.
    [40] S. Ullah, M. A. Khan, M. Farooq, et al. A fractional model for the dynamics of tuberculosis (TB) using Atangana-Baleanu derivative, Discrete Cont. Dyn. S, 13 (2019).
    [41] A. A. Velayati, Tuberculosis in children, International Journal of Mycobacteriology, 5 (2016).
    [42] C. Castillo-Chaves, B. Song, Dynamic models of tuberculosis and their applications, Math. Biosci. Eng., 1 (2004), 361–404.
    [43] Z. M. Odibat, N. T. Shawagfeh, Generalized Taylors formula, Appl. Math. Comput., 186 (2007), 286–293.
    [44] W. Lin, Global existence theory and chaos control of fractional differential equations, J. Math. Anal. Appl., 332 (2007), 709–726.
    [45] P. van den Driessche, J. Watmough, Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission, Math. Biosci., 180 (2002), 29–48.
    [46] K. Diethelm, N. J. Ford, A. D. Freed, A predictor-corrector approach for the numerical solution of fractional differential equations, Nonlinear Dynamics, 29 (2002), 3–22.
    [47] K. Diethelm, N. J. Ford, A. D. Freed, Detailed error analysis for a fractional Adams method, Numerical Algorithms, 36 (2004), 31–52.
    [48] C. P. Bhunu, Mathematical analysis of a three-strain tuberculosis transmission model, Appl. Math. Model., 35 (2011), 4647–4660.
    [49] S. Athithan, M. Ghosh, Optimal control of tuberculosis with case detection and treatment, World Journal of Modelling and Simulation, 11 (2015), 111–122.
    [50] Word Health Organization, Factsheet on the World Tuberculosis Report 2017, WHO, 2017. Available from: https://www.who.int/en/news-room/fact-sheets/detail/tuberculosis.
    [51] J. J. Tewa, S. Bowong, B. Mewoli, Mathematical analysis of two-patch model for the dynamical transmission of tuberculosis, Appl. Math. Model., 36 (2012), 2466–2485.
    [52] M. Toufik, A. Atangana, New numerical approximation of fractional derivative with non-local and non-singular kernel: application to chaotic models, The European Physical Journal Plus, 132 (2017), 444.
    [53] L. Xu, H. Hub, F. Qinc, Ultimate boundedness of impulsive fractional differential equations, Appl. Math. Lett., 62 (2016), 110–117.
    [54] L. Xu, J. Li, S. S. Ge, Impuls ivestabilization of fractional differential systems, ISA T., 70 (2017), 125–131.
    [55] L. Xu, X. Chu, H. Hu, Exponential ultimate boundedness of non-autonomous fractional differential systems with time delay and impulses, Appl. Math. Lett., 99 (2020), 106000.
  • This article has been cited by:

    1. Zehui Ma, Rahmatjan Imin, A High‐Precision Meshless Method for Time‐Fractional Mixed Diffusion and Wave Equations, 2025, 126, 0029-5981, 10.1002/nme.70020
  • Reader Comments
  • © 2020 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)
通讯作者: 陈斌, bchen63@163.com
  • 1. 

    沈阳化工大学材料科学与工程学院 沈阳 110142

  1. 本站搜索
  2. 百度学术搜索
  3. 万方数据库搜索
  4. CNKI搜索

Metrics

Article views(8260) PDF downloads(973) Cited by(56)

Figures and Tables

Figures(11)  /  Tables(2)

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog