Numerical Solutions of Fractional System of Partial Differential Equations By Haar Wavelets

F.Bulut,Ö.Oruçand A.Esen

Numerical Solutions of Fractional System of Partial Differential Equations By Haar Wavelets

F.Bulut1,2,Ö.Oruç3and A.Esen3

In this paper,time fractional one dimensional coupled KdV and coupled modified KdV equations are solved numerically by Haar wavelet method.Proposed method is new in the sense that it doesn’t use fractional order Haar operational matrices.In the proposed method L1 discretization formula is used for time discretization where fractional derivatives are Caputo derivative and spatial discretization is made by Haar wavelets.L2and L∞error norms for various initial and boundary conditions are used for testing accuracy of the proposed method when exact solutions are known.Numerical results which produced by the proposed method for the problems under consideration confirm the feasibility of Haar wavelet method combined with L1 discretization formula.

Haar wavelet method,Fractional coupled KdV equation,Fractional coupled MKdV equation,Linearization,Numerical solution.

1 Introduction

In this paper we will consider time fractional coupled KdV(FCKdV)equation

where a and b are constants and time fractional coupled modified KdV(FCMKdV)equation

where

is the fractional derivative in the Caputo’s sense[Kilbas,Srivastava,and Trujillo(2006);Podlubny(1999)],here n is an integer.When α=β=1,Eq.(1)corresponds to classical coupled KdV equation which was first introduced by Hirota and Satsuma(1981)in 1981 and describes interaction of two long waves with different dispersion relations.When α =β =1,Eq.(2)corresponds to integer order coupled modified KdV equation[Fan(2000,2001,2002)].

In recent years,fractional differential equations have been broadly used in different branches of physics and engineering[West,Bolognab,and Grigolini(2003)].For example;nonlinear oscillation of earthquake[He(1999)],diffusion in a certain type of porous medium[Podlubny(1999)],phenomenons occurred in electromagnetics,acoustics electrochemistry and material science[Podlubny(1999);Chechkin,Gorenflo,Sokolov,and Gonchar(2003);Miller and Ross(1993);Samko,Kilbas,and Marichev(1993)],can be well modeled by the aid of fractional derivative.Many researchers focus their attention to fractional differential equation due to suitable applications of fractional differential equations in different fields of engineering and science[Agrawal(2002)].Even though there are some techniques for getting analytical solutions of the fractional differential equations[Kilbas,S-rivastava,and Trujillo(2006);Agrawal(2002)],usually their exact solutions are not known.Therefore various numerical techniques devised by many researchers for obtaining approximate solutions of fractional differential equations.Finite difference method[Sun,Chen,and Li(2012);Meerschaert,Scheffler,and Tadjeran(2006);Yuste(2006);Yuste and Acedo(2005);Quintana-Murillo and Yuste(2011);Sweilam,Khader,and Mahdy(2012);Celik and Duman(2012)],finite element method[Sun,Chen,and Sze(3 Sep 2011);Esen,Ucar,Yagmurlu,and Tasbozan(2013)],second kind Chebyshev wavelets method with differential operator matrix and the product operation matrix(Chen,Sun,Li,and Fu,2013),Bernstein polynomials operational matrixes method[Chen,Liu,Li,and Sun(2014)],homotopy analysis method[Hashim and Abdulaziz(2009)],generalized differential transform method[Momani and Odibat(2007);Odibat and Momani(2008);L.J.Cun and H.G.Lin(2010)],adomian decomposition method[Hosseini(2006);EI-Kalla(2008);YongandHong-Li(2008)]and variational iteration method[Odibat(2010)]are some of these techniques.

