Haar wavelet Picard method for fractional nonlinear partial differential equations

Applied Mathematics and Computation 264 (2015) 310–322

Contents lists available at ScienceDirect

Applied Mathematics and Computation journal homepage: www.elsevier.com/locate/amc

Haar wavelet Picard method for fractional nonlinear partial differential equations Umer Saeed a,∗, Mujeeb ur Rehman b a b

Department of Mathematics, Faculty of Basic and Applied Sciences, University of Poonch Rawalakot, Rawalakot, Pakistan School of Natural Sciences, National University of Sciences and Technology, Sector H-12, Islamabad, Pakistan

a r t i c l e

i n f o

Keywords: Fractional differential equations Haar wavelet Operational matrices Picard iteration

a b s t r a c t In this article, we present a solution method for fractional nonlinear partial differential equation. The proposed technique utilizes the Haar wavelets operational matrices method in conjunction with Picard technique. The operational matrices are derived and utilized for the solution of fractional nonlinear partial differential equations. Convergence analysis for the proposed technique has also been given. Numerical examples are provided to illustrate the efficiency and accuracy of the technique. © 2015 Elsevier Inc. All rights reserved.

1. Introduction Fractional calculus is a branch of mathematics which deals with derivatives and integrals of non-integer orders. In recent years, numerous applications of fractional order ordinary and partial differential equations have appeared in many areas of physics and engineering. There have been found a number of works, especially in hereditary solid mechanics and in viscoelasticity theory, where fractional order derivatives are used for a better description of material properties. This is the main advantage of fractional derivatives in comparison with classical integer order models in which such effects are neglected. For most of fractional order differential equations, exact solutions are not known. Therefore different numerical methods have been applied for providing approximate solutions. Some of these techniques include, the Adomian decomposition method (ADM) [1], the homotopy perturbation method (HPM) [2], the variational iteration method (VIM) [3], the generalized differential transform method (DTM) [4] and wavelet methods [5,6]. There are many kinds of wavelets. The simplest orthonormal wavelet with compact support is the Haar wavelet. Haar wavelet is a rectangular function. A good feature of the Haar wavelets is the possibility to integrate them analytically arbitrary times. A drawback of the Haar wavelets is their discontinuity. Since the derivatives do not exist in the points of discontinuity, it is not possible to apply the Haar wavelets for solving differential equation directly. There are two possibilities to handle this situation. One way is to regularize the Haar wavelets with interpolating splines. This approach has been introduced by Cattani [7]. But this approach complicates the solution process and loses the simplicity of the Haar wavelets. The other way is to make use of the integral method, which was proposed by Chen and Hsiao [8]. They expand the highest derivative in the differential equation into Haar series, and lower derivatives are obtained through integrations. A detailed study of Haar wavelet is given in [9]. The wavelet algorithms for solving partial differential equations are based on the Galerkin techniques or on the collocation method. Haar wavelet algorithm is based on the collocation method. Numerical solution of linear fractional partial differential equations by Haar wavelet method has been given in [10]. Convergence analysis of Haar wavelet method is also discussed. Lepik ∗

Corresponding author. Tel.: +92 333 5700144. E-mail address: [email protected] (U. Saeed).

http://dx.doi.org/10.1016/j.amc.2015.04.096 0096-3003/© 2015 Elsevier Inc. All rights reserved.

U. Saeed, M. ur Rehman / Applied Mathematics and Computation 264 (2015) 310–322

311

[11] organized a method by using two dimensional Haar wavelets for solving partial differential equations, he considered the diffusion and Poisson equations as test problems. An efficient numerical method for the solution of nonlinear partial differential equations based on the Haar wavelets approach is introduced by Celik [6,12], and tested in the case of generalized Burgers–Huxley equation and magnetohydrodynamic flow equations. Hariharan et al. [13,14] developed a solution for Fisher’s equation, Nowell– whitehead equation, Cahn–Allen equation, FitzHugh–Nagumo equation, Burger’s equation and the Burgers–Fisher equation by Haar wavelet method. In this paper, we have constructed the operational matrix of the Haar wavelets and developed a method by combining the Haar wavelet method and Picard technique for numerical solutions of nonlinear fractional partial differential equations. We discretize the nonlinear fractional partial differential equation by Picard technique and then convert the obtained discretized equation into a Sylvester equation by the Haar wavelet technique to get the solution. The convergence analysis of the method is also investigated. We consider a parabolic problem with a quadratic term and Burgers equation as test problems. The exact solution is known for these problems, so the comparisons are possible. These problems are solved to show the applicability of the Haar wavelet method with Picard technique. 2. The uniform Haar wavelets The Haar functions contains just one wavelet during some subinterval of time, and remains zero elsewhere and are orthogonal. The uniform Haar wavelets are useful for the treatment of solution of differential equations which has no abrupt behavior. The lth uniform Haar wavelet hl (x), x  [0, 1) is defined as [8]:

⎧ a(l) ≤ x < b(l); ⎨1, hl (x) = −1, b(l) ≤ x < c(l); ⎩ 0, otherwise,

(2.1)

k k+1 j j j , b(l) = k+0.5 where a(l) = m m , c(l) = m , l = 2 + k + 1, j = 0, 1, 2, . . . , J is dilation parameter, m = 2 and k = 0, 1, 2, . . . , 2 − 1 is translation parameter. When k = 0, j = 0 we have l = 2, which is the minimal value of l and the maximal value of l is 2M where M = 2J , J is maximal level of resolution. In particular h1 (x)  χ [0, 1) (x), where χ [0, 1) (x) is characteristic function on interval [0, 1) , is the Haar scaling function. For the uniform Haar wavelet, the wavelet-collocation method is applied. The collocation points for the uniform Haar wavelets are usually taken as xj = j+0.5 2M , where j = 1, 2, . . . , 2M.

