Fitted reproducing kernel Hilbert space method for the solutions of some certain classes of time-fractional partial differential equations subject to initial and Neumann boundary conditions

Computers and Mathematics with Applications (

)

–

Contents lists available at ScienceDirect

Computers and Mathematics with Applications journal homepage: www.elsevier.com/locate/camwa

Fitted reproducing kernel Hilbert space method for the solutions of some certain classes of time-fractional partial differential equations subject to initial and Neumann boundary conditions Omar Abu Arqub Department of Mathematics, Faculty of Science, Al-Balqa Applied University, Salt 19117, Jordan

article

info

Article history: Available online xxxx Keywords: Reproducing kernel Hilbert space method Time-fractional partial differential equations Neumann boundary conditions (Heat, cable, anomalous subdiffusion, reaction subdiffusion, Fokker–Planck, Fisher’s, and Newell–Whitehead) equations

abstract Latterly, many problems arising in different fields of science and engineering can be reduced, by applying some appropriate discretization, to a series of time-fractional partial differential equations. Unlike the normal case derivative, the differential order in such equations is with a fractional order, which will lead to new challenges for numerical simulation. The purpose of this analysis is to introduce the reproducing kernel Hilbert space method for treating classes of time-fractional partial differential equations subject to Neumann boundary conditions with parameters derivative arising in fluid-mechanics, chemical reactions, elasticity, anomalous diffusion, and population growth models. The method provides appropriate representation of the solutions in convergent series formula with accurately computable components. Numerical experiments with different order derivatives degree are performed to support the theoretical analyses which are acquired by interrupting the n-term of the exact solutions. Finally, the obtained outcomes showed that the proposed method is competitive in terms of the quality of the solutions found and is very valid for solving such time-fractional Neumann problems. © 2016 Elsevier Ltd. All rights reserved.

1. Introduction The one-dimensional nonlinear fractional models, which are representative in the time-fractional partial differential equations (PDEs), have been studied to describe numerous realism matters successfully not only in physics, but also in engineering, biology, economics, and other sciences [1–5]. Such models are utilized extensively by many experts to explain their complicated structures easily, simplified the controlling design without any loss of hereditary behaviors as well as create nature issues closely understandable for these phenomena. Consequently, fractional derivatives provide more accurate models of realism problems than integer-order derivatives; they are actually found to be a suitable tool to describe certain physical and engineering problems including reaction diffusion models, dynamical mathematical models, electrical circuits models, signal processing models, and so on [6–13]. Developing analytical and numerical methods for the solutions of time-fractional PDEs is a very important task. Indeed, it is difficult to obtain exact solutions form in general for most cases. Therefore, attempts have been made to propose analytical methods that approximate the exact solutions of such equations [6–35]. The purpose of this analysis is to investigate and implement a computational iterative method, the reproducing kernel Hilbert space method (RKHSM), in finding approximate solutions for various certain classes of Neumann time-fractional E-mail address: [email protected]. http://dx.doi.org/10.1016/j.camwa.2016.11.032 0898-1221/© 2016 Elsevier Ltd. All rights reserved.

2

O. Abu Arqub / Computers and Mathematics with Applications (

)

–

PDEs with parameters derivative in the sense of Riemann–Liouville and Caputo fractional derivatives. More specifically, we consider two types of time-fractional models in the fractional operator form. (I) The following general form of nonlinear second-order time-fractional PDE with constant coefficients:

     1−β  2 1−β [u (x, t )] + f (x, t , u (x, t )) , ∂t u (x, t ) = Aα 0 Dt1−α + Bβ 0 Dt ∂x2 u (x, t ) − Cα 0 Dt1−α + Dβ 0 Dt

(1)

subject to the following initial and Neumann boundary conditions: u (x, 0) = ω (x) ,

(2)

∂x u (L, t ) = υ2 (t ) ,

∂x u (0, t ) = υ1 (t ) ,

where 0 < α, β < 1, 0 ≤ t ≤ T ∈ R, 0 ≤ x ≤ L ∈ R, Aα , Bβ , Cα , Dβ are nonnegative real constants, f is continuous real1−γ

valued function, u is an unknown function to be determined. Here, 0 Dt derivatives operator of order 1 − γ of a function u (x, t ) and defined as 1−γ u 0 Dt

( x, t ) =

1

Γ (γ )

∂t

t



(t − τ )γ −1 u (x, τ ) dτ ,

denote the Riemann–Liouville time-fractional

0 < τ < t , 0 < γ < 1.

(3)

0

To specify more, the time-fractional PDE of Eq. (1) consists of the following well-known certain equations as special cases:

• If Bβ = Cα = Dβ = 0 and f (x, t , u (x, t )) = g (x, t ), then we obtain the time-fractional heat equation. The heat equation is derived from Fourier’s law and conservation of energy, it is used in describing the distribution of heat or variation in temperature in a given region over time [14,15]. • If Bβ = Cα = 0, then we obtain the time-fractional cable equation. The cable equation is derived from the cable equation for electrodiffusion in smooth homogeneous cylinders, it is occurred due to anomalous diffusion and is used in modeling the ion electrodiffusion at the neurons [16,17]. • If Cα = Dβ = 0 and f (x, t , u (x, t )) = g (x, t ), then we obtain the time-fractional modified anomalous subdiffusion equation. The modified anomalous subdiffusion equation is derived from the neural cell adhesion molecules, it is used for describing processes that become less anomalous as time progresses by the inclusion of a second fractional time derivative acting on the diffusion term [18,19]. (II) The following general form of nonlinear second-order time-fractional PDE with variable coefficients:

  ∂tαα u (x, t ) = ∂x22 p (x, u (x, t )) − ∂x q (x, u (x, t )) [u (x, t )] + f (x, t , u (x, t )) ,

(4)

subject to the following initial and Neumann boundary conditions: u (x, 0) = ω1 (x) ,

∂x u (0, t ) = υ1 (t ) ,

(5)

∂x u (L, t ) = υ2 (t ) ,

where 0 < α < 1, 0 ≤ t ≤ T ∈ R, 0 ≤ x ≤ L ∈ R, f is continuous real-valued functions, u is an unknown function to be determined. Here, ∂tαα = ∂ α /∂ t α denote the Caputo time-fractional derivatives operator of order α of a function u (x, t ) and defined as

∂tαα u (x, t ) =

1

Γ (1 − α)

t



(t − τ )−α ∂τ u (x, τ ) dτ ,

0 < τ < t , 0 < α < 1.

(6)

0

To specify more, the time-fractional PDE of Eq. (4) consists of the following well-known certain equations as special cases:

• If p (x, u (x, t )) = 1, q (x, u (x, t )) = 0, and f (x, t , u (x, t )) = −Eα u (x, t ) + g (x, t ), then we obtain the timefractional reaction subdiffusion equation. The reaction subdiffusion equation appears in many different areas of chemical reactions, such as exciton quenching, recombination of charge carriers or radiation defects in solids, and predator–pray relationships in ecology [20,21]. • If f (x, t , u (x, t )) = 0, then we obtain the time-fractional Fokker–Planck equation [22,23]. The Fokker–Planck equation arises in many phenomena in plasma and polymer physics, population dynamics, neurosciences, nonlinear hydrodynamics, pattern formation, and psychology [24–27]. • If p (x, u (x, t )) = 1, q (x, u (x, t )) = 0, and f (x, t , u (x, t )) = Fα u (x, t ) (1 − Gα un (x, t )) + g (x, t ), then we obtain the time-fractional Fisher’s equation when n = 1 and the time-fractional Newell–Whitehead equation when n = 2. The Fisher’s equation is used to describe the population growth models [28,29], whilst, the Newell–Whitehead equation arises in fluid dynamics model and capillary–gravity waves [30,31]. The RKHSM is a numerical, as well as, analytical technique for solving a large variety of ordinary and PDEs associated to different kinds of order derivatives degree, and usually provides the solutions in terms of rapidly convergent series with components that can be elegantly computed. The advantages of the utilized approach lie in the following; firstly, it

O. Abu Arqub / Computers and Mathematics with Applications (

)

–

3

can produce good globally smooth numerical solutions, and with ability to solve many fractional differential systems with complex constraints conditions, which are difficult to solve; secondly, the numerical solutions and their derivatives converge uniformly to the exact solutions and their derivatives, respectively; thirdly, the numerical solutions and all their derivatives are applicable for each arbitrary point in the given domain; fourthly, the method is mesh-free, needs no time discretization, and easily implemented. The theory of reproducing kernel was used for the first time at the beginning of the 20th century as a novel solver for the boundary value problems of harmonic and biharmonic function types. This theory, which is representative in the RKHSM, has been successfully applied to various important applications in numerical analysis, computational mathematics, image processing, machine learning, probability and statistics, and finance [36–39]. The RKHSM is a useful framework for constructing numerical solutions of great interest to applied sciences. In the recent years, based on this theory, extensive work has been proposed and discussed for the numerical solutions of several integral and differential operators side by side with their theories. It has been applied in the numerical solutions of PDEs [40–43], in the numerical solutions of fractional PDEs [32–35], in the numerical solutions of differential equations of fractional order [44–46], in the numerical solutions of differential algebraic systems [47,48], in the numerical solutions of fuzzy differential equations [49–51], in the numerical solutions of integral equations [52,53], in the numerical solutions of integrodifferential equations [54,55], in the numerical solutions of singularly perturbed differential equations [56–58], and so on. In fact, the application of the RKHSM in the field of fractional PDEs is not new; in [32] the authors have been used the RKHSM for solving the time-fractional telegraph equation. Authors in [33] have been applied the RKHSM to the solution of certain class of time–space-fractional PDEs subject to integral boundary conditions. The RKHSM is carried out in [34] for solving the time–space-fractional advection–dispersion equation. Recently, the RKHSM for solving the time–spacefractional modified anomalous subdiffusion equation is proposed in [35]. The structure of the present paper is as follows. In the next section, we utilized some necessary definitions and results from the reproducing kernel theory. In Section 3, two extended inner product spaces and two extended reproducing kernel functions are constructed. The formulation of the considered problems is presented in the spaces W (Ω ) and H (Ω ) as shown in Section 4. Essential results for the constructed method is given shortly in Section 5. Representation of the numerical solutions is described in Section 6 based upon the reproducing kernel theory. In Section 7, convergent theorem and error behavior of the solutions are presented. Numerical algorithm and numerical outcomes are discussed as utilized in Section 8. A final section provides brief conclusions. 2. Background for the reproducing kernel theory In linear algebra, an inner product space is a vector space with an additional structure called an inner product. This additional structure associates each pair of vectors in the space with a scalar quantity known as the inner product of the vectors. Inner products allow the rigorous introduction of intuitive geometrical notions such as the length of a vector or the means of defining orthogonality between vectors. In this section, several inner product spaces and several reproducing kernel functions are constructed on the finite domains [0, T ] and [0, L], respectively. Before any further discussions, some fundamental concepts from the reproducing kernel theory are needed to describe the RKHSM. For illustration, a Hilbert space which possesses a reproducing kernel is called a RKHS. Definition 1 ([40]). The space W21 [0, T ] is defined as W21 [0, T ] = {z = z (t ) : z is one-variable absolutely continuous realvalued function on [0, T ] and z ′ ∈ L2 [0, T ]}. Whilst, the inner product and the norm of W21 [0, T ] are given, respectively, as