Recently,for getting numerical solutions of differential equations,a lot of techniques related to Haar wavelet based methods are developed.We can put in order some of studies specially devoted to fractional differential equations as follows.Fractional Volterra and Fredholm integral equations are solved by Lepik(2009).Wu(2009),solved some fractional partial differential equations by fractional order Haar wavelet operational matrix.Li and Zhao(2010),solved Bagley-Torvik,Ricatti and composite fractional oscillation equations with Haar wavelet operational matrix of the fractional order.Mujeeb[Rehman and Khan(2013)]et al.,employed Haar wavelets for obtaining numerical solutions of boundary value problems for linear fractional partial differential equations.Ray and Patra(2013)solved fractional order nonlinear oscillatory Van der Pol system with Haar wavelet operational matrix.Nonlinear fractional ordinary differential equations are solved by Haar wavelet and a quasilinearization technique by Saeed and Rehman(2013).Wang,Ma,and Meng(2014)solved some fractional partial differential equations by Haar wavelet operational matrix.Neutron point kinetics equation is solved with Haar wavelet operational method by Patra and Ray(2014b,a).Variable coefficient fractional differential equations are solved by M.Yi and Huang(2014)with Haar wavelet operational matrix of fractional order.Most of studies aforementioned use Haar wavelet operational matrix of fractional order for obtaining numerical solutions of fractional differential equations.As an alternative approach,in this study we propose a relatively new and simple technique for obtaining numerical solutions of fractional differential equations which consists of Haar wavelets and L1 discretization formula.

In this paper,to get numerical solutions of systems(1)and(2),we have applied Haar wavelet method with L1 discretization formula.The paper organized as follows;In Section 2,Haar wavelets are introduced.Time and spatial discretizations are described in Section 3.Numerical results are given in Section 4 and finally the paper concluded in Section 5.

2 Haar wavelets

Haar wavelets were first introduced by Hungarian mathematician Alfred Haar in 1910 which are the simplest of possible wavelets with compact support.They are piecewise constant functions which form an orthonormal system on the interval[0,1)in the space of square integrable functions.Because of their discontinuity they can not be directly used in solution process of differential equations.To alleviate this situation,Chen and Hsiao(1997)used the integral method,in which the highest derivative appeared in the differential equation is expanded into Haar series.In recent years,Haar wavelets are used widely to obtain numerical solutions of differential equations and are favored by researchers because of their simplicity and desirable computational features.The ith Haar wavelet is defined as follows for x∈[0,1]

where m=2j,j=0,1,...,J and k=0,1,...,m-1 is dilation parameter and translation parameter,respectively.The index of hiin Eq.(3)can be found by relation i=m+k+1.For the lowest values of m=1,k=0,we have i=2 and the greatest value of i will be i=2M=2J+1;where J is the maximum resolution of the wavelet.For i=1 we have Haar scaling function

Any function u(x)∈L2[0,1)can be expanded into Haar series as

where cican be found by

In practice,for approximating a function u(x)∈L2[0,1),finite terms of Haar series are needed,hence one may write

where M=2jand T denotes transpose.

While using Haar wavelet method for solving any order partial differential equation one needs to following integrals in solution process.

The first three integrals can be calculated from Eq.(3)as follows;

3 Discretization scheme for the equations

3.1 Time discretization for FCKdV

In this section we start with the discretization of the fractional time derivative that appears in Eq.(1)in the Caputo’s sense by L1 formula[Oldham and Spanier(1974)]

where

We also take time averages of the other terms of Eq.(1),as follows

We apply Rubin Graves linearization[Rubin and Graves(1975)]formula un+1(ux)n+un(ux)n+1-(uux)nto nonlinear term(uux)n+1.Similar linearization is made for(vvx)n+1and(uvx)n+1terms.Hence we get

with initial conditions

and boundary conditions

where un+1and vn+1are the solutions of the Eq.(8)at the(n+1)th time step,Δt denotes the step size in time so that tn=n×Δt and

3.2 Time discretization for FCMKdV

Again the Caputo derivative appeared in Eq.(2)is discretized by L1 formula as made in subsection and time averages used for some of the terms of Eq.(2)as follows