2.1. Fractional integral of the uniform Haar wavelets The Riemann–Liouville fractional integral of the Haar scaling function and the uniform Haar wavelets are given as

pα ,1 (x) = Iaα(1)h1 (x) =

1

 (α)



x

a(1)

(x − s)α−1 ds, α > 0,

(2.2)

⎧ x ⎪ ⎪ (x − s)α−1 ds, a(l) ≤ x < b(l); ⎪ ⎪ ⎪ a(l) ⎪ ⎪ ⎪   x 1 ⎨ b(l) (x − s)α−1 ds − (x − s)α−1 ds, b(l) ≤ x < c(l); pα ,l (x) = Iaα hl (x) =  (α) ⎪ a(l) b(l) ⎪ ⎪ ⎪  c(l) ⎪ b(l) ⎪ ⎪ ⎪ ⎩ (x − s)α−1 ds − (x − s)α−1 ds, x ≥ c(l). a(l)

(2.3)

b(l)

Eqs. (2.2) and (2.3) imply (1))α pα ,1 (x) = (x−a  (α +1) ,

and

(2.4)

⎧ ⎨(x − a(l))α , 1 (x − a(l))α − 2(x − b(l))α , pα ,l (x) =  (α + 1) ⎩(x − a(l))α − 2(x − b(l))α + (x − c(l))α ,

a(l) ≤ x < b(l); b(l) ≤ x < c(l); x ≥ c(l).

(2.5)

The uniform grids are useful for the treatment of differential equations whose solutions have smooth behavior. 3. Function approximations and Haar matrices Any function y  L2 [0, 1] can be represented in terms of the non-uniform Haar series as

y(x) =

∞ 

bl hl (x),

(3.1)

l=1

where bl are the Haar wavelet coefficients given by bl =

1 0

y(x)hl (x)dx.

312

U. Saeed, M. ur Rehman / Applied Mathematics and Computation 264 (2015) 310–322

The function y(x) can be approximated by the truncated Haar wavelet series as

y(x) ≈ yM (x) =

2M 

bl hl (x),

l = 2j + k + 1,

j = 0, 1, 2, . . . , J,

k = 0, 1, 2, . . . , 2j − 1.

(3.2)

l=1

The wavelets coefficients bl are determined in such a way that the integral square error E

 E=

1

y(x) −

0

2M 

2 bl hl (x)

dx

(3.3)

l=1

is minimized. In order to find the numerical approximations of a function, we put the Haar wavelets into a discrete form. For this purpose, we utilized the collocation method. The collocation points for the Haar wavelets are taken as xc (i) = i+0.5 2M , where i = 1, 2, . . . , 2M. In discrete form, Eq. (3.3) is written as

EM = x



2M 

y(xc (i)) −

i=1

2M 

2 bl hl (xc (i))

.

l=1

Any function of two variables u(x, t)  L2 [a, b] × [a, b] can be approximated as

u(x, t) ≈

2M 2M  

cl,i hl (x)hi (t) = HT (x)CH(t),

l=1 i=1

where C is 2M × 2M coefficient matrix which can be determined by the inner product cl, i = hl (x), u(x, t), hi (t) 3.1. Haar matrix Chen and Hsiao [8] established an operational matrix for integration via Haar wavelets, and a procedure for applying the matrix to analyze lumped and distributed- parameters dynamic systems is formulated. The highest derivatives appearing in the differential equations are first expanded into Haar series. The lower order derivatives and the solution functions can then be obtained quite easily by using Haar operational matrix of integration. The discrete form of (3.2) is

yM (xc (i)) =

2M l=1

bl hl (xc (i)), i = 1, 2, . . . , 2M.

(3.4)

We can represent (3.4) in vector form as

y = bH,

(3.5)

where b = [b1 b2  b2M ] and y = [y1 y2  y2M ] are 2M dimensional row vectors and

⎡

H2M×2M

h1 (xc (1)) ⎢ h2 (xc (1)) ⎢ =⎢ .. ⎣ .

h2M (xc (1))

h1 (xc (2)) · · · h2 (xc (2)) · · · .. .. . . h2M (xc (2)) · · ·

⎤ h1 (xc (2M)) h2 (xc (2M)) ⎥ ⎥ ⎥. .. ⎦ . h2M (xc (2M))

In particular, for J = 2, we get 2M = 8 and the Haar matrix is given as

H8×8

⎡ 1 1 1 ⎢1 1 1 ⎢ ⎢1 1 −1 ⎢ ⎢0 0 0 =⎢ ⎢1 −1 0 ⎢ ⎢0 0 1 ⎢ ⎣0 0 0 0 0 0

1 1 1 −1 −1 0 0 1 0 0 −1 0 0 1 0 0

⎤ 1 1 1 −1 −1 −1⎥ ⎥ 0 0 0⎥ ⎥ 1 −1 −1⎥ ⎥. 0 0 0⎥ ⎥ 0 0 0⎥ ⎥ −1 0 0⎦ 0 1 −1

The Haar coefficients bl can be determined by matrix inversion

b = yH−1 , where H−1 is the inverse of H. Eq. (3.6) gives the Haar coefficients bl which are used in (3.2) to get the solution y(x).

(3.6)

U. Saeed, M. ur Rehman / Applied Mathematics and Computation 264 (2015) 310–322

313

3.2. Haar wavelet operational matrix of fractional integration Haar matrix H is obtained by using the collocation points in (2.1), H(l, i) = hl (xc (i)). Similarly, we can obtain the fractional order integration matrix P of Haar function by substituting the collocation points in Eqs. (2.4) and (2.5), P(l, i) = pα , l (xc (i)), as

⎡

P2M×2M

pα ,1 (xc (1)) pα ,1 (xc (2)) · · · ⎢ pα ,2 (xc (1)) pα ,2 (xc (2)) · · · ⎢ =⎢ .. .. .. ⎣ . . . pα ,2M (xc (1)) pα ,2M (xc (2)) · · ·

⎤ pα ,1 (xc (2M)) pα ,2 (xc (2M)) ⎥ ⎥ ⎥. .. ⎦ . pα ,2M (xc (2M))

In particular, we fix J = 2, α = 0.75, we get 2M = 8 and the Haar wavelet operational matrix of fractional integration is

P8×8

⎡ 0.1360 ⎢0.1360 ⎢ ⎢ ⎢0.1360 ⎢ ⎢ 0 ⎢ =⎢ ⎢0.1360 ⎢ ⎢ 0 ⎢ ⎢ ⎣ 0 0

⎤ 0.3100 0.4548 0.5853 0.7067 0.8215 0.9312 1.0367 0.3100 0.4548 0.5853 0.4347 0.2014 0.0216 −0.1340⎥ ⎥ ⎥ 0.3100 0.1828 −0.0347 −0.0668 −0.0391 −0.0275 −0.0210⎥ ⎥ 0 0 0 0.1360 0.3100 0.1828 −0.0347⎥ ⎥ ⎥. 0.0380 −0.0293 −0.0142 −0.0091 −0.0066 −0.0051 −0.0042⎥ ⎥ 0 0.1360 0.0380 −0.0293 −0.0142 −0.0091 −0.0066⎥ ⎥ ⎥ 0 0 0 0.1360 0.0380 −0.0293 −0.0142⎦ 0 0 0 0 0 0.1360 0.0380

3.3. Haar wavelet operational matrix of fractional integration for boundary value problems We drive another operational matrix of fractional integration to solve the fractional boundary value problems. Let η > 0 and g : [0, η] → R be a continuous function and assume that Haar functions have [0, η] as compact support, then

g(x)I0α h1 (η) = g(x)

 η 0

(η − s)α−1 ds,

(3.7)

vα ,η,1 = g(x)Cα ,1 , and

g(x)I0α hl (η) = g(x)



b(l)

a(l)

(η − s)α−1 ds −

 c(l) b(l)



(η − s)α−1 ds ,

(3.8)

vα ,η,l = g(x)Cα ,l , ηα where vα ,η,1 = g(x)I0α h1 (η), vα ,η,l = g(x)I0α hl (η), Cα ,1 =  (α and Cα ,l =  (α1+1) [(η − a(l))α − 2(η − b(l))α + (η − c(l))α ]. Also l = +1) j j 2 + k + 1, j = 0, 1, 2, . . . , J and k = 0, 1, 2, . . . , 2 − 1. Let xc (i) = η i−0.5 2M , i = 1, 2, . . . , 2M and define a matrix V by using the collocation points, xc , in (3.7) and (3.8), we get

⎡

α ,η

V2M×2M

g(xc (2))I0α h1 (η) · · · g(xc (1))I0α h1 (η) ⎢ g(x (1))Iα h (η) g(xc (2))I0α h2 (η) · · · c ⎢ 0 2 ⎢ =⎢ .. .. .. ⎢ . . . ⎣ g(xc (1))I0α h2M (η) g(xc (2))I0α h2M (η) · · ·

⎤ g(xc (2M))I0α h1 (η) g(xc (2M))I0α h2 (η) ⎥ ⎥ ⎥ ⎥. .. ⎥ . ⎦ α g(xc (2M))I0 h2M (η)

In particular, for η = 1, g(x) = x2 , α = 1.25, J = 2, we get

1.25,1 V8×8

⎡ ⎤ 0.0034 0.0310 0.0862 0.1689 0.2793 0.4172 0.5827 0.7757 ⎢0.0005 0.0049 0.0137 0.0269 0.0444 0.0664 0.0927 0.1234⎥ ⎢ ⎥ ⎢ ⎥ ⎢0.0001 0.0008 0.0021 0.0041 0.0069 0.0102 0.0143 0.0190⎥ ⎢ ⎥ ⎢0.0002 0.0021 0.0058 0.0113 0.0187 0.0279 0.0390 0.0519⎥ ⎢ ⎥ =⎢ ⎥. ⎢0.0000 0.0002 0.0005 0.0009 0.0015 0.0023 0.0032 0.0042⎥ ⎢ ⎥ ⎢0.0000 0.0002 0.0006 0.0012 0.0019 0.0029 0.0041 0.0054⎥ ⎢ ⎥ ⎢ ⎥ ⎣0.0000 0.0003 0.0009 0.0017 0.0029 0.0043 0.0060 0.0080⎦ 0.0001 0.0009 0.0024 0.0048 0.0079 0.0117 0.0164 0.0218

The Haar matrices H, P and V are constructed to solve the fractional order initial and boundary value problems.

314

U. Saeed, M. ur Rehman / Applied Mathematics and Computation 264 (2015) 310–322

4. Procedure of implementation In this section, we describe the procedure of implementing the method for fractional nonlinear partial differential equation. The first step is to convert the fractional nonlinear partial differential equation into discrete fractional partial differential equation by Picard technique. The second step is to solve the obtained discretized fractional partial differential equation by the Haar wavelet operational matrix method. Consider the following fractional nonlinear partial differential equation,

∂αu ∂γ u ∂βu − a ( x ) + b ( x ) u + d(x)up = f (t, x), ∂ tα ∂ xγ ∂ xβ 0 < α ≤ 2, 1 < β ≤ 2, 0 < γ ≤ 1, p > 1,

(4.1)

with the initial and boundary conditions,

∂u | = g2 (x), 0 < x < 1, ∂ t t=0 u(0, t) = Y1 (t), u(1, t) = Y2 (t), t ≥ 0.

u(x, 0) = g1 (x),

Applying the Picard iteration [15] to Eq. (4.1), we get





∂ α ur+1 ∂ γ ur ∂ β ur+1 − a(x) = f t, x, ur , , α β ∂t ∂ xγ ∂x 0 < α ≤ 2, 1 < β ≤ 2, 0 < γ ≤ 1, r ≥ 0,

(4.2)

with the initial and boundary conditions,

∂ ur+1 | = g2 (x), 0 < x < 1, ∂ t t=0

ur+1 (x, 0) = g1 (x), ur+1 (0, t) = Y1 (t), γ f (t, x, ur , ∂∂ xuγ r

ur+1 (1, t) = Y2 (t),

γ f (t, x) − b(x)ur ∂∂ xuγ r

t ≥ 0,

where )= − d(x)ur Applying the Haar wavelet method to Eq. (4.2), we approximate the higher order term by the Haar wavelet series as p

2M 2M  ∂ β ur+1  r+1 = cl,i hl (x)hi (t) = HT (x)C r+1 H(t). β ∂x l=1 i=1

(4.3)

β

Apply the fractional integral operator Ix on Eq. (4.3)

ur+1 (x, t) =

2M 2M  

β

r+1 cl,i Ix hl (x)hi (t) + p(t)x + q(t),

(4.4)

l=1 i=1

where p(t) and q(t) are functions of t. Using the boundary conditions and Eqs. (2.2) and (2.3), we get

q(t) = Y1 (t),

p(t) = −

2M 2M  

 r+1 cl,i

l=1 i=1

1

 (β)



1



(1 − s)β −1 h (s)ds l

0

hi (t) + Y2 (t) − Y1 (t),

or

p(t) = −

2M 2M  

  β r+1 cl,i I1 hl (x) hi (t) + Y2 (t) − Y1 (t).

l=1 i=1

Eq. (4.4) can be written as

ur+1 (x, t) =

2M  2M  l=1 i=1

r+1 cl,i pβ ,l (x)hi (t) − x

 2M 2M 

   β r+1 I1 hl (x) hi (t) + x(Y2 (t) − Y1 (t)) + Y1 (t). cl,i

γ Take the Caputo fractional derivative ∂∂xγ , of order γ , to Eq. (4.5)





2M 2M 2M  2M     x1−γ x1−γ ∂ γ ur+1  β r+1 r+1 I h = c p ( x ) h ( t ) − c h ( x ) ( t ) + (Y (t) − Y1 (t)). i i l β − γ ,l 1 l,i l,i ∂ xγ  (2 − γ ) l=1 i=1  (2 − γ ) 2 l=1 i=1

For simplicity, let

(4.5)

l=1 i=1

(4.6)

U. Saeed, M. ur Rehman / Applied Mathematics and Computation 264 (2015) 310–322

  ∂ γ ur S(x, t) = f t, x, ur , , γ ∂x =

2M 2M  

315

(4.7)

ml,i hl (x)hi (t),

l=1 i=1

where ml, i = hl (x), S(x, t), hi (t) . Use Eqs. (4.3) and (4.7) in Eq. (4.2), to obtain 2M 2M 2M  2M    ∂ α ur+1 r+1 = a ( x ) c h ( x ) h ( t ) + ml,i hl (x)hi (t). i l l,i ∂ tα l=1 i=1 l=1 i=1

(4.8)

Apply fractional integral operator Itα to (4.8) and use the initial conditions to obtain,

ur+1 (x, t) = a(x)

2M 2M  

r+1 cl,i hl (x)pα ,i (t) +

l=1 i=1

2M 2M  

ml,i hl (x)pα ,i (t) + g2 (x)t + g1 (x).

(4.9)

l=1 i=1

Let K(x, t) = −g2 (x)t − g1 (x) + x(Y2 (t) − Y1 (t)) + Y1 (t). From Eqs. (4.5) and (4.9), 2M 2M  

r+1 cl,i pβ ,l

(x)hi (t) − x

l=1 i=1

 2M 2M 



r+1 cl,i Cβ ,l hi

(t) + K (x, t)

l=1 i=1

= a(x)

2M 2M  

r+1 cl,i hl (x)pα ,i (t) +

l=1 i=1

2M 2M  

ml,i hl (x)pα ,i (t).

(4.10)

l=1 i=1

β

where Cβ ,l = I1 hl (x) and is given in Section 3. In discrete form, Eq. (4.10) can be written as 2M 2M  

r+1 cl,i pβ ,l

(xc (j))hi (tc (n)) − xc (j)

l=1 i=1

 2M 2M 



r+1 cl,i Cβ ,l hi

(tc (n)) + K (xc (j), tc (n))

(4.11)

l=1 i=1

= a(xc (j))

2M 2M  

r+1 cl,i hl (xc (j))pα ,i (tc (n)) +

l=1 i=1

2M 2M  

ml,i hl (xc (j))pα ,i (tc (n)),

l=1 i=1

2j−1 where xc (j) = 2j−1 4M , tc (n) = 4M , j, n = 1, 2, . . . , 2M, are collocation points. In matrix form, Eq. (4.11) can be written as

β

Px C r+1 H − Vβ ,1,f (x)C r+1 H − AHT C r+1 Pαt − HT M Pαt + K = 0,

(4.12)

β where H is the 2M × 2M Haar matrix, Vβ ,1,f (x) = f (x)I1 HT is the 2M × 2M fractional integration matrix for boundary value problem β β T and P = I H and Pα = Iα H are the 2M × 2M matrices of fractional integration of the Haar functions. These matrices are derived x

x

t

t

in Section 3. Also K = K(xc (j), tc (j)), j = 1, 2, . . . , 2M is the 2M × 2M matrix determined at the collocation points, M is 2M × 2M coefficient matrix determined by the inner products ml, i = hl (x), S(x, t), hi (t) , and f(x) = x. Let Q  (AHT )−1 is the 2M × 2M matrix, where A is the diagonal matrix and is given by

⎡ 0 ··· a(x1 ) ⎢ 0 ( x ) · ·· a 2 ⎢ A=⎢ . . ⎢ . .. .. ⎣ . . 0 0 ···

⎤ 0 0 ⎥ ⎥ .. ⎥ ⎥ . ⎦ a(x2M )

(4.12) can be written as

  β Q Px − Vβ ,1,f (x) C r+1 − C r+1 Pαt H−1 + Q (K − HT M Pαt )H−1 = 0,

(4.13)

which is the Sylvester equation. Solve Eq. (4.13) for Cr + 1 , which is 2M × 2M coefficient matrix, and substituting Cr + 1 in Eq. (4.5) or (4.9), we get solution ur + 1 (x, t) at the collocation points. In particular, given an initial approximation u0 (x, t), we get a linear fractional partial differential equation in u1 (x, t) by substituting r = 0 in Eq. (4.2), which is solved by above procedure to get u1 (x, t) at the collocation points. Similarly for r = 1, we obtain u2 (x, t) and so on. 5. Convergence analysis Let us assume that ur + 1 (x, t) be a continuously differentiable function on [0, 1) × [0, 1), and | 1) × [0, 1), K > 0.