⟨z1 (t ), z2 (t )⟩W 1 = z1 (0) z2 (0) + 2

T



z1′ (t )z2′ (t )dt ,

(7)

0

and ∥z1 ∥2W 1 = ⟨z1 (t ) , z1 (t )⟩W 1 , where z1 , z2 ∈ W21 [0, T ]. 2

2

Theorem 1 ([40]). The Hilbert space W21 [0, T ] is a complete reproducing kernel with the reproducing kernel function Rs{1} (t ) =



1 + t, 1 + s,

t ≤ s, t > s.

(8)

 1 [0, L] = {z = z (x) : z is one-variable In similar fashion, if [0, L] is the domain space in the x direction, then the space W 2 ′ 2  1 [0, L] can be given absolutely continuous real-valued function on [0, L] and z ∈ L [0, L]}. Whilst, the inner product of W 2 as ⟨z1 (x), z2 (x)⟩W  1 = z1 (0) z2 (0) + 2

L 0

{1}

z1′ (x)z2′ (x)dx with the reproducing kernel function  Ry (x) = 1 + x if x ≤ y and

{1}  Ry (x) = 1 + y if x > y.

Definition 2 ([41]). The space W22 [0, T ] is defined as W22 [0, T ] = {z = z (t ) : z , z ′ are one-variable absolutely continuous real-valued functions on [0, T ], z ′′ ∈ L2 [0, T ], and z (0) = 0}. Whilst, the inner product and the norm of W22 [0, T ] are given,

4

O. Abu Arqub / Computers and Mathematics with Applications (

)

–

respectively, as

⟨z1 (x), z2 (x)⟩W 2 = 2

1 

(i)

(i)

z1 (0) z2 (0) +

T



z1′′ (t ) z2′′ (t ) dt ,

(9)

0

i =0

and ∥z1 ∥2W 2 = ⟨z1 (t ) , z1 (t )⟩W 2 , where z1 , z2 ∈ W22 [0, T ]. 2

2

Theorem 2 ([41]). The Hilbert space W22 [0, T ] is a complete reproducing kernel with the reproducing kernel function

R{s2} (t ) =

 1 1  st + st 2 − t 3 ,

t ≤ s,

 st + s2 t − s3 ,

t > s.

2 1

6 1

2

6

(10)

Definition 3. The space W23 [0, L] is defined as W23 [0, L] = {z = z (x) : z , z ′ , z ′′ are one-variable absolutely continuous real-valued functions on [0, L], z ′′′ ∈ L2 [0, L], and z ′ (0) = z ′ (L) = 0}. Whilst, the inner product and the norm of W23 [0, L] are given, respectively, as

⟨z1 (x), z2 (x)⟩W 3 = 2

and ∥z1 ∥2

W23

1 

(i)

(i)

z1 (0) z2 (0) + z1 (L) z2 (L) +

L



z1′′′ (x) z2′′′ (x) dx,

(11)

0

i =0

= ⟨z1 (x) , z1 (x)⟩W 3 , where z1 , z2 ∈ W23 [0, L]. 2

Theorem 3. The Hilbert space W23 [0, L] is a complete reproducing kernel with the reproducing kernel function

{3}

Ry

1

  

5 (x) = 120L 1  

120L5

   (L − y)3 x3 6x2 y2 + 3Lxy (x − 5y) + L2 10y2 − 5xy + x2 ,

x ≤ y,

,

x > y.

(12)

(L − x)3 y 6x2 y2 + 3Lxy (y − 5x) + L 10x2 − 5xy + y  3

 2

 2

Proof. The proof of the completeness and the reproducing property of W23 [0, L] is similar to the proof in [42]. Let us find {3}

out the expression form of Ry (x) in W23 [0, L]. By using the tabular integration formula, one can get



z (x) , R{y3} (x)



W23

=

2 



z (i) (0) ∂xi R{y3} (0) + (−1)i+1 ∂x55−−ii Ry{3} (0)



i =0

 L 2  (−1)i z (i) (L) ∂x55−−ii R{y3} (L) − z (x) ∂x66 Ry{3} (x) dx.

+

(13)

0

i =0

{3}

Hence, if Ry (x) satisfy the following generalized differential equations:

∂xi i R{y3} (0) = ∂xi i R{y3} (L) ,

i = 0, 1, 2,

∂

(x) = −δ (y − x) ,

∂

(y + 0) = ∂

6 {3} R x6 y i {3} R xi y

i {3} R xi y

δ dirac-delta function,

(y − 0) ,

(14)

i = 0, 1, . . . , 4,

∂x55 R{y3} (y + 0) − ∂x55 R{y3} (y − 0) = −1, {3}

{3}

then the expression form of Ry (x) can be written as Ry (x) =

5

i =0

{3}

ai (y) xi , if x ≤ y and Ry (x) =

5 {3}

i =0

bi (y) xi if x > y.

Through the last descriptions and by using Maple 13 software package, the unknown coefficients of Ry (x) can be obtained directly as utilized in Eq. (12). This completes the proof. 

3. Two extended inner product spaces To solve the given Neumann time-fractional PDEs using the RKHSM, we need to introduce two extended inner product spaces H (Ω ) and W (Ω ) and two extended reproducing kernel functions Rt (s) and rt (s). Anyhow, the space W (Ω ) must be satisfied the constraints conditions of the given equation.

O. Abu Arqub / Computers and Mathematics with Applications (

)

–

5

i +j Henceforth and not to conflict unless stated otherwise, we denote the following symbols: Ω = [0, L] ⊗ [0, T ], ∂xi t j =  i i  j j ∂ /∂ x ∂ /∂ t , whenever i, j = 1, 2, and so on. Next, we combined the results in the previous section to generate two

new Hilbert spaces. Definition 4. The inner product Hilbert space W (Ω ) can be defined as W (Ω ) = {u : ∂x22 ∂t u is complete continuous

function in Ω , ∂x33 ∂t22 u ∈ L2 (Ω ), and u (x, 0) = ∂x u (0, t ) = ∂x u (L, t ) = 0}. The inner product and the norm of W (Ω ) are building, respectively, as

⟨u1 (x, t ) , u2 (x, t )⟩W =

1  ⟨∂tjj u1 (x, 0) , ∂tjj u2 (x, 0)⟩W 3 2

j =0

+

 T  1 0

∂ ∂

2 j u t 2 xj 2

(0, t ) ∂ ∂

(0, t ) + ∂

2 u t2 1

(L, t ) ∂

2 u t2 2

(L, t ) dt

j =0 T



 2 j u t 2 xj 1

L



∂x33 ∂t22 u1 (x, t ) ∂x33 ∂t22 u2 (x, t ) dxdt ,

+ 0

0

(15)

and ∥u1 ∥2W = ⟨u1 (x, t ) , u1 (x, t )⟩W , where u1 , u2 ∈ W (Ω ). Theorem 4. The Hilbert space W (Ω ) is a complete reproducing kernel with the reproducing kernel function R(y,s) (x, t ) = R{y3} (x) R{s2} (t ) ,

(16) { }

= u (y, s) and R(y,s) (x, t ) = R(x,t ) (y, s), where Ry3 (x) and {2} Rs (t ) are the reproducing kernel functions of the spaces W23 [0, L] and W22 [0, T ], respectively.

such that for any u (x, t ) ∈ W (Ω ), we have u (x, t ) , R(y,s) (x, t )





W

Proof. By using properties of the inner product spaces W23 [0, L] and W22 [0, T ] with respect to the differentials dx and dt, respectively, one can find

⟨u (x, t ) , R{y3} (x) R{s2} (t )⟩W =

1  ⟨∂tjj u (x, 0) , ∂tjj R{y3} (x) R{s2} (0)⟩W 3 2

i=0 T





+

1 

0

(0, t ) ∂ ∂

0

(L, t ) ∂

2 {3} R t2 y

(L) Rs (t ) dt {2}

1  ⟨∂tjj u (x, 0) , R{y3} (x) ∂tjj Rs{2} (0)⟩W 3 2

T





+

1 

0

∂ ∂

(0, t ) ∂

2 {2} R t2 s

0

0

2 u t2

(L, t ) Ry (L) ∂ {3}

∂tjj ⟨u (x, 0) , R{y3} (x)⟩W 3 ∂tjj Rs{2} (0) 2

T



∂t22 R{s2} (t ) ∂t22

+ 0

 1 

∂xj j u (0, t ) ∂xj j R{y3} (0) + u (L, t ) Ry{3} (L)

j =0



L



∂

+ 0

i =0

(t ) ∂xj Ry (0) + ∂ {3}

∂x33 ∂t22 u (x, t ) ∂x33 R{y3} (x) ∂t22 R{s2} (t ) dxdt ,

i=0

1 

j

L



+ 1 

 2 j u t 2 xj

j =0 T



=

2 u t2

∂x33 ∂t22 u (x, t ) ∂x33 ∂t22 R{y3} (x) R{s2} (t ) dxdt ,

0

i=0

=

(0) Rs (t ) + ∂ {2}

L



+

=

∂ ∂

2 j {3} R t 2 xj y

j =0 T



 2 j u t 2 xj

j

3 u x3

(x, t ) ∂

3 {3} R x3 y

j

(x) dx dt

∂t j u (y, 0) ∂t j Rs (0) + {2}

T

 0

∂t22 Rs{2} (t ) ∂t22 ⟨u (x, t ) , Ry{3} (x)⟩W 3 dt 2

2 {2} R t2 s

(t ) dt

6

O. Abu Arqub / Computers and Mathematics with Applications (

=

1 

∂t22 Rs{2} (t ) ∂t22 u (y, t ) dt

0

i =0

–

T



∂tjj u (y, 0) ∂tjj R{s2} (0) +

)