In the above equation by using linearizationswe get

with initial conditions

and boundary conditions

where un+1and vn+1are the solutions of the Eq.(10)at the(n+1)th time step,Δt denotes the step size in time so that tn=n×Δt and

3.3 Space discretization by Haar wavelets

In this subsection we describe discretization of spatial derivatives that appeared in Eqs.(1),(2).We start with the highest derivative by Haar wavelets.To do so we assume

Integrating Eq.(12)with respect to x from 0 to x,we get the following equation

In Eq.(13),(uxx)n+1(0)is unknown so to find it,we need to integrate Eq.(13)from 0 to 1.After that,using boundary conditions(9)we get

Substituting(14)into Eq.(13)results in the following equation

Now,if we integrate Eq.(15)from 0 to x we get

In Eqs.(14),(15)and(16),(ux)n+1(0)term is unknown.So to find(ux)n+1(0)term we integrate Eq.(16)from 0 to 1 and use boundary conditions(9).Therefore we have

Now by plugging the calculated value of(ux)n+1(0)into Eq.(16)we obtain

Finally,integrating(17)from 0 to x,we obtain

If we summarize,we have

Similarly,we have

For FCKdV system,if we substitute Eqs.(19),(20)into Eq.(8)and discretize the results at the collocation pointswe obtain following system of equations

where

ciand diare wavelet coefficients.By solving Eq.(21),wavelet coefficients ciand dican be calculated successively.Then by plugging these coefficients into the Eqs.(19)and(20),the numerical solutions can be constructed.

For FCMKdV system,if we substitute Eqs.(19),(20)into Eq.(10)and discretize the results at the collocation pointswe obtain following system of equations

where

ciand diare wavelet coefficients.By solving Eq.(22),wavelet coefficients ciand dican be calculated successively.Again by plugging these coefficients into the Eqs.(19)and(20),the numerical solutions can be constructed.

3.4 Error analysis

To analyze the convergence of the method,we use the asymptotic expansion of Eq.(18)as given in Kumar and Pandit(2014),the resulting equation is as follows

Lemma 1.Suppose that u(x)∈L2(R)withandThen

Lemma2.Let u(x)∈L2(R)be a continuous function defined in(0,1).Then the error norm at J th level satisfies the following inequality

Theorem:Suppose that u(x)is exact and u2M(x)is approximate solution of the Eq.(18),then

Proof.See Kumar and Pandit(2014)

Similar approach is valid for v2M(x).It is clear from above equation that the level of resolution J of the Haar wavelet is inversely proportional to the error bound.Therefore the accuracy of the method increases as we increase the level of resolution J.

4 Numerical Examples

Free software package GNU Octave is used in this study for numerical computations and Matplotlib package is used[Hunter(2007)]for graphical outputs.In order to check accuracy of the proposed methods we considered the error norms L2and L∞defined by

4.1 Example 1

Firstly,we consider following initial conditions for the Eq.(1)

and the boundary conditions

This problem have the following exact solution[Hirota and Satsuma(1981)]for α=β=1.

where

We solve the problem for the various values of 2M in the interval-25≤x≤25 and tabulated the L2,L∞error norms when α = β =1 in Table 1.We can easily see from the table that when we increase the collocation points the error norms decrease and they are sufficiently small for t=0.1.To see the accuracy of the Haar wavelet method we depicted the evolution of numerical solutions of u and v at α = β =1 in Fig.1.also for α =1/2,β =1/3 in Fig.2.The comparison of these figures shows that these results are almost same and this verifies the accuracy of Haar wavelet method.

Table 1:Numerical results for Δt=0.005,λ =0.5,a=-0.125,b=-3 and α=β=1 at t=0.1

4.2 Example 2

Secondly,we consider the Eq.(1)with following initial conditions

Figure 1:Exact solutions for α = β =1,Δt=0.005 and 2M=256 at t=0.1

