IBN AL- HAITHAM J. FO R PURE & APPL. SC I. VO L.24 (1) 2011 Explicit Finite Difference Approximation for the Two- Dimensional Fractional Dispersion Equation I. I. Gorial Department of Mathematics, Ibn Al–Haitham College Education, Baghdad University Received in Feb,7,2010 Accepted in May,13,2010 Abstract In this paper, we introduce and discuss an algorithm for the numerical solution of two- dimensional fractional dispersion equation. The algorithm for the numerical solution of this equation is based on explicit finite difference approximation. Consistency, conditional stability, and convergence of this numerical method are described. Finally, numerical example is presented to show the dispersion behavior according to the order of the fractional derivative and we demonstrate that our explicit finite difference approximation is a computationally efficient method for solving two-dimensional fractional dispersion equation. Key words: Fractional derivative, Two-dimensional probem, Explicit Euler method, fractional dispersion equation, S tability, Convergence Introduction The space fractional dispersion equation is obtained from the classical dispersion equation by replacing the second space derivative by a fractional derivative. Numerical methods associated with integer-order differential equations have been treated extensively in the literature. On other hand, studies of the numerical methods and error estimates of fractional order differential equations are quite limited to date [1, 2, 3]. Many works by researchers from various fields of science and engineering deal with dynamical systems described by fractional partial differential equations, which have been used to represent many natural processes in physics [4], finance [5,6], and hydrology [7,8]. In this paper, we find the numerical solution of the two-dimensional fractional dispersion equation of the form: ),,( ),,( ),( ),,( ),( ),,( tyxq y tyxu yxb x tyxu yxa t tyxu              …………(1) subject to the initial condition u (x,y,0) = f(x,y), for x0  x  xR and y0  y yR …………………(2) and the boundary conditions u (x0,y,t) = 0, for y0  y yR and 0 t  Τ u (x,y0,t) = 0, for x0  x  xR and 0 t  Τ …………(3) u (xR,y,t) = g(y,t), for y0  y yR and 0 t  Τ u (x,yR,t) = k(x,t), for x0  x  xR and 0 t  Τ IBN AL- HAITHAM J. FO R PURE & APPL. SC I. VO L.24 (1) 2011 where a, b and f are known functions of x and y, g is a known function of y and t, k is aknwon function of x and t.  and  are given fractional number. q is a knwon function of x, y and t. We use a variation on the classical explicit Euler method. We prove this method by using a novel shifted version of the usual Grunwaled finite difference an approximation for the non-local fractional derivative operator. Explicit Finite Difference Approximation for Solving the Two-Dimensional Fractional Dispersion Equation In this section, we propose explicit finite difference approximation for solving the initial and boundary value problem two-dimensional fractional dispersion equation (1)-(3). The finite difference method starts by dividing the x-interval [x0, xR] into n subintervals to get the grid points xi = x0 + ix, where   nxxx R 0 and i=0,1,…,n. Also we divide the y- interval [y0, yR] into m subintervals to get the grid points yj = y0 + jy, where   myyy R 0 and j=0,1,…,m. Also, the t-interval [0,T] is divided into M subintervals to get the grid points ts = st, s = 0, …,M, where MTt  . Now, we evaluate eq.(1) at (x i, y j, t s) and we use the explicit finite difference approximation to get ),,( ),,( ),( ),,( ),( ),,(),,( 1 sji sji ji sji ji sjisji tyxq y tyxu yxb x tyxu yxa t tyxutyxu              …..(4) Then use the shifted Grunwald estimate to the  , - the fractional derivative, [9]: )(),,)1(( )( 1),,( 0 , xOtyxkxug xx tyxu M k k          …………….(5) )(),)1(,( )( 1),,( 0 , yOtykyxug yy tyxu M k k          to reduce eq.(4) as in the following form ),,( 1 ),( 1 ),( 1 0 1,,,1 1 0 , , 1 , sji j k s kjikji s jki i k kji s ji s ji tyxqug y yxbug x yxa t uu                          s ji s kji j k k ji i k s jkik ji s ji s ji qug y b ug x a t uu ,1, 1 0 , , 1 0 ,1, ,, 1 ,                 , 1,...,1  ni , Msmj ,...,0,1,...,1  ……………. (6) where ),,(, sji s ji tyxuu  , ),(, jiji yxaa  , ),(, jiji yxbb  , ),,(, sji s ji tyxqq  , ! )1()1( )1(, k k g k k      , k=0,1,2,… and ! )1()1( )1(, k k g k k      , k=0,1,2,… The resulting equation can be explicitly solved for 1 , s jiu to give s ji s ji s kji j k kji i k s jkikji s ji utqug y t bug x t au ,,1, 1 0 ,, 1 0 ,1,, 1 ,                 ...............(7) Also form the initial condition and boundary conditions one can get IBN AL- HAITHAM J. FO R PURE & APPL. SC I. VO L.24 (1) 2011 jiji fu , 0 ,  , i=0,…, n 0,0 s ju , j=0,…, m and s=1,…,M 00, s iu , i=0,…, n and s=1,…,M s j s jR gu , , j=0,…, m and s=1,…,M s i s Ri ku , , i=0,…, n and s=1,…,M where ),,(, sjiji tyxff  , ),( sj s j tygg  and ),( si s i txk  Stability Analysis of the Explicit Finite Difference Approximation Define the following fractional partial difference operators:       1 0 ,1, , ,, i k s jkik jis jix ug x a u  and s kji j k k jis jiy ug y b u 1, 1 0 , , ,,        which is an )( xO  approximation to the  th fractional derivative and )( yO  approximation to the th fractional derivative term. Then eq.(7) may be written in the operator form s ji s jiyx s ji tquttu ,,,, 1 , )1(    ................(8) eq.(8) may be written in form s ji s jiyx s ji tquttu ,,,, 1 , )1)(1(    ..................(9) where Ts mn s m s m s n sss n sss uuuuuuuuuU ],,,,,,,,,,,,[ 1,11,21,12,12,22,11,11,21,1   To solve the problem for each fixed yj to obtain an intermediate solution  jiu , from 1 ,,,, )1(   s ji s jijix utqut  ......................(10) Then solve for each fixed xi  ji s jiy uut ,,, )1(  .......................(11) Now, we must prove each one-dimensional explicit system defined by the linear difference eqs. (10) and (11) is conditionally stable for all 1 < < 2, 1 < < 2. Theorem: The explicit system defined by the linear difference eqs.(10) and (11) with 1< < 2, 1 < < 2 is conditionally stable if max 1 ax t     and max 1 by t     Proof: At each grid point yk, for 1,,1  mk  , the system of equation defined by eq.(10) can be written in the explicit matrix form s kkk s k QtUCU  1 where IBN AL- HAITHAM J. FO R PURE & APPL. SC I. VO L.24 (1) 2011 ,],,,[ 1 ,1 1 ,2 1 ,1 1 Ts kn s k s k s k uuuU      ,],,,[ ,1,2,1 T knkkk uuuU      Ts kn s k s k s k qqqQ ],,,[ ,1,2,1   kC is the matrix of coefficients, and is the sum of a lower triangular matrix and a super diagonal matrix at the grid point y k, where the matrix entries along the ith row are defined from eq.(10). For example, for i = 1 the equation becomes s kkkkkkk s k tqugugugu ,1,20,,1,11,,1,02,,1 1 ,1 )1(     for i = 2 we have s kkkkkkkkk s k tqugugugugu ,2,30,,2,21,,2,12,,2,03,,2 1 ,2 )1(     and for i = n - 1 we get          knknknknknkn s kn ugugugu ,11,,1,22,,1,0,,1 1 ,1 )1(    s knknkn tqug ,1,0,,1     Where the coefficients  x t a jiki    ,, , Therefore the resulting matrix entries jiC , for 1,,1  ni  and 1,,1  nj  by are defined by          1,, 0,, 2,, 1,, , 1 jki ki ki ki ji g g g g C         for for for for 1 1 1     ij ij ij ij According to the Greshgorin theorem [9], the eigenvalues of the matrix C lie in the union of the circles centered at iic , with radius     n il l lii cr 0 , . Here we have   kikiii gc ,1,,, 11  and          n il l i il l kilikilii gcr 0 1 0 ,1,,,   and therefore 1,  iii rc . We also have  kikikiiii rc ,,,, 211         x t a ki,21          x t amax21 Therefore, for the spectral radius of the matrix C to be at most one, it suffices to have 121 max           x t a        1max x t a      1max   x t a max 1 ax t     Same method above, resulting the system of equation defined by eq.(11) is then defined by , k s kk UUS where ,],,,[ 1,2,1, Ts mk s k s k s k uuuU   ,],,,[ 1,2,1, T mkkkk uuuU      IBN AL- HAITHAM J. FO R PURE & APPL. SC I. VO L.24 (1) 2011 kS is the matrix of coefficients, and is the sum of a lower triangular matrix and a super diagonal matrix at the grid point xk for 1,,1  nk  . Therefore the resulting matrix entries kS for 1,,2,1  mi  and 1,,1  mj  by are defined by          1,, 0,, 2,, 1,, , 1 jik ik ik ik ji g g g g S         for for for for 1 1 1     ij ij ij ij where the coefficients   y t b jiki    ,, So, and in the same way, according to the Greshgorin theorem [9], to get max 1 by t     Consistency and Convergent Analysis of the Explicit Finite Difference Approximation We note that the three difference operators used in eq.(6) each have a local truncation error with ),( tO  ),( xO  and )( yO  respectively. The )( tO  , for the time derivative term, is obtained from the classical Taylor's expansion. The )( xO  and )( yO  for the local truncation error of the fractional derivative terms was proved in [10]. Therefore, the explicit finite difference approximation is consistency. Theorem above shows that s jiyx u ,,, )(   converges to the mixed fractional derivative linearly, as )()( yOxO  . Therefore, the local truncation error of the explicit Euler method eq.(8) is )()()( yOxOtO  . This consistency of the explicit finite difference approximation together with the above result on conditional stability implies that the explicit finite difference approximation is convergent and this convergence is )( tyxO  . Numerical Example In this section, numerical example is presented which confirm our theoretical results. Example : Consider the two-dimensional fractional dispersion equation: tt exyeyx y tyxuy x tyxu x t tyxu 22225.0 6.1 6.16.1 5.1 5.1 ),,( 2 )4.1(),,( )5.0( ),,(          subject to the initial condition u (x,y,0) = xy2 , 0  x  0.5, 0 < y < 0.5 and the boundary conditions u (0,y,t) = 0, 0 < y < 0.5, 0  t 0.025 IBN AL- HAITHAM J. FOR PURE & APPL. SCI. VO L.24 (1) 2011 u (x,0,t) = 0, 0 < x < 0.5, 0  t 0.025 u (,0.5,y,t) = 0.5 e 2t , 0 < y < 0.5, 0  t 0.025 u (x,, 0.5,t) = 0.25e2tx , 0 < x < 0.5, 0  t 0.025 This fractional dispersion equation together with the above initial and boundary condition is constructed such that the exact solution is u(x,y,t) = e2txy2. Table (1) and (2) show the numerical solution using the explicit finite difference approximation. From table (1) and (2), it can be seen that thereisa good agreement between the numerical solution and exact solution. 6. Conclusions In this paper, a numerical method for solving the two-dimensional fractional dispersion equation has been described and demonstrated. The explicit difference approximation is proved to be conditionally stable and converges. Furthermore numerical example is presented to show that good agreement between the numerical solution and exact solution has been noted. References 1. Huang, F. and Liu, F. (2005), Anziam Journal, 46: 1-14. 2. Shen, S. and Lin,F. (2005),ANZIAM J., 46(E): 871-887. 3. Meerschaert, M.M.; Scheffler, H.P. and Tadjeran,C. (2006), J. Comput. Phys., 211:249–261. 4. Metzler, R. and Klafter, J. (2004), J. Physics, A 37: R161-R208. 5. Sabatelli, L.; Keating, S.; Dudley,J. and Richmond, P. (2002), Eur. Phys. J., B 27: 273–275. 6. Scalas,E.; Gorenflo, R. and Mainardi, F. (2000), Phys., A 284: 376–384. 7. Benson, D.; Wheatcraft, S. and Meerschaert, M. (2000), Water Resour. Res., 36:1403–1412. 8. Schumer,R.; Benson, D.A.; Meerschaert, M.M. and Baeumer, B. (2003), Water Resour. Res., 39: 1022–1032. 9. Isaacson, E. and Keller, H.B. (1966), Wiley, New York. 10.Meerschaert, M.M. and Tadjeran,C. (2004), J. Comput. Appl. Math., 172: 65–77. Table (1) The numerical solution of example by using the explicit finite difference approximation for 1.0,1.0  yx and 0125.0t Numerical Solution Exact Solution Error 9.730E-4 1.02532 E -3 6.72814 E -2 7.876E-3 8.20252 E -3 3.26521 E -4 0.027 2.76835 E -2 6.83508 E -4 0.064 6.56202 E -2 1.62017 E -3 9.795E-4 1.05127 E -3 7.17711 E -5 7.759E-3 8.41017 E -3 6.51169 E -4 0.027 2.83843 E -2 1.38432 E 3 0.036 6.72814 E -2 3.12814 E -2 Table (2) The numerical solution of example by using the explicit finite difference approximation for 125.0,125.0  yx and 0125.0t Numerical Solution Exact Solution Error 1.908E-3 2.00257 E -3 9.45686 E -5 0.015 1.60205 E -2 1.02055 E -3 0.052 5.40693 E -2 2.06935 E -3 1.929E-3 2.05326 E -3 1.24264 E -4 0.015 1.64261 E -2 1.42611 E -3 0.036 5.54381 E -2 1.94381 E -2 2011) 1( 24مجلة ابن الھیثم للعلوم الصرفة والتطبیقیة المجلد تتشتال دلةلمعا الفروق المنتهیة الصریحةّ تقریب ذات البعدین الكسریة ایمان ایشو كولایر جامعة بغداد ،ابن الهیثم -كلیة التربیة ،قسم الریاضیات 2010 شباط 7استلم البحث في 2010 ایار 13قبل البحث في الخالصة وارزمیــة الحـل العــددي وان خ. البعـدین ذي الكســریة تتشـتال ةلمعادلــفـي هـذا البحــث قـدمنا وناقشــنا خوارزمیـة للحــل العـددي شـرطي، والتقـارب للطریقـة السـتقرار االو ،تسـاق اال اكمـا ناقشـن. الصـریحة الفـروق المنتهیـة تقریـب قائمة على اساس دلةالتلك المع .العددیة وبینـا أن تقریـب الفـروق المنتهیـة الصـریحة ةحسـب مرتبـة االشـتقاق الكسـری التشـتتسلوك ظهرلی اعددی اخیرا قدمنا مثاال .البعدین الكسریة ذي التشتت هو طریقة فعالة حسابیا لمعادلة .التقارب االستقرار، التشتت الكسریةمعادلة ، طریقة اوبلر الصریحة، مسألة ذي بعدین، كسریة مشتقة -:الكلمات المفتاحیة