= ⟨u (y, t ) , R{s2} (t )⟩W 2 = u (y, s) .

(17)

2

Thus, u (x, t ) , R(y,s) (x, t )



R(x,t ) (ξ , ζ ) , R(y,s) (ξ , ζ )





 W

= u (y, s). Whilst on the other hand, R(y,s) (x, t ) = = R(x,t ) (y, s). So, the proof of the theorem is complete. 

W



 R(y,s) (ξ , ζ ) , R(x,t ) (ξ , ζ ) W

=

Definition 5. The inner product Hilbert space H (Ω ) can be defined as H (Ω ) = {u : u is complete continuous function in Ω and ∂x ∂t u ∈ L2 (Ω )}. The inner product and the norm of H (Ω ) are building, respectively, as

⟨u1 (x, t ) , u2 (x, t )⟩H = ⟨u1 (x, 0) , u2 (x, 0)⟩W 1 2   T ∂t u1 (0, t ) ∂t u2 (0, t ) dt + +

T 0

0



L

∂xt2 u1 (x, t ) ∂xt2 u2 (x, t ) dxdt ,

(18)

0

and ∥u1 ∥2W = ⟨u1 (x, t ) , u1 (x, t )⟩H , where u1 , u2 ∈ H (Ω ). Theorem 5. The Hilbert space H (Ω ) is a complete reproducing kernel with the reproducing kernel function r(y,s) (x, t ) =  R{y1} (x) R{s1} (t ) ,

(19)

such that for any u (x, t ) ∈ H (Ω ), we have u (x, t ) , r(y,s) (x, t )



{1}

 H

{1}

= u (y, s) and r(y,s) (x, t ) = r(x,t ) (y, s), where Ry (x) and

 1 [0, L] and W 1 [0, T ], respectively. Rs (t ) are the reproducing kernel functions of spaces W 2 2 Proof. Similar to the proof of Theorem 4.



4. Formulation of the considered time-fractional PDE Throughout the remaining sections, we will focus our constructions and results on the first type only in order not to increase the length of the paper without the loss of generality for the remaining type and its results. Actually, in the same manner, we can employ the RKHSM to construct the analytic and the approximate solutions. Problem formulation is normally the most important part of the process. It is the selection of design variables, constraints, models of the discipline/design, and the appropriate fractional linear operator. Anyhow, in order to apply the RKHSM on Eqs. (1) and (2), we must homogenize the nonhomogeneous constraints conditions by suitable transformations. This normalizing is needed to put the mentioned conditions into the space W (Ω ) to make an agreement feasible between the constructed inner product spaces and the corresponding kernel functions which are on Ω , for the convenience, we still denote the solution of the new equation by u (x, t ). So, let

     1−β  2 1−β [u (x, t )] + δ (x, t , u (x, t )) , ∂t u (x, t ) = Aα 0 Dt1−α + Bβ 0 Dt ∂x2 u (x, t ) − Cα 0 Dt1−α + Dβ 0 Dt

(20)

subject to the following initial and Neumann boundary conditions: u (x, 0) = 0,

∂x u (0, t ) = 0,

(21)

∂x u (L, t ) = 0.

For the conduct of proceedings in the RKHSM formulation on the extended spaces H (Ω ) and W (Ω ), we will now define the fractional differential linear operator F as follows:

F : W (Ω ) → H (Ω ) ,

(22)

such that



1−β

F u (x, t ) := ∂t u (x, t ) − Aα 0 Dt1−α + Bβ 0 Dt



   1−β [u (x, t )] . ∂x22 u (x, t ) + Cα 0 Dt1−α + Dβ 0 Dt

(23)

Thus, the time-fractional PDE to be solved is governed by the following equivalent functional equation:

F u (x, t ) = δ (x, t , u (x, t )) ,

in which u ∈ W (Ω ) and δ ∈ H (Ω ) .

(24)

{(xi , ti )}∞ i =1

Next, we construct an orthogonal function systems of W (Ω ) as follows: choose a countable dense subset in Ω , define ϕi (x, t ) = r(xi ,ti ) (x, t ) and ψi (x, t ) = F ∗ ϕi (x, t ), where F ∗ is the adjoint operator of F and r is the reproducing kernel function of H (Ω ).  ∞ Algorithm 1. The orthonormal function systems ψ i (x, t ) i=1 of W (Ω ) can be derived from the process of the Gram–Schmidt orthogonalization of {ψi (x, t )}∞ i=1 as follows.

O. Abu Arqub / Computers and Mathematics with Applications (

)

–

7

 Step 1: Set cip (t ) = ψi (t ) , ψ p (t ) W . Step 2: For i = 1, 2, . . . and k = 1, 2, . . . , i set 

1 µ11 = ∥ψ1 ∥− W ; 

µii = ∥ψ ∥ − 2 i W

i−1 

−0.5 ci2,p

(t )

;

p=1

 µik = − ∥ψ ∥ − 2 i W

i −1 

−0.5 ci2,p

(t )

i−1 

ci2,p (t ) µpk .

p=k

p=1

Step 3: For i = 1, 2, . . . set

ψ i (x, t ) =

i 

µik ψk (x, t ) .

(25)

k=1

Indeed, in order to apply the RKHSM, we divide the finite domain Ω into a p × q mesh point with the space step size L in the x direction of [0, L] and the time step size 1t = Tq in the t direction of [0, T ], respectively, in which p and q p are positive integers. Anyhow the grid points (xl , tm ) in the space–time domain Ω are defined simultaneously as

1x =

(xl , tm ) = (l1x, m1t ) ,

l = 0, 1, 2, . . . , p, m = 0, 1, 2, . . . , q.

(26)

5. Essential results for the constructed method In this section, we will see what the influence choice of the continuous linear operator F in order to show its

∞ boundedness. ∞After that, we will show the sequence {ψi (x, t )}i=1 is a complete function system, whilst the sequence R(xi ,ti ) (x, t ) i=1 is linearly independent in the space W (Ω ).

Lemma 1. F : W (Ω ) → H (Ω ) is a bounded linear operator. Proof. The linearity part is obvious, for the boundedness part, we need to prove that ∥F u (x, t )∥2H

≤ M ∥u∥2W ,  i j M > 0. By the reproducing property of R(y,s) (x, t ), we have u (x, t ) = u (◦, •) , R(x,t ) (◦, •) W and ∂xi ∂t j F u (x, t ) =   j j u (◦, •) , ∂xi i ∂t j F R(x,t ) (◦, •) , i, j = 0, 1. Thus, from the continuity of R(x,t ) (◦, •) and ∂xi i ∂t j F u (x, t ), and the Schwarz 

W

inequality, one can write

        i j      j j ∂xi ∂t j F u (x, t ) =  u (◦, •) , ∂xi i ∂t j F R(x,t ) (◦, •)  ≤ ∥u∥W ∂xi i ∂t j F R(x,t ) (◦, •)

W

W

≤ Mi,j ∥u∥W .

(27)

From the definition of the inner product and the norm of W (Ω ), we have

∥F u (x, t )∥2H = ⟨F u (x, t ) , F u (x, t )⟩H = ⟨F u (x, 0) , F u (x, 0)⟩W 1 + 2

= (F u (0, 0))2 +

T



(∂t F u (0, t ))2 dt + 0

(∂x F u (x, 0))2 dx + 0

T



L



+ 0

L



0

L



T





2 ∂xt2 F u (x, t ) dxdt

0

T



(∂t F u (0, t ))2 dt 0

 2 2 ∂xt F u (x, t ) dxdt

0

≤ M02,0 ∥u∥2W +

L

 0

M12,0 ∥u∥2W dx +

T

 0

M02,1 ∥u∥2W dt +

T

 0

L

 0

M12,1 ∥u∥2W dxdt

  = M02,0 + M12,0 L + M02,1 T + M12,1 LT ∥u∥2W . In other words, ∥F u (x, t )∥

2 W21

≤ M ∥ u∥

2 W,

where M =

M02,0

+

M12,0 L

(28)

+

M02,1 T

+

M12,1 LT .



∞ Theorem 6. Suppose that {(xi , ti )}∞ i=1 is dense on Ω , then the sequence {ψi (x, t )}i=1 is a complete function system in W (Ω ) with

    ψi (x, t ) = ∂s − Aα 0 Ds1−α + Bβ 0 Ds1−β ∂y22 + Cα 0 Ds1−α + Dβ

1−β 0 Ds



R(y,s) (x, t ) 

 

(y,s)=(xi ,ti )

.

(29)

8

O. Abu Arqub / Computers and Mathematics with Applications (

)

–

Proof. In this proof, the subscript (y, s) by the operator F , denoted by F(y,s) , indicates that the operator F applies to the function of (y, s). Indeed, it is easy to see that.

  ψi (x, t ) = F ∗ ϕi (x, t ) = F ∗ ϕi (y, s) , R(x,t ) (y, s) W   = ϕi (y, s) , F(y,s) R(x,t ) (y, s) H  = F(y,s) R(x,t ) (y, s)(y,s)=(x ,t ) i i  = F(y,s) R(y,s) (x, t )(y,s)=(x ,t ) i i        = ∂s − Aα 0 Ds1−α + Bβ 0 Ds1−β ∂y22 + Cα 0 Ds1−α + Dβ 0 Ds1−β R(y,s) (x, t ) 

(y,s)=(xi ,ti )

.

(30)

Clearly, ψi (x, t ) ∈ W (Ω ). Now, for each fixed u ∈ W (Ω ), let ⟨u (x, t ) , ψi (x, t )⟩W = 0 i = 1, 2, . . .. Then, ⟨u (x, t ) , ψi (x, t )⟩W = ⟨u (x, t ) , F ∗ ϕi (x, t )⟩W = ⟨F u (x, t ) , ϕi (t )⟩H = F u (xi , ti ) = 0. Whilst, {(xi , ti )}∞ i=1 is dense on Ω , we must have F u (x, t ) = 0 from the existence of F −1 , it follows that u = 0. The proof is complete.   ∞ Theorem 7. Suppose that {(xi , ti )}∞ i=1 is dense on Ω . Then the sequence R(xi ,ti ) (x, t ) i=1 is linearly independent in W (Ω ).  m Proof. It is enough to show that R(xi ,ti ) (x, t ) i=1 is linearly independent for each m ≥ 1. In fact, if {ci }m i=1 satisfies  m 1, l = k, for each l = 1, 2, . . . , m, then it follows i=1 ci R(xi ,ti ) (x, t ) = 0, taking hk (x, t ) ∈ W (Ω ) such that hk (xl , tl ) = 0, l ̸= k, that:

 0 =