Figure 2:Numerical solutions for α =1/2,β =1/3,Δt=0.005 and 2M=256 at t=0.1

Boundary conditions are taken from exact solutions,computer simulations for this problem are done with parameters k=1,a=1 and b=-1 in the interval-10≤x≤10.L2and L∞error norms for α =β =1 at different time levels are tabulated in Table 2.We see from the table that as t gets larger so as the error norms.We depicted the exact solutions of the problem for α =1,β =1 in Fig.3.,and to see the evolution of numerical solutions of u and v for fractional derivatives we also depicted them for α =1/2,β =1/3 in Fig.4.Comparison of these figures also proves that the present method gives results with high accuracy.

Table 2:Numerical results for 2M=128,Δt=0.005,k=1,a=1,b=-1 and α=β=1 at different times.

Figure 3:Exact solutions for α = β =1,Δt=0.005 and 2M=256 at t=0.1

4.3 Example 3

Finally,we consider Eq.(2)with initial conditions

The exact solutions of this problem for α =β =1 given by Fan(2000,2001,2002)as follows

Figure 4:Numerical solutions for α =1/2,β =1/3,Δt=0.005 and 2M=256 at t=0.1

Figure 5:Exact solutions for α = β =1,Δt=0.001 and 2M=256 at t=0.05

where

The boundary conditions taken from exact solutions.Numerical simulation of this problem are done for b=1,λ=0.1 and k=1/3 in the interval-15≤x≤15 for growing times,L2and L∞error norms are tabulated when α = β =1 in Table 3.We see from the table the error norms are sufficiently small.We depicted the exact solutions for Δt=0.001,α =1,β =1 at t=0.05 in Fig.5.also we depicted the numerical solutions for α =1/2,β =2/3 in Fig.6.Comparison of these figures proves that high accuracy results can be achieved by using Haar wavelet method.

Table 3:Numerical results for Δt=0.005,k=1/3,λ =0.1,b=1,2M=128 and α=β=1.

Figure 6:Numerical solutions for α =1/2,β =2/3,Δt=0.001 and 2M=256 at t=0.05

5 Conclusion

In this paper,we applied Haar wavelet method integrated with L1 discretization formula in the Caputo’s sense to time fractional nonlinear coupled partial differential equations with various initial and boundary conditions.The main idea of the proposed method is using L1 formula for fractional time derivatives and Haar wavelets for spatial derivatives.This approach gives a system of algebraic equations which can be solved easily with suitable methods.Also the proposed method is novel in the sense that it doesn’t use fractional order Haar operational matrices which is a more complex method when compared with the proposed method.The numerical results are quite good in the cases where the exact solutions are known.The proposed method can easily handle similar fractional system of partial differential equations.

Agrawal,O.(2002):Solution for a fractional diffusion-wave equation defined in a bounded domain.Nonlinear Dynamics,vol.29,pp.145-155.

Celik,C.;Duman,M.(2012):Crank-nicolson method for the fractional diffusion equation with the riesz fractional derivative.Journal of Computational Physics,vol.231,pp.1743-1750.

Chechkin,A.;Gorenflo,R.;Sokolov,I.;Gonchar,V.Y.(2003): Distributed order time fractional diffusion equation. Fractional Calculus and Applied Analysis,vol.6,pp.259-280.

Chen,C.;Hsiao,C.(1997): Haar wavelet method for solving lumped and distributed parameter systems.Control Theory and Applications IEE Proceedings,vol.144,pp.87-94.

Chen,Y.;Liu,L.;Li,X.;Sun,Y.(2014):Numerical solution for the variable order time fractional diffusion equation with bernstein polynomials.Cmes-Computer Modeling in Engineering&Sciences,vol.97(1),pp.81-100.

Chen,Y.;Sun,L.;Li,X.;Fu,X.(2013):Numerical solution of nonlinear fractional integral differential equations by using the second kind chebyshev wavelets.Cmes-Computer Modeling in Engineering&Sciences,vol.90(5),pp.359-378.