∂ ur+1 (x,t) | ≤ K, for all (x, t)  [0, ∂x

316

U. Saeed, M. ur Rehman / Applied Mathematics and Computation 264 (2015) 310–322 Haar solution at 9th iteration and J=4

Haar solution at α=0.5, β=1.5 Haar solution at α=0.75, β=1.75 Haar solution at α=0.95, β=1.95 Exact solution at α=1, β=2

3.5 3

u(x,t)

2.5 2 1.5 1 0.5 0 0 0.5

1

0

0.2

x

0.4

0.6

0.8

1

t

Fig. 1. Solutions by the Haar wavelet Picard technique at 9th iteration and level of resolution J = 4, for different values of α and β .

Suppose uM r+1 (x, t) be the Haar wavelet approximation for the function u(x, t) at the (r+1)th iteration, that is,

ur+1 (x, t) ≈ uM r+1 (x, t ) =

2M 2M  

r+1 clp hl (x)hp (t).

(5.1)

l=1 p=1

L2 error norm for the Haar wavelet approximation of u(x, t) at the (r+1)th iteration is given as 2 ur+1 (x, t) − uM r+1 (x, t ) =



1 0



1 0

2 (ur+1 (x, t) − uM r+1 (x, t )) dxdt.

Convergence analysis of Haar wavelet method for linear partial differential equations has been discussed in [10]. By following the similar procedure, we get the convergence of developed method as 2 ur+1 (x, t) − uM r+1 (x, t ) ≤

K2 1 , 3 (2M)2

(5.2)

where M = 2J and Eq. (5.2) implies that error between the exact and approximate solution at the (r+1)th iteration is inversely proportional to the maximal level of resolution, J. This implies that uM r+1 (x, t) converges to ur + 1 (x, t) as J → . Since ur + 1 (x, t) is obtained at (r+1)th iteration of Picard technique then according to the convergence analysis of Picard technique [15] which states that ur + 1 (x, t) converges to u(x, t) as r approaches to infinity, if there is convergence at all. This suggests that solution by Haar wavelet Picard technique, uM r+1 (x, t), converges to u(x, t) as J and r approaches to infinity. 6. Application In this section we utilize the proposed method for the solution of fractional nonlinear parabolic equation and Burgers’ equation. The obtained solutions are compared with the exact solution to show the applicability of Haar wavelet Picard method. Problem 1: Consider the fractional nonlinear parabolic problem

∂αu ∂βu − − u2 = et sin(π x), 0 ≤ x ≤ 1, t ≥ 0, 0 < α ≤ 1, 1 < β ≤ 2, ∂ t α ∂ xβ

(6.1)

subject to the initial and boundary conditions

u(x, 0) = sin(π x), u(0, t) = 0, u(1, t) = 0. The exact solution of Eq. (6.1), when α = 1 and β = 2, is given as

u(x, t) = et sin(π x).

(6.2)

Consider u0 (x, t) = sin (π x) as an initial approximation and apply the Haar wavelet Picard technique, as described in Section 4, to Eq. (6.1), we get the approximate solution of the fractional nonlinear parabolic problem. Fig. 1 shows that solutions by Haar wavelet Picard technique of Eq. (6.1) at different values of α , β converges to the exact solution at α = 1, β = 2, when α and β approaches to 1 and 2 respectively. According to Figs. 2 and 3, we get more accurate results while increasing iteration r or level of resolution J or both, as given in convergence analysis.

U. Saeed, M. ur Rehman / Applied Mathematics and Computation 264 (2015) 310–322

st

Haar solution at 1 iteration

Haar solution at 2

0 1

2 0 1

1

0.5 0 0

0 1

0 0

0 0

|u(x,t)−u2(x,t)| 0.4

0.2

t

0 1 0.5 0

rd

Haar solution at 4 iteration

0 1

0 1

1

0.5

2

t

0 0

0 1

0 0

0 0

0.4

0.2

t

0 1

1

0.8

0.6

0.01

0.5 0

th

Haar solution at 6 iteration

0 1

0 1

1

0.5 t

0 0

0 1

0 0

|u(x,t)−u6(x,t)|

0 1 0.5 0

0

0.4

0.2

t

1

0.8

0.6

4 2 0 1 0.5 0 t

x

th

Haar solution at 7 iteration

0.2

0.4

0.6

x

4 u(x,t)

7

0

Exact solution

4 u (x,t)

x

Absolute Error

−3

0.01

0.5 0 0

t

x

x 10

0.005

2 0 1

2 0 1

1

0.5 0 0

1

0.5

0.5 t

0.5 0 0

t

x

x

Absolute Error

−3

x 10 2

7

|u(x,t)−u (x,t)|

|u(x,t)−u5(x,t)|

Absolute Error

0 0

1

0.5

0.5 t

x

2 0 1

1

0.5

0.5 t

x

2

1

0.5

0.5

4 u(x,t)

2

1 0 1 0.5 0 t

0

1

Exact solution

4 u6(x,t)

u(x,t)

4

0.8

0.6

x

th

Exact solution

2

0.4

0.2

0

t

x

Haar solution at 5 iteration

u5(x,t)

0.02

|u(x,t)−u4(x,t)|

|u(x,t)−u3(x,t)|

0.5

4

x

Absolute Error

0 1 0

0.5 0 0

t

x

Absolute Error 0.05

0

1

0.5

0.5 t

x

2 0 1

1

0.5

0.5 t

x

4

2

1

0.5

0.5

1

Exact solution

4 u4(x,t)

u(x,t)

4

0.8

0.6

x

th

Exact solution

2

0.4

0.2

0

t

x

Haar solution at 3 iteration

u3(x,t)

0.1

1

0.8

0.6

0.2

u(x,t)

|u(x,t)−u1(x,t)|

0.5

4

x

Absolute Error

0 1 0

0.5 0 0

t

x

Absolute Error 0.5

0

1

0.5

0.5 t

x

2 0 1

1

0.5

0.5 t

x

2

1

0.5

0.5 t

Exact solution 4 u(x,t)

2

iteration

4 u2(x,t)

4 u(x,t)

4 u1(x,t)

nd

Exact solution

317

0.2

0.4

0.6

0.8

1

x

Fig. 2. Comparison of exact and Haar wavelet solution for α = 1 and β = 2, at J = 4 and different iterations.

0.8

1

318

U. Saeed, M. ur Rehman / Applied Mathematics and Computation 264 (2015) 310–322

Fig. 3. Comparison of exact and Haar wavelet solution for α = 1 and β = 2, at 8th iteration and different J.

Problem 2: Consider the fractional Burgers’ equation

∂αu ∂u ∂βu = +u , 0 ≤ x ≤ 1, t ≥ 0, 0 < α ≤ 1, 1 < β ≤ 2, α ∂t ∂ x ∂ xβ

(6.3)

subject to the initial and boundary conditions

u(x, 0) = 2x, u(0, t) = 0,

u(1, t) =

2 . 1 + 2t

The exact solution of Eq. (6.3), when α = 1 and β = 2, is given as

u(x, t) =

2x . 1 + 2t

(6.4) ∂u

Consider u0 (x, t) = 2x, ∂ x0 = 2, as a initial approximations and apply the Haar wavelet Picard technique to Eq. (6.3), we get the following results: We plot the solutions by Haar wavelet Picard technique by fixing α = 1, β = 2 and J = 7, and compare the obtained solution at different iterations with the exact solution, as shown in Fig. 4. It shows that absolute error reduces while increasing iterations. We also utilize the method for solving fractional Burgers’ Eq. (6.3) at different values of α and β . Fig. 5 shows that obtained solution converges to the exact solution when α and β approaches to 1 and 2 respectively. Problem 3: The Burger–Fisher equation has important applications in various fields of financial mathematics, gas dynamic, traffic flow, applied mathematics and physics applications. This equation shows a prototypical model for describing the

U. Saeed, M. ur Rehman / Applied Mathematics and Computation 264 (2015) 310–322

319

Fig. 4. Comparison of exact solution and solutions by Haar wavelet Picard technique at different iterations by fixing J = 7 and for α = 1 and β = 2.

interaction between the reaction mechanism, convection effect, and diffusion transport [16]. Consider the fractional time derivative generalized Burger–Fisher equation,

∂ α u ∂ 2u ∂u + bu(uγ − 1) = 0, 0 ≤ x ≤ 1, t ≥ 0, 0 < α ≤ 1 − + auγ ∂ t α ∂ x2 ∂x subject to the initial and boundary conditions

(6.5)

320

U. Saeed, M. ur Rehman / Applied Mathematics and Computation 264 (2015) 310–322 Haar solution at 6th iteration and J=3 Haar solution at α=0.4, β=1.4 Haar solution at α=0.6, β=1.6 Haar solution at α=0.8, β=1.8 Exact solution at α=1, β=2

2

u(x,t)

1.5

1

0.5

0 1 1 0.8

0.5

0.6 0.4 t

0

0.2 0

x

Fig. 5. Solution by Haar wavelet Picard technique at different values of α , β and exact solution at α = 1, β = 2. Table 1 Comparison of solution y the Haar wavelet Picard method uHaar , at 6th iteration, J = 5 and α = 1, with exact solution uExact , solution by reduced differential transform method and variational iteration method, by fixing a = 0.001, b = 0.001 and γ = 1. 6th iteration

J=5

x

t

uHaar

uExact

EHaar

ERDTM [16]

EVIM [16]

0.01

0.02 0.04 0.06 0.08 0.02 0.04 0.06 0.08 0.02 0.04 0.06 0.08

0.50000375125000 0.50000875250000 0.50001375375000 0.50001875499999 0.50000000125000 0.50000500250000 0.50001000375000 0.50001500500000 0.49999500125000 0.50000000250000 0.50000500375000 0.50001000500000

0.50000375125000 0.50000875250000 0.50001375375000 0.50001875499999 0.50000000125000 0.50000500250000 0.50001000375000 0.50001500500000 0.49999500125000 0.50000000250000 0.50000500375000 0.50001000500000

2.22045e − 016 1.11022e − 016 3.33067e − 016 2.22045e − 016 1.11022e − 016 0 1.11022e − 016 3.33067e − 016 2.77556e − 016 1.11022e − 016 0 2.22045e − 016

0.49994e − 05 0.49994e − 05 1.49994e − 05 1.99994e − 05 0.49975e − 05 0.99975e − 05 1.49975e − 05 1.99975e − 05 0.49950e − 05 0.99950e − 05 1.49950e − 05 1.99950e − 05

2.50311e − 03 2.50811e − 03 2.51312e − 03 2.51812e − 03 9.99620e − 03 1.00012e − 02 1.00062e − 02 1.00112e − 02 1.99794e − 02 1.99844e − 02 1.99894e − 02 1.99944e − 02

0.04

0.08

u(x, 0) =

u(0, t) =

u(1, t) =



 γ1  1 1 aγ x , − tanh 2 2 2(1 + γ )



   2   γ1 1 aγ 1 a + b(1 + γ )2 − tanh , − t 2 2 2(1 + γ ) a(1 + γ )



   2   γ1 1 aγ 1 a + b(1 + γ )2 − tanh , 1− t 2 2 2(1 + γ ) a(1 + γ )

where a, b andγ are constants. Generalized Burger equation is obtained when b = 0. Apply the proposed scheme, as described in Section (4), to Eq. (6.5). The obtained results are shown in Tables 1, 2 and 3. We fix the order of the differential Eq. (6.5), α = 1, level of resolution, J = 5 and consider a = 0.001, b = 0.001, γ = 1. The comparison of the results by proposed scheme at sixth iteration with the exact solution, reduced differential transform method and variational iteration method are shown in Table 1. Here uHaar represent the solution by Haar wavelet Picard technique, uExact represent the exact solution. EHaar , ERDTM and EVIM represent the absolute error by Haar wavelet Picard technique, reduced differential transform method and variational iteration method respectively. We also solved the generalized Burger–Fisher Eq. (6.5) for a = 0.001, b = 0.001 and γ = 2, and compare the solution by Haar wavelet Picard technique with solution obtained by reduced differential transform method and exact solution as shown in Table 2. Numerical solution of the generalized Burger–Fisher equation for α = 1, J = 6, γ = 1 and at different values of a, b are shown in Table 3. We compared the solution by Haar wavelet Picard technique with homotopy perturbation method. Where EHaar and EHPM denotes the absolute error by Haar wavelet Picard technique and homotopy perturbation method respectively.

U. Saeed, M. ur Rehman / Applied Mathematics and Computation 264 (2015) 310–322

321

Table 2 Comparison of solution by the Haar wavelet Picard technique uHaar , at 6th iteration, level of resolution J = 5 and α = 1, with exact solution uExact and solution by reduced differential transform method, by fixing a = 0.001, b = 0.001 and γ = 2. 6th iteration

J=5

t

x

uHaar

uExact

EHaar

ERDTM [16]

0.01

0.02 0.04 0.06 0.08 0.02 0.04 0.06 0.08 0.02 0.04 0.06 0.08

0.70710796008575 0.70710560306262 0.70710324603111 0.70710088899223 0.70711856777058 0.70711621078328 0.70711385378751 0.70711149678430 0.70713271109673 0.70713035415672 0.70712799720816 0.70712564025213

0.70710796008970 0.70710560306710 0.70710324603664 0.70710088899833 0.70711856777268 0.70711621078544 0.70711385379034 0.70711149678738 0.70713271110241 0.70713035416232 0.70712799721437 0.70712564025857

3.9578e-012 4.4785e-012 5.5282e-012 6.0935e-012 2.0928e-012 2.1563e-012 2.8214e-012 3.0833e-012 5.6780e-012 5.6025e-012 6.2091e-012 6.4343e-012

4.7133e-06 9.4271e-06 1 .4142e-05 1 .8855e-05 4.7117e-06 9.4260e-06 1 .4140e-05 1 .8854e-05 4.7104e-06 9.4241e-06 1 .4138e-05 1.8852e-05

0.04

0.08

Table 3 Comparison of solution by the Haar wavelet Picard technique, at 5th iteration, level of resolution J = 6 and α = 1, with exact solution and solution by homotopy perturbation method at γ = 1 and different values of a and b. J=6

5th iteration

J=6

J=6

x

t

EHPM [17] a = b = 0.01

EHaar a = b = 0.01

EHPM [17] a = b = 0.1

EHaar a = b = 0.1

EHPM [17] a = b = 0.5

EHaar a = b = 0.5

0.1

0.2 0.4 0.6 0.8 0.2 0.4 0.6 0.8 0.2 0.4 0.6 0.8

6.28000e − 11 5.08000e − 11 1.63800e − 10 3.02500e − 10 3.50000e − 10 6.02000e − 10 1.65600e − 9 2.71000e − 9 6.69900e − 9 2.69600e − 9 1.31100e − 9 5.32000e − 9

9.49241e − 014 3.78475e − 013 5.83755e − 013 4.37983e − 013 2.67564e − 014 2.71116e − 013 4.80727e − 013 4.04565e − 013 6.55587e − 014 4.70735e − 014 2.31148e − 013 2.23821e − 013

4.32620e − 8 1.08832e − 7 1.74575e − 7 2.40128e − 7 3.85160e − 7 6.65330e − 7 1.71580e − 6 2.76581e − 6 7.28031e − 6 3.08014e − 6 1.12091e − 6 5.32152e − 6

1.02372e − 010 4.05313e − 010 6.21942e − 010 4.62845e − 010 3.59754e − 011 3.05358e − 010 5.30327e − 010 4.45039e − 010 6.50776e − 011 6.13558e − 011 2.59859e − 010 2.47570e − 010

6.17680e − 6 1.60295e − 5 2.58023e − 5 3.54471e − 5 7.87747e − 5 7.89518e − 5 2.36284e − 4 3.92443e − 4 1.24466e − 3 6.22452e − 4 2.80910e − 6 6.28040e − 4

2.01866e − 008 7.60004e − 008 1.11933e − 007 8.52507e − 008 1.60003e − 008 8.94412e − 008 1.51149e − 007 1.57278e − 007 8.96199e − 009 2.03433e − 008 6.30765e − 008 6.62272e − 008

0.4

0.8

7. Conclusion We have constructed the Haar wavelets matrix, Haar wavelet operational matrix of fractional order integration and Haar wavelets operational matrix of fractional order integration for boundary value problems. We have successfully utilized the combination of Haar wavelet operational matrices method and Picard technique for the solutions of fractional nonlinear partial differential equations. The Haar wavelet Picard method gives excellent results when applied to different fractional nonlinear partial differential equations. Absolute error by the Haar wavelet Picard method reduces while increasing iterations or level of resolution or both. The solution of the fractional nonlinear partial differential equation converges to the solution of integer order partial differential equation, when α , β approach to integer values, as shown in Figs. 1 and 5. The comparison in Tables 1–3 show that our results are more accurate as compared to reduced differential transform method [16], variational iteration method [16] and homotopy perturbation method [17] for the Burger–Fisher equation. Different types of non-linearities can easily be handled by the Haar wavelet Picard technique. Acknowledgment We are grateful to the anonymous reviewers for their valuable comments which led to the improvement of the manuscript. References [1] S.S. Ray, R.K. Bera, An approximate solution of a nonlinear fractional differential equation by adomian decomposition method, Appl. Math. Comput. 167 (2005) 561–571. [2] S. Momani, Z. Odibat, Homotopy perturbation method for nonlinear partial differential equations of fractional order, Phys. Lett. A 365 (2007) 345– 350. [3] R. Y. Molliq, M.S.M. Noorani, I. Hashim, Variational iteration method for fractional heat- and wave–like equations, Nonlinear Anal. Real World Appl. 10 (2009) 1854–1869.

322

U. Saeed, M. ur Rehman / Applied Mathematics and Computation 264 (2015) 310–322

[4] V.S. Erturk, S. Momani, Z. Odibat, Application of generalized differential transform method to multi-order fractional differential equations, Commun. Nonlinear Sci. Numer. Simul. 13 (2008) 1642–1654. [5] U. Saeed, M. Rehman, Haar wavelet–quasilinearization technique for fractional nonlinear differential equations, Appl. Math. Comput. 220 (2013) 630– 648. [6] I. Celik, Haar wavelet method for solving generalized burgers–huxley equation, Arab J. Math. Sci. 18 (2012) 25–37. [7] C. Cattani, Haar wavelet spline, J. Interdiscipl. Math. 4 (2001) 35–47. [8] C.F. Chen, C.H. Hsiao, Haar wavelet method for solving lumped and distributed-parameter systems, IEE Proc. Control Theory Applic. 144 (1997). [9] G. Hariharan, K. Kannan, An overview of haar wavelet method for solving differential and integral equations, World Appl. Sci. J. 23 (2013) 1–14. [10] L. Wang, Y. Ma, Z. Meng, Haar wavelet method for solving fractional partial differential equations numerically, Appl. Math. Comput. 227 (2014) 66–76. [11] Ü. Lepik, Solving PDEs with the aid of two-dimensional Haar wavelets, Comp. Math. Appl. 61 (2011) 1873–1879. [12] I. Celik, Haar wavelet approximation for magnetohydrodynamic flow equations, Appl. Math. Model. 37 (2013) 3894–3902. [13] G. Hariharan, K. Kannan, K.R. Sharma, Haar wavelet method for solving Fisher’s equation, Appl. Math. Comput. 211 (2009) 284–292. [14] G. Hariharan, K. Kannan, Haar wavelet method for solving some nonlinear parabolic equations, J. Math. Chem. 48 (2010) 1044–1061. [15] R.E. Bellman, R.E. Kalaba, Quasilinearization and Nonlinear Boundary-value Problems, American Elsevier Publishing Company, 1965. [16] D. Kocacoban, A.B. Koc, A. Kurnaz, Y. Keskin, A better approximation to the solution of Burger–Fisher equation, Proc. World Congress Eng. I (2011). [17] M.M. Rashidi, D.D. Ganji, S. Dinarvand, Explicit analytical solutions of the generalized Burger and Burger–Fisher equations by homotopy perturbation method, Numer. Meth. Partial Differ. Equ. 25 (2009) 409–417.