hk (x, t ) ,

m 

 ci R(xi ,ti ) (x, t )

i=1

=

W

m  

 ci hk (x, t ) , R(xi ,ti ) (x, t ) W

i =1

=

m 

ci hk (xi , ti ) = ck ,

k = 1 , 2 , . . . , m.

(31)

i =1

Hence, {ψi (x, t )}∞ i=1 is linearly independent in W (Ω ). So, the proof of the theorem is complete.



6. Representation of the numerical solutions In this section, we will show how to solve Eqs. (24) and (21) by using the RKHSM in detail. After that, we will present the expansion forms of the exact and the numerical solutions in the spaces H (Ω ) and W (Ω ). Through the remainder sections and not to conflict unless stated otherwise, we will use the following brief symbols: δ = δ (x, t , u (x, t )), δk := δ (xk , tk , u (xk , tk )), and δkn := δ (xk , tk , un (xk , tk )), whenever k = 1, 2, . . . , ∞. Theorem 8. Suppose that {(xi , ti )}∞ i=1 is dense on Ω , if u ∈ W (Ω ) is the solution of Eqs. (24) and (21), then u (x, t ) =

∞  i 

µik δk ψ¯ i (x, t ) .

(32)

i=1 k=1

Proof. Since, ⟨u (x, t ) , ϕi (x, t )⟩W = u (xi , ti ) for each u ∈ W (Ω ), whilst,

∞ i ¯ i=1 k=1 µik δk ψi (x, t ) is the Fourier series  ∞ expansion about ψ i (x, t ) i=1 , then it is a convergent series in the sense of ∥·∥W . Thus, u (x, t ) =

∞  

¯ i (x, t ) u (x, t ) , ψ

 W

ψ¯ i (x, t )

i =1

=

∞ 

 u (x, t ) ,

i =1

=

∞  i 

i  k=1

 µik ψk (x, t )

ψ¯ i (x, t ) W

  µik u (x, t ) , F ∗ ϕk (x, t ) W ψ¯ i (x, t )

i=1 k=1

=

∞  i  i=1 k=1

  µik F u (x, t ) , r(xk ,tk ) (x, t ) H ψ¯ i (x, t )

O. Abu Arqub / Computers and Mathematics with Applications (

∞  i 

=

)

–

9

µik F u (xk , tk ) ψ¯ i (x, t )

i=1 k=1

∞  i 

=

µik δk ψ¯ i (x, t ) .

(33)

i=1 k=1

Therefore, the form of Eq. (32) is the exact solution of Eqs. (24) and (21). This completes the proof. Anyhow, because W (Ω ) is a Hilbert space, the series u (x, t ) = series un (x, t ) =

n  i 

µik δk ψ¯ i (x, t ) ,

n = pq =

L

T

1x 1t

∞ i i=1

k=1



µik δk ψ¯ i (x, t ) < ∞. Therefore, the truncating

,

(34)

i=1 k=1

is convergent in the sense of the norm of W (Ω ) and the numerical solution un (x, t ) of Eqs. (24) and (21) can be calculated directly by Eq. (34). To analyze the most comprehensive computations, the basis of our RKHSM for solving Eqs. (24) and (21) is summarized below for the exact and the numerical solutions depending on the internal structure of the function δ . Case 1: If δ is linear, then the exact and the numerical solutions can be obtained directly from Eqs. (32) and (34), respectively. Case 2: If δ is nonlinear, then the exact and the numerical solutions can be obtained iteratively by using the following process. According to Eq. (32), the representation form of the exact solution can be written as u (x, t ) =

∞ 

¯ i ( x, t ) , Ai ψ

i=1

Ai =

i 

µik δk .

(35)

k=1

For numerical computations, put (x1 , t1 ) = (0, 0), then from the initial condition of Eq. (21), the value of u (x1 , t1 ) is known. Set u0 (x1 , t1 ) = u (x1 , t1 ) and define the n-term numerical solution of u (x, t ) as follows: un (x, t ) =

n 

¯ i (x, t ) , Ai ψ

i=1

Ai =

i 

µik δk , n = pq =

L

T

1x 1t

.

(36)

k=1

Anyhow, in the iterative process of Eq. (36), we can guarantee that the numerical solution un satisfies the constraints conditions of Eq. (21). 7. Convergence analysis of the method Convergence allows approximations to be computed accurately on an individual realization basis. For some usages, such detailed pathwise information is required. In this section, the n-truncation numerical solution is proved to converge uniformly to the exact solution as n → ∞. Lemma 2. If ∥un−1 − u∥W → 0, (xn , tn ) → (y, s) as n → ∞, ∥un−1 ∥W is bounded, and δ is continuous, then δnn−1 → δ as n → ∞. Proof. By the reproducing property of R(y,s) (x, t ), we have u (x, t )



i j R xi t j (x,t )

u (◦, •) , ∂ ∂

(◦, •)

one can find

   i j  ∂xi ∂t j u (x, t ) =

 W

=



u (◦, •) , R(x,t ) (◦, •) i j u xi t j

, i, j = 0, 1. From the continuity of R(x,t ) (◦, •) and ∂ ∂

j

 W

and ∂xi i ∂t j u (x, t )

=

(x, t ), and the Schwarz inequality,

     j  u (◦, •) , ∂xi i ∂t j R(x,t ) (◦, •)  W     ≤ ∥u∥W ∂xi i ∂tjj R(x,t ) (◦, •) ≤ Ni,j ∥u∥W . (37) W     From un−1 (y, s) = un−1 (◦, •) , R(y,s) (◦, •) W , u (y, s) = u (◦, •) , R(y,s) (◦, •) W , and ▽u (◦, •) = ∂x u (◦, •) e1 + ∂t u (◦, •) e2 , where e1 = ⟨1, 0⟩ and e2 ⟨= 0, 1⟩, we get    |un−1 (y, s) − u (y, s)| =  un−1 (x, t ) − u (x, t ) , R(y,s) (x, t ) W    ≤ ∥un−1 − u∥W R(y,s) (x, t )W ≤ N ∥un−1 − u∥W . (38)

10

O. Abu Arqub / Computers and Mathematics with Applications (

|▽un−1 (ζ , ξ )| = ≤



)

–

(∂x un−1 (ζ , ξ ))2 + (∂t un−1 (ζ , ξ ))2



N12,0 + N02,1 ∥un−1 ∥W .

(39)

On the other aspect as well, by merging Eqs. (38) and (39), yields that

|un−1 (xn , tn ) − u (y, s)| = |un−1 (xn , tn ) − un−1 (y, s) + un−1 (y, s) − u (y, s)| ≤ |un−1 (xn , tn ) − un−1 (y, s)| + |un−1 (y, s) − u (y, s)| ≤ |▽un−1 (ζ , ξ )| |(xn , tn ) − (y, s)| + |un−1 (y, s) − u (y, s)|  ≤ N12,0 + N02,1 ∥un−1 ∥W |(xn , tn ) − (y, s)| + N ∥un−1 − u∥W ,

(40)

where ζ lies between xn , y and ξ lie between tn , s. From ∥un−1 − u∥W → 0 and (xn , tn ) → (y, s) as n → ∞, and the boundedness of ∥un−1 ∥W , it follows that |un−1 (xn , tn ) − u (y, s)| → 0 as n → ∞. As a result, by the means of the continuation of δ , it implies that δnn−1 → δ as n → ∞.  Theorem 9. Suppose that ∥un ∥W is bounded in Eq. (36), {(xi , ti )}∞ i=1 is dense on Ω , and Eqs. (24) and (21) have a unique solution. Then the n-term numerical solution un (x, t ) converges to the exact solution u (x, t ).

¯ n+1 (x, t ). Proof. Firstly, we will prove the convergence of un (x, t ). From Eq. (36), we infer that un+1 (x, t ) = un (x, t )+ An+1 ψ ∞ 2 2 2 2 2 2 ¯ From the orthogonality of ψi (x, t ) i=1 it follows that ∥un+1 ∥W = ∥un ∥W + An+1 = ∥un−1 ∥W + An + An+1 = · · · =

 ∥u0 ∥2W + in=+11 A2i . In other formulation, it holds that ∥un+1 ∥W ≥ ∥un ∥W . Due to the condition that ∥un ∥W is bounded,  2 ∞  2 2 ∥un ∥W is convergent and there exists a positive constant c such that ∞ i=1 Ai = c. This implies that Ai i=1 ∈ l . On the other hand, since (um (t ) − um−1 (t ))⊥(um−1 (t ) − um−2 (t ))⊥ . . . ⊥(un+1 (t ) − un (t )), it follows for ζ > n that: ∥um − un ∥2W = ∥um − um−1 + um−1 − · · · + un+1 − un ∥2W = ∥um − um−1 ∥2 + ∥um−1 − um−2 ∥2 + · · · + ∥un+1 − un ∥2W .

(41)

2 Furthermore, ∥um − um−1 ∥2W = A2m . Consequently, as n, m → ∞, we have ∥um − un ∥2W = l=n+1 Al → 0. Considering the completeness of W (Ω ), there exist u ∈ W (Ω ) such that un (x, t ) → u (x, t ) as n → ∞. Secondly, solution of Eqs. (24) and (21). Taking the limits in Eq. (36), one can get ∞we will prove that u (x, t ) is the  ∞ ¯ ¯ u (x, t ) = A ψ x , t . Because F u x , t = ( ) ( ) i i i =1 i=1 Ai Lψi (x, t ), we get,

m

F u (xl , tl ) =

∞ 

¯ i (x, t ) , ϕl (x, t ) Ai F ψ 

 H

i =1

=

∞ 

¯ i (x, t ) , F ∗ ϕl (x, t ) Ai ψ 

 W

i =1

=

∞ 

¯ i (x, t ) , ψl (x, t ) Ai ψ 

 W

.

(42)

i =1

Multiplying both sides of Eq. (42) by µil and taking the summation i 

µil F u(xl , tl′ ) =

l =1

∞ 

 ¯ i ( x, t ) , Ai ψ

i =1

=

∞ 

i 

i

l =1 ,

it follows that:

 µil ψl (x, t )

l =1

W

 ¯ i (x, t ) , ψ¯ l (x, t ) = Al . Ai ψ W 

(43)

i =1