EI-Kalla,I.(2008): Convergence of the adomian method applied to a class of nonlinear integral equations.Applied Mathematics and Computation,vol.21,pp.372-376.

Esen,A.;Ucar,Y.;Yagmurlu,N.;Tasbozan,O.(2013):A galerkin finite element method to solve fractional diffusion and fractional diffusion-wave equations.Mathematical Modelling and Analysis,vol.18(2),pp.260-273.

Fan,E.(2000):Extended tanh-function method and its applications to nonlinear equations.Physics Letters A.,vol.277,pp.212-218.

Fan,E.(2001): Soliton solutions for a generalized hirota-satsuma coupled kdv equation and a coupled mkdv equation.Physics Letters A,vol.282,pp.18-22.

Fan,E.(2002):Solving kadomtsev-petviashvili equation via a new decompositionand darboux transformation. Communications in Theoretical Physics.(Beijing,China),vol.37,pp.145-148.

Hashim,I.;Abdulaziz,O.(2009):Homotopy analysis method for fractional ivps.Communications in Nonlinear Science and Numerical Simulation,vol.14(3),pp.674-684.

He,J.H.(1999):Some applications of nonlinear fractional differential equations and their approximations.Bulletin of Science,Technology&Society,vol.15(2),pp.86-90.

Hirota,R.;Satsuma,J.(1981):Soliton solutions of a coupled korteweg-de vries equation.Physics Letters A,vol.85(8-9),pp.407-408.

Hosseini,M.(2006): Adomian decomposition method for solution of nonlinear differential algebraic equations.Applied Mathematics and Computation,vol.181,pp.1737-1744.

Hunter,J.D.(2007): Matplotlib:A 2d graphics environment.Computing In Science&Engineering,vol.9(3),pp.90-95.

Kilbas,A.A.;Srivastava,H.;Trujillo,J.(2006): Theory and Applications of Fractional Differential Equations.Elsevier,Amsterdam.

Kumar,M.;Pandit,S.(2014):A composite numerical scheme for the numerical simulation of coupled burgers’equation.Computer Physics Communications,vol.185(3),pp.809-817.

Lepik,U.(2009): Solving fractional integral equations by the haar wavelet method.Applied Mathematics and Computation,vol.214,pp.468-478.

Li,Y.;Zhao,W.(2010): Haar wavelet operational matrix of fractional order integration and its applications in solving the fractional order differential equations.Applied Mathematics and Computation,vol.216,pp.2276-2285.

L.J.Cun;H.G.Lin(2010):New approximate solution for time-fractional coupled kdv equations by generalised differential transform method.Chinese Physics B,vol.19(11),pp.110203.

Lu,H.;Wang,M.(1999): Exact soliton solutions of some nonlinear physical models.Physics Letters A,vol.255,pp.249-252.

Meerschaert,M.;Scheffler,H.;Tadjeran,C.(2006):Finite difference methods for two-dimensional fractional dispersion equation. Journal of Computational Physics,vol.211,pp.249-261.

Miller,K.;Ross,B.(1993): An Introduction to the Fractional Calculus and Fractional Differential Equations.Wiley,New York.

Momani,S.;Odibat,Z.(2007): Generalized differential transform method for solving a space and time-fractional diffusion-wave equation.Physics Letters A,vol.370,pp.379-387.

M.Yi;Huang,J.(2014): Wavelet operational matrix method for solving fractional differential equations with variable coefficients.Applied Mathematics and Computation,vol.230,pp.383-394.

Odibat,Z.(2010): A study on the convergence of variational iteration method.Mathematical and Computer Modelling,vol.51,pp.1181-1192.

Odibat,Z.;Momani,S.(2008): Generalized differential transform method:application to differential equations of fractional order.Applied Mathematics and Computation,vol.197,pp.467-477.