In view of Eq. (36), we have F u (xl , tl ) = δll−1 . For the conduct of proceedings in the proof, since {(xi , ti )}∞ i=1 is dense  ∞   on Ω , then for each (y, s) ∈ Ω , there exists a subsequence xnj , tnj , such that xnj , tnj → (y, s) as j → ∞. j =1

But since, we have known that F u xnj , tnj





n −1

= δnjj

. Hence, let j → ∞, by the continuity of δ and Lemma 2, we have

¯ i (x, t ) ∈ W , then u (x, t ) satisfies the constraints F u (y, s) = δ (y, s, u (y, s)). Hence, u (x, t ) satisfies Eq. (24). Also, since ψ conditions of Eq. (21). So, the proof of the theorem is complete.  j

Theorem 10. The partial derivatives of the numerical solution ∂tii ∂xj un (x, t ) are converging uniformly to the partial derivatives i j u t i xj

of the exact solution ∂ ∂

(x, t ), whenever i = 0, 1 and j = 0, 1, 2 as n → ∞.

O. Abu Arqub / Computers and Mathematics with Applications (

)

–

11

Proof. Since W (Ω ) is a Hilbert space, from Eqs. (32) and (34), it follows that, ∥u − un ∥W → 0 as n → ∞. Again, since

      i j    j j ∂t i ∂xj u (x, t ) − ∂tii ∂xj un (x, t ) =  u (y, s) − un (y, s) , ∂tii ∂xj F R(x,t ) (y, s)  W     ≤ ∥u − un ∥W ∂tii ∂xj j F R(x,t ) (y, s) ≤ Mi,j ∥u − un ∥W . (44) W     j j j j Thus, ∂tii ∂xj u (x, t ) − ∂tii ∂xj un (x, t ) → 0 as n → ∞ or ∂tii ∂xj un (x, t ) converge to ∂tii ∂xj u (x, t ) uniformly. This completes the proof.



For the error behavior, if εn (x, t ) = u (x, t ) − un (x, t ), then using Eqs. (32) and (34), one can write ∥εn ∥2W

∞

i=n+1

Because



i =1

2



2

∥ = i=n+2 µik δk . Clearly, εn+1 ≥ ¯ k=1 µik δk ψi (x, t ) is convergent series in ∥·∥W , then ∥εn ∥W → 0 as n → ∞.

µik δk ∞ i i k=1

and ∥ε

∞

2 n +1 W

i k=1

εn , and consequently {εn }∞ n =1

=

is decreasing.

8. Numerical algorithm and numerical outcomes To demonstrate the simplicity and effectiveness of the RKHSM, numerical solutions for some certain classes of Neumann time-fractional PDEs with respect to the Riemann–Liouville and the Caputo fractional derivatives are constructed in this section. The results reveal that the method is highly accurate, rapidly converge, and convenient to handle various physical and engineering problems in fractional calculus. To allocating more, the solvability analysis of Neumann time-fractional PDEs presented in this work with parameters derivative in the sense of the Riemann–Liouville and the Caputo fractional derivatives have been discussed by many authors, using different numerical and approximate methods. To mention a few, in [14] the authors have discussed the compact finite difference method for solving the time-fractional heat equation, whilst, the compact finite difference with Laplace transform method is carried out in [15] for solving the same equation. The implicit compact difference method has been applied to solve the time-fractional cable equation as described in [17]. The implicit difference method has been applied to solve the timefractional modified anomalous subdiffusion equation as utilized in [18], in the meantime, the compact difference method is carried out in [19] for solving the same equation. In [20] the authors have been discussed the compact difference method for solving the time-fractional reaction subdiffusion equation, whilst, the finite difference/element approaches are carried out in [21] for solving the same equation. In [26] the authors have been used the variational iteration method and the Adomian decomposition method for solving the time-fractional Fokker–Planck equation, in the meantime, the Galerkin finite element method is carried out in [27] for solving the same equation. The Adomian decomposition method has been used to solve the time-fractional Fisher’s equation as presented in [28], whilst, the quadratic spline functions method is carried out in [29] for solving the same equation. In [31] the authors have proposed the homotopy analysis method for solving the time-fractional Newell–Whitehead equation. Algorithm 2. To approximate the solution un (x, t ) of u (x, t ) for Eqs. (24) and (21), we do the following steps: Step 1: Choose n = pq =

L

T

1x 1t

collocation points in the finite domain Ω ;

Step 2: Set ψi (xi , ti ) = F(y,s) R(y,s) (xi , ti )(y,s)=(x ,t ) ; Step Step Step Step Step

3: 4: 5: 6: 7:



i i

Obtain the orthogonalization coefficients µik using Algorithm 1; i Set ψ i (xi , ti ) = k=1 µik ψk (xi , ti ) for i = 1, 2, . . . , n; Choose an initial approximationu0 (x1 , t1 ); Set i = 1; i k−1 Set Ai = ; k=1 µik δk

¯ Step 8: Set ui (x, t ) = k=1 Ak ψk (xk , tk ); Step 9: If i < n, then set i = i + 1 and go to step 7, else stop. i

Using RKHSM, taking (xl , tm ) =





l pL , m Tq , l = 0, 1, 2, . . . , p, m = 0, 1, 2, . . . , q as a two-dimensional partition of

the domain Ω in un (xl , tm ) of Eq. (36), generating the reproducing kernel functions R(y,s) (x, t ) , r(y,s) (x, t ), and applying Algorithm 2 throughout the numerical computations, some graphical results, tabulate data, and numerical comparison are presented and discussed quantitatively at some selected grid points on Ω to illustrate the numerical solutions for the following Neumann time-fractional PDEs. In the process of computation, all the symbolic and numerical computations are performed by using Maple 13 software package. 8.1. Numerical solution of time-fractional heat equation Here, we present an example on fractional heat PDE on a finite domain Ω = [0, 1] ⊗ [0, 1] to experiment the efficiency and applicability of the RKHSM.

12

O. Abu Arqub / Computers and Mathematics with Applications (

)

–

Table 1 Absolute errors of approximating the solution in Example 1 using RKHSM. t

α = 0.7

0

x

α = 0.8

α = 0.9

α=1

0.25 0.5 0.75 1

−4

2.30525887 × 10 2.72190425 × 10−4 4.49405731 × 10−4 7.92880972 × 10−4

−6

3.02968747 × 10 4.86353611 × 10−5 5.62195895 × 10−5 6.43749963 × 10−5

−7

3.52378874 × 10 5.92724581 × 10−6 5.65878716 × 10−5 4.44093021 × 10−6

1.33937295 × 10−7 3.97161323 × 10−7 2.50871339 × 10−6 1.05895446 × 10−6

0.25

0.25 0.5 0.75 1

9.36658357 × 10−4 1.12537228 × 10−4 1.74271330 × 10−3 5.67127464 × 10−4

1.66003981 × 10−6 1.83934011 × 10−5 3.78552002 × 10−4 2.61389747 × 10−5

2.34224182 × 10−7 3.44578276 × 10−6 1.54642236 × 10−5 5.58617092 × 10−6

1.94333904 × 10−7 1.58613144 × 10−7 2.95074593 × 10−6 4.42962295 × 10−7

0.5

0.25 0.5 0.75 1

9.89807447 × 10−4 3.34607969 × 10−4 3.98483283 × 10−4 7.93340963 × 10−4

2.52006257 × 10−6 3.55636654 × 10−5 2.82845905 × 10−5 5.24832441 × 10−5

8.07961372 × 10−6 5.86970786 × 10−6 2.14574307 × 10−5 2.60316764 × 10−5

1.54843758 × 10−7 5.54247608 × 10−6 2.23875734 × 10−6 6.71016353 × 10−6

0.75

0.25 0.5 0.75 1

1.18018508 × 10−4 2.67385286 × 10−4 1.81472957 × 10−4 3.23475156 × 10−3

9.39044763 × 10−5 3.82845905 × 10−5 3.08075152 × 10−5 2.82845905 × 10−4

4.56096589 × 10−6 6.67169716 × 10−6 2.93774914 × 10−5 2.09543496 × 10−5

6.37433259 × 10−7 4.36180144 × 10−7 6.55401471 × 10−6 2.91220100 × 10−6

1

0.25 0.5 0.75 1

1.98164203 × 10−3 5.60468393 × 10−3 2.09142267 × 10−4 1.00359011 × 10−4

4.56114192 × 10−4 2.21521908 × 10−4 2.42257454 × 10−5 8.30328245 × 10−5

2.97780416 × 10−6 6.26208443 × 10−6 2.34696009 × 10−5 3.10862446 × 10−6

5.84300796 × 10−7 2.41151316 × 10−7 1.84027293 × 10−6 9.00661661 × 10−7

Example 1. Consider the following linear time-fractional PDE:

  ∂t u (x, t ) = 0 Dt1−α ∂x22 u (x, t ) + g (x, t ) ,

0 ≤ x ≤ 1, 0 ≤ t ≤ 1, 0 < α < 1,

(45)

subject to the following initial and Neumann boundary conditions: u (x, 0) = 0,

∂x u (0, t ) = 0, where g (x, t )

(46)

∂x u (1, t ) = 0,

 = ex (2 + α) x2 (1 − x)2 t 1+α −

u (x, t ) = ex x2 (1 − x)2 t α+2 .

Γ (3+α) Γ (2+2α)





2 − 8x + x2 + 6x3 + x4 t 2α+1 . Here, the exact solution is



In order to demonstrate the agreement between the exact and the RKHSM approximate solutions, Table 1 shows the absolute error of approximate solution of Example 1 obtained at various values of (x, t ) ∈ Ω when α ∈ {0.7, 0.8, 0.9, 1}. From the table, it can be seen that the numerical solutions at each level characteristics α are in good agreement with the exact solutions. The geometric behaviors of the memory and hereditary properties of the RKHSM approximate solutions and their level characteristics are studied next. Anyhow, the comparisons of between the computational values of the RKHSM approximate solutions for different values of α for Example 1 have been depicted on the domain Ω as shown in Fig. 1. It is clear from Fig. 1 that each of the graphs are nearly coinciding and similar in their behaviors with good agreement with the RKHSM approximate solutions when the ordinary derivatives are considered, whilst, each of the graphs are nearly identical and in excellent agreement to each other in terms of the accuracy. As a result, one can note that the RKHSM approximate solutions continuously depend on the time-fractional derivatives assumed. 8.2. Numerical solution of time-fractional cable equation Here, we present an example on fractional cable PDE on a finite domain Ω = [0, 1] ⊗ [0, 1] to experiment the efficiency and applicability of the RKHSM. Example 2. Consider the following nonlinear time-fractional PDE:

  1−β ∂t u (x, t ) = 0 Dt1−α ∂x22 u (x, t ) − 0 Dt [u (x, t )] + u3 (x, t ) − u (x, t ) + g (x, t ) ,

0 ≤ x ≤ 1, 0 ≤ t ≤ 1, 0 < α, β < 1,

(47)

subject to the following initial and Neumann boundary conditions: u (x, 0) = 0,

∂x u (0, t ) = 0,

∂x u (1, t ) = 0,

(48)

O. Abu Arqub / Computers and Mathematics with Applications (

)

–

13

Fig. 1. Comparisons of between the computational values of the RKHSM approximate solutions at various values of (x, t ) ∈ Ω for Example 1: black α = 1; blue: α = 0.75; green: α = 0.5; red: α = 0.25. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Table 2 Absolute errors of approximating the solution in Example 2 using RKHSM. t

(α = 0.7, β = 0.7)

(α = 0.8, β = 0.8)

(α = 0.9, β = 0.9)

(α = 1, β = 1)

0

0.25 0.5 0.75 1

−4

2.46164796 × 10 5.73112761 × 10−4 2.67360782 × 10−4 1.56545297 × 10−3

−5

8.12198171 × 10 1.83064095 × 10−4 2.34560072 × 10−4 7.71487108 × 10−4

−6

1.51768545 × 10 7.43294317 × 10−5 4.61142359 × 10−5 4.88461269 × 10−5

7.35432483 × 10−6 1.33898771 × 10−5 3.63443697 × 10−5 2.01440016 × 10−5

0.25

0.25 0.5 0.75 1

1.05426574 × 10−4 4.35954614 × 10−4 4.46691295 × 10−4 8.90802326 × 10−4

8.28629496 × 10−5 9.01175766 × 10−4 1.72527559 × 10−4 5.03151727 × 10−4

5.36363310 × 10−6 2.04941737 × 10−5 4.40964112 × 10−5 1.06450476 × 10−5

6.43828074 × 10−6 1.25535616 × 10−5 1.11046251 × 10−5 1.75425020 × 10−5

0.5

0.25 0.5 0.75 1

1.63909021 × 10−3 5.50183630 × 10−4 1.73060738 × 10−3 9.70734024 × 10−4

7.79317206 × 10−5 2.34267929 × 10−4 6.57304469 × 10−4 4.09312209 × 10−4

7.22234737 × 10−5 5.35533905 × 10−5 3.26231207 × 10−5 4.86021988 × 10−5

5.91574113 × 10−6 3.76435102 × 10−5 2.62099213 × 10−5 3.96861427 × 10−5

0.75

0.25 0.5 0.75 1

1.19767356 × 10−3 8.93465439 × 10−4 1.01687632 × 10−3 5.62183171 × 10−4

9.16711576 × 10−4 4.34534562 × 10−5 1.32674872 × 10−4 1.64536688 × 10−4

1.56449018 × 10−5 3.82660465 × 10−6 2.91695180 × 10−5 2.45760742 × 10−5

7.23968800 × 10−5 3.64156054 × 10−6 4.40183937 × 10−5 1.92976281 × 10−5

1

0.25 0.5 0.75 1

7.81353986 × 10−4 2.94909674 × 10−3 1.53166722 × 10−4 2.71410028 × 10−4

5.03245856 × 10−5 5.73332234 × 10−4 3.75234064 × 10−4 2.74570789 × 10−4

3.28511520 × 10−5 6.20878250 × 10−5 3.69804944 × 10−5 9.54414697 × 10−5

5.82093136 × 10−6 3.15722807 × 10−5 9.32427167 × 10−5 4.97696223 × 10−5

x



Γ (2+α+β)

where g (x, t ) = sin (3π (1 − x)) t α+β (1 + α + β) + 9π 2 Γ (1+2α+β) t α +

Γ (2+α+β) β t Γ (1+α+2β)

 − t 3+2(α+β) sin2 (3π (1 − x)) + t .

Here, the exact solution is u (x, t ) = t 1+α+β sin (3π (1 − x)). Our next goal is to illustrate some numerical results of the RKHSM approximate solutions of the aforementioned timefractional cable equation in numeric values. Anyhow, the agreement between the exact and the numerical solutions is investigated for Example 2 at various values of (x, t ) ∈ Ω when α, β ∈ {0.7, 0.8, 0.9, 1} by computing the numerical approximation of their exact solutions for the corresponding equivalent fractional equations as shown in Table 2.

14

O. Abu Arqub / Computers and Mathematics with Applications (

)

–

Fig. 2. Comparisons of between the computational values of the RKHSM approximate solutions at various values of (x, t ) ∈ Ω for Example 2: black (α, β) = (1, 1); blue: (α, β) = (0.75, 0.5); green: (α, β) = (0.75, 0.25); red: (α, β) = (0.5, 0.25). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Again, to further show the advantage of the memory and hereditary properties of the RKHSM approximate solutions and their level characteristics, the comparisons of between the computational values of the RKHSM approximate solutions for different values of α, β ∈ {1, 0.75, 0.5, 0.25} for Example 2 have been depicted on the domain Ω as shown in Fig. 2. As a result, the process described by the RKHSM is slightly more skewed to the down than that modeled by the standard cable equation. Whilst, from the graphs, it can be seen that, the RKHSM approximate solutions are stable, convergent, and continuously depend on the time-fractional derivatives used for solving Eqs. (47) and (48). 8.3. Numerical solution of time-fractional modified anomalous subdiffusion equation Here, we utilized an example on fractional modified anomalous subdiffusion PDE on a finite domain Ω = [0, 1] ⊗ [0, 1] to experiment the efficiency and applicability of the RKHSM. Example 3. Consider the following linear time-fractional PDE:

   1−β  2 ∂t u (x, t ) = 0.5 0 Dt1−α + 0.5 0 Dt ∂x2 u (x, t ) + g (x, t ) ,

0 ≤ x ≤ 1, 0 ≤ t ≤ 1, 0 < α, β < 1,

(49)

subject to the following initial and Neumann boundary conditions: u (x, 0) = 0,

    ∂x u (0, t ) = t t α + t β , ∂x u (1, t ) = et t α + t β ,    Γ (2+α) where g (x, t ) = ex (1 + α) t α − Γ (1+2α) t 2α + ex (1 + β) t β −   ex t t α + t β .

(50) Γ (2+β) 2β t Γ (1+2β)



. Here, the exact solution is u (x, t ) =

Now, the influence of memory and hereditary properties of the RKHSM approximate solutions and their level characteristics on the corresponding errors are investigated for different values of α and β . Table 3 gives the relevant data for Example 3 at the middle space x direction and the middle time t direction, whilst Table 4 gives the relevant data at the end space x direction and the end time t direction. Regarding the convergence speed, it is obvious that the difference between the exact and the RKHSM nodal values decreases initially till a maximum level fractional derivative values are reached.

O. Abu Arqub / Computers and Mathematics with Applications (

)

–

15

Table 3 Absolute errors of approximating the solution u (0.5, 0.5) in Example 3 using RKHSM.

(α, β)

0.1

0.2 0.4 0.6 0.8 1

5.50428 × 10 2.80147 × 10−4 6.24043 × 10−4 1.38723 × 10−5 3.11749 × 10−5

0.3 −4

0.5

8.80327 × 10 1.98010 × 10−4 1.41135 × 10−5 9.98154 × 10−5 1.98996 × 10−5

0.7

1.02952 × 10 3.86054 × 10−5 2.87845 × 10−5 3.27808 × 10−5 4.09574 × 10−5

−4

−5

0.9

1

1.41789 × 10 3.27026 × 10−5 6.76824 × 10−5 4.51141 × 10−5 2.69578 × 10−5

1.71758 × 10 2.66027 × 10−5 5.52043 × 10−5 4.70036 × 10−6 5.15447 × 10−6

4.00153 × 10−5 2.93697 × 10−5 3.72442 × 10−5 2.54442 × 10−6 2.91630 × 10−6

0.7

0.9

1

−5

−5

Table 4 Absolute errors of approximating the solution u (1, 1) in Example 3 using RKHSM.

(α, β)

0.1

0.2 0.4 0.6 0.8 1

4.34684 × 10 2.76215 × 10−3 3.73522 × 10−4 1.21602 × 10−3 2.04563 × 10−5

0.3 −3

0.5

6.17377 × 10 1.05346 × 10−3 6.17754 × 10−4 3.46512 × 10−3 2.48518 × 10−5

7.00658 × 10 5.73371 × 10−4 4.62827 × 10−4 3.45841 × 10−4 3.99991 × 10−5

−4

1.81446 × 10 6.49945 × 10−5 4.86235 × 10−4 1.32759 × 10−5 2.24625 × 10−5

−4

−4

6.43268 × 10 3.51936 × 10−5 9.87317 × 10−5 6.31582 × 10−6 1.43924 × 10−5

−4

8.09878 × 10−4 2.57169 × 10−5 1.18442 × 10−5 7.43064 × 10−5 1.48902 × 10−6

Table 5 Absolute errors of approximating the solution in Example 4 using RKHSM. t

α = 0.7

0

x

α = 0.8

α = 0.9

α=1

0.25 0.5 0.75 1

−5

3.09071699 × 10 5.87785252 × 10−4 8.09016994 × 10−5 1.51056516 × 10−4

−6

3.66720799 × 10 3.10433268 × 10−5 4.07687904 × 10−6 3.45101525 × 10−5

−6

2.78766607 × 10 9.14441244 × 10−6 2.78513106 × 10−5 9.99988656 × 10−6

4.73459555 × 10−7 8.03729722 × 10−6 3.86053602 × 10−6 5.90093459 × 10−6

0.25

0.25 0.5 0.75 1

3.23108510 × 10−5 2.99462917 × 10−4 1.32045983 × 10−4 9.98446542 × 10−5

2.14208946 × 10−6 2.64539396 × 10−5 1.27063506 × 10−5 1.95274968 × 10−5

6.01300398 × 10−6 8.39556673 × 10−6 9.24694579 × 10−6 1.67147802 × 10−5

6.34429292 × 10−7 5.27502327 × 10−6 6.43737879 × 10−6 6.84057432 × 10−6

0.5

0.25 0.5 0.75 1

6.34737175 × 10−5 2.57011621 × 10−4 6.85297043 × 10−4 1.90685894 × 10−4