Oldham,K.B.;Spanier,J.(1974): The Fractional Calculus.Academic,New York.

Patra,A.;Ray,S.(2014): A numerical approach based on haar wavelet operational method to solve neutron point kinetics equation involving imposed reactivity insertions.Annals of Nuclear Energy,vol.68,pp.112-117.

Patra,A.;Ray,S.(2014): Numerical simulation based on haar wavelet operational method to solve neutron point kinetics equation involving sinusoidal and pulse reactivity.Annals of Nuclear Energy,vol.73,pp.408-412.

Podlubny,I.(1999): Fractional Differential Equations.Academic,San Diego.

Quintana-Murillo,J.;Yuste,S.(2011): An explicit difference method for solving fractional diffusion and diffusion-wave equations in the caputo form. Journal of Computational and Nonlinear Dynamics,vol.6:021014,pp.http://dx.doi.org/10.1115/1.4002687.

Ray,S.(2012): On haar wavelet operational matrix of general order and its application for the numerical solution of fractional bagley torvik equation.Applied Mathematics and Computation,vol.218,pp.5239-5248.

Ray,S.S.;Patra,A.(2013):Haar wavelet operational methods for the numerical solutions of fractional order nonlinear oscillatory van der pol system. Applied Mathematics and Computation,vol.220,pp.659-667.

Rehman,M.;Khan,R.A.(2013):Numerical solutions to initial and boundary value problems for linear fractional partial differential equations.Applied Mathematical Modelling,vol.37,pp.5233-5244.

Rubin,S.G.;Graves,R.A.(1975): Cubic spline approximation for problems in fluid mechanics.Washington,DC.

Saeed,U.;Rehman,M.(2013): Haar wavelet-quasilinearization technique forfractional nonlinear differential equations. Applied Mathematics and Computation,vol.220,pp.630-648.

Samko,S.;Kilbas,A.;Marichev,O.(1993): Fractional Integrals and Derivatives:Theory and Applications.Gordon and Breach,Yverdon.

Sun,H.;Chen,W.;Li,C.(2012): Finite difference schemes for variable-order time fractional diffusion equation.International Journal of Bifurcation and Chaos,vol.22(4).

Sun,H.;Chen,W.;Sze,K.(3 Sep 2011): A semi-analytical finite element method for a class of time-fractional diffusion equations. arXiv:1109.0641v1,[math-ph].

Sweilam,N.;Khader,M.;Mahdy,A.(2012):Crank-nicolson finite difference method for solving time-fractional diffusion equation.Journal of Fractional Calculus and Applications,vol.2,pp.1-9.

Wang,L.;Ma,Y.;Meng,Z.(2014):Haar wavelet method for solving fractional partial differential equations numerically.Applied Mathematics and Computation,vol.227,pp.66-76.

West,B.;Bolognab,M.;Grigolini,P.(2003): Physics of Fractal Operators.Springer,New York.

Wu,J.L.(2009): A wavelet operational method for solving fractional partial differential equations numerically.Applied Mathematics and Computation,vol.214,pp.31-40.

Yong,C.;Hong-Li,A.(2008): Numerical solutions of a new type of fractional coupled nonlinear equations.Communications in Theoretical Physics,vol.49,pp.839-844.

Yuste,S.(2006):Weighted average finite difference methods for fractional diffusion equations.Journal of Computational Physics,vol.216,pp.264-274.

Yuste,S.;Acedo,L.(2005):An explicit finite difference method and a new von neumann-type stability analysis for fractional diffusion equations.SIAM Journal of Numerical Analysis,vol.42,pp.1862-1874.

1˙Inonu University,Department of Physics,Malatya,TURKEY.

2Corresponding author:E-mail:fatih.bulut@inonu.edu.tr Tel:+904223773745,Fax:+904223410 037.

3˙Inonu University,Department of Mathematics,Malatya,TURKEY.