2.36489808 × 10−6 1.44664795 × 10−5 2.50882946 × 10−5 8.00078170 × 10−5

1.24211576 × 10−6 1.03234526 × 10−5 2.82647062 × 10−5 5.18446732 × 10−6

5.93281102 × 10−7 6.41442955 × 10−6 2.37693411 × 10−6 7.63105954 × 10−7

0.75

0.25 0.5 0.75 1

2.72032600 × 10−5 2.49613574 × 10−4 1.46718025 × 10−4 7.09335263 × 10−4

6.57530599 × 10−5 7.97256221 × 10−6 1.93599327 × 10−5 1.64383548 × 10−5

3.13151289 × 10−6 4.80400646 × 10−5 3.37042362 × 10−5 5.73822739 × 10−6

4.38986266 × 10−7 3.92699159 × 10−6 4.27905261 × 10−6 2.63910727 × 10−6

1

0.25 0.5 0.75 1

2.48925371 × 10−5 1.35533459 × 10−4 2.14053150 × 10−4 6.47018812 × 10−5

9.85746163 × 10−6 1.69433614 × 10−5 1.57032448 × 10−5 3.53089079 × 10−6

1.42538367 × 10−6 1.61395757 × 10−6 2.07980688 × 10−5 6.28809306 × 10−6

2.67527898 × 10−7 2.28146241 × 10−6 1.93982155 × 10−6 1.78987732 × 10−6

8.4. Numerical solution of time-fractional reaction subdiffusion equation In this subsection, we present one example on fractional reaction subdiffusion PDE on a finite domain Ω = [0, 1] ⊗ [0, 1] to experiment the efficiency and applicability of the RKHSM. Example 4. Consider the following linear time-fractional PDE:

∂tαα u (x, t ) = ∂x22 u (x, t ) − u (x, t ) + g (x, t ) ,

0 ≤ x ≤ 1, 0 ≤ t ≤ 1, 0 < α < 1,

(51)

subject to the following initial and Neumann boundary conditions: u (x, 0) = x2 (1 − x)2 ,

∂x u (0, t ) = 0,

∂x u (1, t ) = 0,

(52)

  (3 + α) x2 (1 − x)2 (t + 1)2 + (t + 1)α+2 x4 − 2x3 − 11x2 + 12x − 2 . Here, the exact solution is u (x, t ) = x2 (1 − x) (t + 1)α+2 .

where g (x, t ) =

1 Γ 2 2

We solve this problem with the scheme presented in the current paper with different values of α . Anyhow, the computational results of the scheme at various values of (x, t ) ∈ Ω when α ∈ {0.7, 0.8, 0.9, 1} are reported in Table 5. It is clear from the table that the proposed method produces results with good accuracy and applicability.

16

O. Abu Arqub / Computers and Mathematics with Applications (

)

–

8.5. Numerical solution of time-fractional Fokker–Planck equation In this subsection, we present one example on fractional Fokker–Planck PDE on a finite domain Ω = [0, 1] ⊗ [0, 1] to experiment the efficiency and applicability of the RKHSM. Example 5. Consider the following nonlinear time-fractional PDE:



α



1

∂t α u (x, t ) = ∂ (xu (x, t )) − ∂x 3u (x, t ) − x 2 x2



2

[u (x, t )] ,

0 ≤ x ≤ 1, 0 ≤ t ≤ 1, 0 < α < 1,

(53)

subject to the following initial and Neumann boundary conditions: u (x, 0) = 0,

∂x u (0, t ) = Eα (t α ) ,

(54)

∂x u (1, t ) = Eα (t α ) .

Here, the exact solution is u (x, t ) = xEα (t α ), where Eα (z ) =

zk k=0 Γ (1+α k)

∞

is the Mittag-Leffler function.

From the key feature of the RKHSM, it is possible to pick any point in domain Ω and as well the approximate partial derivatives ∂xυυ ∂t un (x, t ), whenever υ, ω = 1, 2 will be applicable. Next, numerical results of approximating these derivatives for Example 5 at various values of (x, t ) ∈ Ω have been tabulated in Table 6. To be specific, the n-term numerical solutions of ∂xυυ ∂t u (x, t ) are computed as follows:

∂xυυ ∂t un (x, t ) =

n  i=1

¯ i (x, t ) , Ai ∂xυυ ∂t ψ

Ai =

i 

µik δk .

(55)

k=1

Further, by Theorem 10, the approximate partial derivatives ∂xυυ ∂t un (x, t ) converge uniformly to exact partial derivatives ∂xυυ ∂t u (x, t ). 8.6. Numerical solution of time-fractional Newell–Whitehead equation In this subsection, we present one example on fractional Newell–Whitehead PDE on a finite domain Ω = [0, 1] ⊗ [0, 1] to experiment the efficiency and applicability of the RKHSM. Example 6. Consider the following nonlinear time-fractional PDE:

  ∂tαα u (x, t ) = ∂x22 u (x, t ) + u (x, t ) 1 − u2 (x, t ) + g (x, t ) ,

0 ≤ x ≤ 1, 0 ≤ t ≤ 1, 0 < α < 1,

(56)

subject to the following initial and Neumann boundary conditions: u (x, 0) = 0,

∂x u (0, t ) = 0, where g (x, t ) t 1+α cos (π x).

=

∂x u (1, t ) = 0,    t cos (π x) Γ (2 + α) + t α π 2 − t α 1 − t 2+2α cos2 (π x) . Here, the exact solution is u (x, t )

(57)

=

The mathematical behaviors of the error function |εn (x, t )| are discussed next geometrically. To be specific, the n-term absolute errors of approximating un (x, t ) is computed as follows:

  n     1+α  ¯ |εn (x, t )| = t cos (π x) − Ai ψi (x, t ) ,   i=1

Ai =

i 

µik δk .

(58)

k=1

Anyhow, Fig. 3 gives the relevant data of the RKHSM results at various values of (x, t ) ∈ Ω when α ∈ {1, 0.9, 0.8, 0.7} for Example 6. Regarding the convergence speed, it is obvious that the difference between the exact and the RKHSM nodal values decreases initially till a maximum level fractional derivative values are reached. As the plots show, while, the value of x moving a way from the boundaries of [0, 1] and the value of t moving a way from the left boundary of [0, 1], the values of |εn (x, t )| vary along the xt-plane by satisfying the initial condition of Eq. (57). We recall that the accuracy and duration of a simulation depend directly on the size of the steps taken by the solver. Generally, decreasing the step size increases the accuracy of the results, while increasing the time required to simulate the problem. 9. Concluding remarks Neumann time-fractional PDEs provide accessible mathematical models that combine deterministic and fractional order derivative components of dynamic behavior. This article is presented in the RKHSM as a novel solver for some certain classes of PDEs which are (heat, cable, anomalous subdiffusion, reaction subdiffusion, Fokker–Planck, Fisher’s,

O. Abu Arqub / Computers and Mathematics with Applications (

)

–

17

Table 6 Absolute errors of approximating the solution in Example 5 using RKHSM. t

∂t u (x, t )

0

x

∂ x u (x , t )

∂x22 u (x, t )

0.25 0.5 0.75 1

−5

7.15413268 × 10 2.72160696 × 10−5 1.15774442 × 10−4 1.18937141 × 10−4

−5

4.56519640 × 10 1.48688216 × 10−4 1.00587717 × 10−4 3.49507964 × 10−5

2.02863751 × 10−3 6.96464466 × 10−4 3.21693125 × 10−4 2.75191396 × 10−3

0.25

0.25 0.5 0.75 1

2.79638034 × 10−4 1.31104901 × 10−4 1.47014414 × 10−4 9.80579618 × 10−5

1.48285745 × 10−4 8.94006404 × 10−5 6.69502657 × 10−5 2.36806586 × 10−4

8.62219172 × 10−4 5.41236537 × 10−4 1.63349127 × 10−3 4.28685224 × 10−4

0.5

0.25 0.5 0.75 1

1.95036373 × 10−4 5.35175370 × 10−5 1.47014414 × 10−4 5.64303520 × 10−5

5.54766234 × 10−5 6.29427036 × 10−5 2.96806586 × 10−4 8.15758583 × 10−5

1.96671757 × 10−3 1.31272619 × 10−3 8.19234477 × 10−4 2.23863742 × 10−3

0.75

0.25 0.5 0.75 1

2.13253789 × 10−4 8.94426328 × 10−5 3.51168422 × 10−5 2.90880404 × 10−4

2.68063742 × 10−4 7.04589911 × 10−5 9.82088752 × 10−5 6.69502657 × 10−5

4.73511765 × 10−4 1.47014831 × 10−3 7.37672476 × 10−4 2.13224807 × 10−3

1

0.25 0.5 0.75 1

3.23200449 × 10−5 7.21527502 × 10−5 3.99155685 × 10−5 6.88094430 × 10−5

5.42552391 × 10−5 7.75036763 × 10−4 1.12462516 × 10−4 1.05642741 × 10−4

6.01453658 × 10−4 2.90501714 × 10−3 1.96691757 × 10−3 1.14233215 × 10−3

Fig. 3. The absolute value of the error function |εn (x, t )| of the RKHSM approximate solutions for Example 6 when: (a) α = 1, (b) α = 0.9, (c) α = 0.8, and (d) α = 0.7.

and Newell–Whitehead) equations to enlarge its applications range. The solutions are fractional processes that represent diffusive dynamics, a common modeling assumption in many application areas. We included a description of fundamental RKHSM and the concepts of convergence, and error behavior for the RKHSM solvers. Results obtained by the proposed method are found outperform in terms of accuracy, generality, and applicability. It is worth to be pointed out that the RKHSM is still suitable and can be employed for solving other strongly linear and nonlinear time-fractional PDEs. Acknowledgments The author would like to express his gratitude to the unknown referees for carefully reading the paper and their helpful comments.

18

O. Abu Arqub / Computers and Mathematics with Applications (

)

–

References [1] K.B. Oldham, J. Spanier, The Fractional Calculus: Theory and Applications of Differentiation and Integration to Arbitrary Order, Dover Publications, USA, 2006. [2] B. Guo, X. Pu, F. Huang, Fractional Partial Differential Equations and their Numerical Solutions, World Scientific Publishing, Singapore, 2015. [3] S. Das, Functional Fractional Calculus for System Identification and Controls, Springer, USA, 2008. [4] A.A.A. Kilbas, H.M. Srivastava, J.J. Trujillo, Theory and Applications of Fractional Differential Equations, Elsevier, USA, 2006. [5] R. Hilfer, Applications of Fractional Calculus in Physics, World Scientific, USA, 2000. [6] A. El-Ajou, O. Abu Arqub, S. Momani, Approximate analytical solution of the nonlinear fractional KdV-Burgers equation: A new iterative algorithm, J. Comput. Phys. 293 (2015) 81–95. [7] O. Abu Arqub, A. El-Ajou, S. Momani, Constructing and predicting solitary pattern solutions for nonlinear time-fractional dispersive partial differential equations, J. Comput. Phys. 293 (2015) 385–399. [8] L. Zhao, W. Deng, A series of high-order quasi-compact schemes for space fractional diffusion equations based on the superconvergent approximations for fractional derivatives, Numer. Methods Partial Differential Equations 31 (2015) 1345–1381. [9] S.S. Ray, New exact solutions of nonlinear fractional acoustic wave equations in ultrasound, Comput. Math. Appl. 71 (2016) 859–868. [10] H. Chen, S. Lü, W. Chen, Spectral methods for the time fractional diffusion–wave equation in a semi-infinite channel, Comput. Math. Appl. 71 (2016) 1818–1830. [11] L. Yana, F. Yang, The method of approximate particular solutions for the time-fractional diffusion equation with a non-local boundary condition, Comput. Math. Appl. 70 (2015) 254–264. [12] A. El-Ajou, O. Abu Arqub, S. Momani, D. Baleanu, A. Alsaedi, A novel expansion iterative method for solving linear partial differential equations of fractional order, Appl. Math. Comput. 257 (2015) 119–133. [13] M.D. Ortigueira, J.A.T. Machado, Fractional signal processing and applications, Signal Process. 83 (2003) 2285–2286. [14] G.S. Priya, P. Prakash, J.J. Nieto, Z. Kayar, Higher-order numerical scheme for the fractional heat equation with Dirichlet and Neumann boundary conditions, Numer. Heat Transfer B 63 (2013) 540–559. [15] B.A. Jacobs, High-order compact finite difference and Laplace transform method for the solution of time-fractional heat equations with Dirchlet and Neumann boundary conditions, Numer. Methods Partial Differential Equations 32 (2016) 1184–1199. [16] T.A.M. Langlands, B.I. Henry, S.L. Wearne, Fractional cable equation models for anomalous electrodiffusion in nerve cells: finite domain solutions, SIAM J. Appl. Math. 71 (2011) 1168–1203. [17] X. Hua, L. Zhanga, Implicit compact difference schemes for the fractional cable equation, Appl. Math. Model. 36 (2012) 4027–4043. [18] F. Liu, C. Yang, K. Burrage, Numerical method and analytical technique of the modified anomalous subdiffusion equation with a nonlinear source term, J. Comput. Appl. Math. 231 (2009) 160–176. [19] H. Ding, C. Li, High-order compact difference schemes for the modified anomalous subdiffusion equation, Numer. Methods Partial Differential Equations 32 (2016) 213–242. [20] J. Cao, C. Li, Y.Q. Chen, Compact difference method for solving the fractional reaction-subdiffusion equation with Neumann boundary value condition, Int. J. Comput. Math. 92 (2015) 167–180. [21] F. Zeng, C. Li, F. Liu, I. Turner, The use of finite difference/element approaches for solving the time-fractional subdiffusion equation, SIAM J. Sci. Comput. 35 (2013) 2976–3000. [22] I. Sokolov, J. Klafter, Field-induced dispersion in subdiffusion, Phys. Rev. Lett. 97 (2006) 140602. [23] B.I. Henry, T.A.M. Langlands, P.L. Straka, Fractional Fokker–Planck equations for subdiffusion with space- and time-dependent forces, Phys. Rev. Lett. 105 (2010) 170602. [24] C. Buet, S. Dellacherie, R. Sentis, Numerical solution of an ionic Fokker–Planck equation with electronic temperature, SIAM J. Numer. Anal. 39 (2011) 1219–1253. [25] V. Palleschi, F. Sarri, G. Marcozzi, M.R. Torquati, Numerical solution of the Fokker–Planck equation: a fast and accurate algorithm, Phys. Lett. A 146 (1990) 378–386. [26] Z. Odibat, S. Momani, Numerical solution of Fokker–Planck equation with space- and time-fractional derivatives, Phys. Lett. A 369 (2007) 349–358. [27] K.N. Le, W. McLean, K. Mustapha, Numerical solution of the time-fractional Fokker–Planck equation with general forcing, SIAM J. Numer. Anal. 54 (2016) 1763–1784. [28] A.M. Wazwaz, A. Gorguis, An analytic study of Fisher’s equation by using Adomian decomposition method, Appl. Math. Comput. 154 (2004) 609–620. [29] T.S. El-Danaf, A.R. Hadhoud, Computational method for solving space fractional Fisher’s nonlinear equation, Math. Methods Appl. Sci. 37 (2014) 657–662. [30] A. Newell, J. Whitehead, Finite bandwidth, finite amplitude convection, J. Fluid Mech. 38 (1969) 279–303. [31] M.S. Mohamed, M. Faisal, T. Rabeaa, An analytic algorithm for the time-space-fractional Newell-Whitehead equation, Int. Rev. Phys. 6 (2012) 337–343. [32] W. Jiang, Y. Lin, Representation of exact solution for the time fractional telegraph equation in the reproducing kernel space, Commun. Nonlinear Sci. Numer. Simul. 16 (2011) 3639–3645. [33] Y. Wang, M. Du, F. Tan, Z. Li, T. Nie, Using reproducing kernel for solving a class of fractional partial differential equation with non-classical conditions, Appl. Math. Comput. 219 (2013) 5918–5925. [34] W. Jiang, Y. Lin, Approximate solution of the fractional advection–dispersion equation, Comput. Phys. Comm. 181 (2010) 557–561. [35] W. Jiang, Z. Chen, A collocation method based on reproducing kernel for a modified anomalous subdiffusion equation, Numer. Methods Partial Differential Equations 30 (2014) 289–300. [36] M. Cui, Y. Lin, Nonlinear Numerical Analysis in the Reproducing Kernel Space, Nova Science, New York, NY, USA, 2009. [37] A. Berlinet, C.T. Agnan, Reproducing Kernel Hilbert Space in Probability and Statistics, Kluwer Academic Publishers, Boston, Mass, USA, 2004. [38] A. Daniel, Reproducing Kernel Spaces and Applications, Springer, Basel, Switzerland, 2003. [39] H.L. Weinert, Reproducing Kernel Hilbert Spaces: Applications in Statistical Signal Processing, Hutchinson Ross, 1982. [40] Y. Lin, M. Cui, L. Yang, Representation of the exact solution for a kind of nonlinear partial differential equations, Appl. Math. Lett. 19 (2006) 808–813. [41] Y. Zhoua, M. Cui, Y. Lin, Numerical algorithm for parabolic problems with non-classical conditions, J. Comput. Appl. Math. 230 (2009) 770–780. [42] L.H. Yang, Y. Lin, Reproducing kernel methods for solving linear initial–boundary-value problems, Electron. J. Differential Equations 2008 (2008) 1–11. [43] W.Y. Wang, B. Han, M. Yamamoto, Inverse heat problem of determining time-dependent source parameter in reproducing kernel space, Nonlinear Anal. RWA 14 (2013) 875–887. [44] F.Z. Geng, M. Cui, A reproducing kernel method for solving nonlocal fractional boundary value problems, Appl. Math. Lett. 25 (2012) 818–823. [45] X. Li, B. Wu, Approximate analytical solutions of nonlocal fractional boundary value problems, Appl. Math. Model. 39 (2015) 1717–1724. [46] O. Abu Arqub, B. Maayah, Solutions of Bagley-Torvik and Painlevé equations of fractional order using iterative reproducing kernel algorithm, Neural Comput. Appl. (2016) 1–15. http://dx.doi.org/10.1007/s00521-016-2484-4. [47] O. Abu Arqub, The reproducing kernel algorithm for handling differential algebraic systems of ordinary differential equations, Math. Methods Appl. Sci. 39 (2016) 4549–4562. [48] O. Abu Arqub, Approximate solutions of DASs with nonclassical boundary conditions using novel reproducing kernel algorithm, Fund. Inform. 146 (2016) 231–254. [49] O. Abu Arqub, M. Al-Smadi, S. Momani, T. Hayat, Numerical solutions of fuzzy differential equations using reproducing kernel Hilbert space method, Soft Comput. 20 (2016) 3283–3302. [50] O. Abu Arqub, M. Al-Smadi, S. Momani, T. Hayat, Application of reproducing kernel algorithm for solving second-order, two-point fuzzy boundary value problems, Soft Comput. (2016) 1–16. http://dx.doi.org/10.1007/s00500-016-2262-3.

O. Abu Arqub / Computers and Mathematics with Applications (

)

–

19

[51] O. Abu Arqub, Adaptation of reproducing kernel algorithm for solving fuzzy Fredholm-Volterra integrodifferential equations, Neural Comput. Appl. (2015) 1–20. http://dx.doi.org/10.1007/s00521-015-2110-x. [52] W. Jiang, Z. Chen, Solving a system of linear Volterra integral equations using the new reproducing kernel method, Appl. Math. Comput. 219 (2013) 10225–10230. [53] L.H. Yang, J.H. Shen, Y. Wang, The reproducing kernel method for solving the system of the linear Volterra integral equations with variable coefficients, J. Comput. Appl. Math. 236 (2012) 2398–2405. [54] O. Abu Arqub, M. Al-Smadi, N. Shawagfeh, Solving Fredholm integro-differential equations using reproducing kernel Hilbert space method, Appl. Math. Comput. 219 (2013) 8938–8948. [55] O. Abu Arqub, M. Al-Smadi, Numerical algorithm for solving two-point, second-order periodic boundary value problems for mixed integro-differential equations, Appl. Math. Comput. 243 (2014) 911–922. [56] F.Z. Geng, S.P. Qian, Reproducing kernel method for singularly perturbed turning point problems having twin boundary layers, Appl. Math. Lett. 26 (2013) 998–1004. [57] F.Z. Geng, S.P. Qian, S. Li, A numerical method for singularly perturbed turning point problems with an interior layer, J. Comput. Appl. Math. 255 (2014) 97–105. [58] F.Z. Geng, S.P. Qian, Modified reproducing kernel method for singularly perturbed boundary value problems with a delay, Appl. Math. Model. 39 (2015) 5592–5597.