# A New Numerical Approach for Solving 1D Fractional Diffusion-Wave Equation

Authors: Umair Ali; Muhammad Asim Khan; Mostafa M. A. Khater; A. A. Mousa; Raghda A. M. Attia
Journal: Journal of Function Spaces
 (2021)

---

## Abstract

Fractional derivative is nonlocal, which is more suitable to simulate physical phenomena and provides more accurate models of physical systems such as earthquake vibration and polymers. The present study suggested a new numerical approach for the fractional diffusion-wave equation (FDWE). The fractional-order derivative is in the Riemann-Liouville (R-L) sense. Discussed the theoretical analysis of stability, consistency, and convergence. The numerical examples demonstrate that the method is more workable and excellently holds the theoretical analysis, showing the scheme’s feasibility.


---

## Body

## 1. Introduction

Nowadays, fractional calculus gains attention because of many applications in different sciences, such as physical and chemical phenomena [1–13]. Therefore, many numerical methods are used for the fractional differential equations [14–21]. The fractional diffusion-wave equations (FDWE) have many applications, for example, continuous-time random walks, unification of diffusion, and wave propagation phenomenon [22, 23]. To find the solution of FDWE, different researchers used different numerical methods; for example, Hu and Zhanga [24] implemented a higher-order difference approach on the system of the FDWE and proved its stability, solvability, and convergence. Their numerical solutions reported that the high-order difference method is more efficient than the Crank-Nicolson method. Yang et al. [25] first converted FDWE into the integrodifferential equation and then apply the fractional multisteps method to solve the FDWE. They discretized the transform equation’s integration term and Lubich’s fractional multistep approach for space derivative and central difference schemes, respectively. Also, they show that the proposed method is convergent stable and solvable for γ∈1,1.71832. Brunner et al. [26] considered the numerical simulation of the 2-D time FDWE. They introduced artificial boundary conditions to transform the problem into the initial boundary value problem (IBVP) and produced a stable and convergence scheme using finite difference approximation. Sweilam et al. [27] solved two-sided space FDWE using a finite difference method. They studied the consistency and stability of both explicit and implicit methods. Also, they concluded that finite difference methods are relatively easy to use and can apply to different fractional differential equations easily. Later, Sweilam et al. [28] solved FDWE with the help of a weighted average finite difference method, which is an extension of the weighted average for the ordinary DWE. The noninteger order derivative is discretized using RL derivative and found that the suggested approach is stable. Liu et al. [29] used the novel finite difference scheme to solve the IBVPs for the FDWE. They reported that the proposed technique is unconditionally stable and convergent, having the order of convergence Oτ3−γ+h2, where h and t represent the time step and mesh size, respectively. Dehghan and Abbaszadeh [30] applied A finite difference method to solve one and 2-D FDWE both space and tempered fractional forms. For fractional space derivative, they used Riemann-Liouville’s fractional derivative. They show using error estimate that the suggested scheme is convergence and unconditionally stable.Sun and Wu [31] researched the fully discrete scheme to solve the solution of the BVP of FDWE. They determined that the proposed numerical approach is solvable, stable, and L∞ convergent with convergence order Oτ3−γ+h2 using the energy method. Furthermore, Du et al. [32] constructed a higher-order approximation to find the solution of FDWE with the order Oτ4−γ+h2, where the fractional derivative is approximated using Caputo fractional derivative. The numerical results show that the Caputo fractional derivative gives convergence Oτ4−γ, which can be many apply to other fractional differential equations. Zhang et al. [33] proposed two new schemes, namely, the Crank-Nicolson (C-N), alternating direction implicit (ADI), and compact ADI to get a solution of 2-D FDWE. The close ADI scheme used the ADI approach for the time, and the relative difference technique is used for the spatial discretizations. Furthermore, they reported that the compact ADI approach is solvable, unconditionally stable, and convergent. The experiments confirmed that the ADI approach is computationally efficient in terms of CPU timing compared to the C-N ADI scheme for different values of γ without affecting the results. Huang et al. [34] used two approximations for the initial-boundary value time FDWEs, which is based on its integrodifferential equations. The proposed schemes were convergent with first and second order in temporal and spatial directions, respectively. Wei [35] used the semidiscrete strategy with the help of the finite difference method and fully implicit discrete scheme using the local Galerkin method for space time-fractional order derivatives, respectively, for the FDWE. The proposed approach results in a fully discrete scheme, which is convergent and unconditionally stable. Ali and Abdullah [36] formulated a numerical scheme for the fractional-order wave equation. They provided numerical solutions for the feasibility of the proposed scheme.The above-cited literature means that fractional calculus is still a novel topic and needs more efficient numerical techniques to investigate the more feasible fractional differential equations (FDEs). This study is aimed at examining the highly accurate and robust numerical approach for fractional-order diffusion-wave equation. So far, the attempt has not been made to discretize the R-L integral and implement in the Riemann-Liouville fractional order derivative to approximate FDWE.The remaining paper is organized as follows: Section2 discussed the related preliminaries. Section 3 explains the methodology of the proposed scheme. Section 4 discusses the stability analysis. In Section 5, consistency and convergence are discussed. Section 6 provides numerical results, and finally, the conclusion is presented in Section 7.Here, we consider the FDWE as following [36]:
(1)∂2Yx,t∂t2=Dt2−q∂2Yx,t∂x2+Fx,t,(2)qε1,2,xε0,L,tε0,T,where the conditions are
(3)Yx,0=v1x,Ytx,0=v2x,Y0,t=YL,t=0,where Fx,t is denoted the source term and Dt2−q represented the R-L of order q lying between 1 and 2.


## 2. Preliminaries

The R-L fractional order derivative is defined as follows:(4)Dt2−qYx,t=1Γqd2dt2∫0tYx,zt−z1−qdz=d2dt2I0qYx,t.Here, the R-L integral can be defined as:(5)I0qYx,t=1Γq∫0tYx,zt−z1−qdz.To approximate the Riemann-Liouville integral in equation (5), following the approach in [37]:
(6)I0qYx,t=1Γq∫0tntn−zq−1Yx,zdz,using the Jumarie property [38] as:
(7)=1qΓq∫0tnYx,zdzq,=1Γ1+q∑k=0n−1∫tktk+1Yx,zdzq,=1Γ1+q∑k=0n−1Yx,tn−k∫tktk+1z0dzq,applying the Jumarie property ∫0tzmdzn=Γ1+mΓ1+n/Γ1+n+mtm+n,
(8)=τqΓ1+q∑k=0n−1Yx,tn−kk+1q−kq,=τqΓ1+q∑k=0n−1CkqYx,y,tn−k,and Ckq=k+1q−kq,k=0,1,2,⋯,n−1.Lemma 1.
Theq0<q<1 order R-L fractional integral of Yx,t in 0,T defined as:
(9)I0qYxi,tn=τqΓ1+q∑k=0n−1CkqYxi,tn−k.Lemma 2.
The coefficientCkqk=0,1,⋯ satisfies the below properties [37]:
(i)
C0q=1,Ckq>0,k=0,1,⋯.(ii)
Ck−1q>Ckq,k=1,2,⋯.(iii)
Here,∃ the +ve constant C1>0. as τ≤C1Ckqτq,k=1,2,⋯.(iv)
∑k=0nCkqτq=n+1q≤Tq


## 3. The Proposed Scheme

The IDS for FDWE is developed in this section. The fractional order part is replaced by Lemma1, and space derivative is approximated by central difference approximation. The step for space is xi=iΔx, and the step for time is tk=kτ, where 1≤i≤M−1, Δx=L/M and 0≤k≤N, τ=T/N, respectively. Using equation (4) in (1) at grid point Yxi,tk, as follows
(10)∂2Yxi,tk∂t2=d2dt2I0qδx2Yxi,tkΔx2+Fxi,tk.Eliminating the second order time derivative by applying backward difference approximation w.r.t time and then integrating fromtk to tk+1 and using trapezoidal rule for source term, we obtained
(11)Yik+1−2Yik+Yik−1=1Δx2I0qδx2Yik+1−2δx2Yik+δx2Yik−1+τ22Fik+1+Fik,substituting Lemma 1 in place of RL fractional integral in equation (11), the simplified approximated difference scheme is obtained for fractional diffusion-wave equation, as follows
(12)Yik+1−2Yik+Yik−1=S1Yi+1k+1−2Yik+1+Yi−1k+1−S1Yi+1k−2Yik+Yi−1k+S1∑j=0k−1Cj+1q−CjqYi+1k−j−2Yik−j+Yi−1k−j−S1∑j=1k−1Cjq−Cj−1qYi+1k−j−2Yik−j+Yi−1k−j+τ22Fik+1+Fik,we know that Ytx,0=v2x, therefore
(13)Yi−1=Yi1+2τv2xi,(14)Yi0=v1xi,Y0k=YMk=0,where S1=τq/Γ1+qΔx2,i=1⋯M−1,k=0⋯N−1.


## 4. Stability

The Von Neumann technique is using to find the stability of the suggested scheme. Letyik represent the exact solution for equation (12), we get
(15)yik+1−2yik+yik−1=S1yi+1k+1−2yik+1+yi−1k+1−S1yi+1k−2yik+yi−1k+S1∑j=0k−1Cj+1q−Cjqyi+1k−j−2yik−j+yi−1k−j−S1∑j=1k−1Cjq−Cj−1qyi+1k−j−2yik−j+yi−1k−j.The error is defined asEik=yik−Yik. And Eik satisfies equation (15) as follows:
(16)Eik+1−2Eik+Eik−1=S1Ei+1k+1−2Eik+1+Ei−1k+1−S1Ei+1k−2Eik+Ei−1k+S1∑j=0k−1Cj+1q−CjqEi+1k−j−2Eik−j+Ei−1k−j−S1∑j=1k−1Cjq−Cj−1qEi+1k−j−2Eik−j+Ei−1k−j.Consider the growth factor in the form of a one Fourier mode as(17)Eik=λke−1iσΔx,where σ and Δx are the mode number and step size, respectively. Equation (17) can be the solution of the above error equation (16), as
(18)λk+1e−1iσΔx−2λke−1iσΔx+λk−1e−1iσΔx=S1λk+1e−1i+1σΔx−2λk+1e−1iσΔx+λk+1e−1i−1σΔx−S1λke−1i+1σΔx−2λke−1iσΔx+λke−1i−1σΔx+S1∑j=0k−1Cj+1q−Cjqλk−je−1i+1σΔx−2λk−je−1iσΔx+λk−je−1i−1σΔx−S1∑j=1k−1Cjq−Cj−1qλk−je−1i+1σΔx−2λk−je−1iσΔx+λk−je−1i−1σΔx.Dividing both sides bye−1iσΔx, then replace e−1σΔx+e−−1iσΔx=2−4S1sin2σΔx/2, we obtained
(19)λk+1=11+4S1sin2σΔx/22+4S1sin2σΔx2λk−λk−1−4S1sin2σΔx2∑j=0k−1Cj+1q−Cjqλk−j−∑j=1k−1Cjq−Cj−1qλk−j.Proposition 3.
Supposeλk+1,k=0,1,⋯,N−1 is the solution of equation (19), then we need to prove that
(20)∣λk+1∣≤λ0.Proof.
Let takek=0 in equation (19), to prove the proposition by induction method, as
(21)λ1=λ02+4S1sin2σΔx/2−λ−11+4S1sin2σΔx/2.By utilizing equations (13) and (17), we get λ−1=λ1 and putting in equation (21). After simplification, the above equation becomes
(22)λ1=λ02+4S1sin2σΔx/22+4S1sin2σΔx/2,obtained the following relation
(23)λ1≤λ0.Let∣λ1∣≤∣λ0∣ holds for k=0,1,⋯,n−1.Using equations (19) and (23) and Lemma 2, we have
(24)λk+1=11+4S1sin2σΔx/22+4S1sin2σΔx2λk−λk−1−4S1sin2σΔx2∑j=0k−1Cj+1q−Cjqλk−j−∑j=1k−1Cjq−Cj−1qλk−j,λk+1≤11+4S1sin2σΔx/22+4S1sin2σΔx2λ0−λ0−4S1sin2σΔx2∑j=0k−1Cj+1q−Cjqλ0−∑j=1k−1Cjq−Cj−1qλ0,≤11+4S1sin2σΔx/22+4S1sin2σΔx2−1−4S1sin2σΔx2∑j=0k−1Cj+1q−Cjq−∑j=1k−1Cjq−Cj−1qλ0,≤11+4S1sin2σΔx/21+4S1sin2σΔx2−4S1sin2σΔx2C1q−Ck−1q−C1q−Ck−2qλ0,≤11+4S1sin2σΔx/21+4S1sin2σΔx2−4S1sin2σΔx2Ck−2q−Ck−1qλ0.From Lemma2, the value 0<Ck−2q−Ck−1q<1, so it is clear that the value
(25)0<1+4S1sin2σΔx/2−4S1sin2σΔx/2Ck−2q−Ck−1q1+4S1sin2σΔx/2<1,λk+1≤λ0.Here,λk+1≤λ0, also Eik+1≤Ei0, so it can written as Eik+12≤Ei02. It shows that the suggested approach is unconditionally stable.


## 5. Consistency

To find the consistency analysis, supposew is the closed form solution, and W is the approximated solution, and the function YW=0 is the approximated scheme for the proposed equation at the mesh point xi,tk. Then, YW=Tik denoted the local truncation error at xi,tk.Theorem 4.
The local transcation errorTx,t for the suggested scheme is Tik=OΔt2+OΔx2.Proof.
The obtained difference scheme in equation (12), can also be written as:
(26)Tik+1=Yik+1−2Yik+Yik−1−S1∑j=0kCjqui+1k−j+1−2uik−j+1+ui−1k−j+1−2∑j=0k−1Cjqui+1k−j−2uik−j+ui−1k−j+∑j=0k−2Cjqui+1k−j−1−2uik−j−1+ui−1k−j−1,using Taylor series, we obtained
(27)Tik+1=uik+Δt∂u∂tik+Δt22!∂2u∂t2ik+Δt33!∂3u∂t3ik+⋯−2uik−uik−Δt∂u∂tik+Δt22!∂2u∂t2ik−Δt33!∂3u∂t3ik+⋯−S1∑j=0kCjquik−j+1+Δx∂u∂xik−j+1+Δx22!∂2u∂x2ik−j+1+Δx33!∂3u∂x3ik−j+1+⋯−2uik−j+1+uik−j+1−Δx∂u∂xik−j+1+Δx22!∂2u∂x2ik−j+1−Δx33!∂3u∂x3ik−j+1+⋯−2S1∑j=0k−1Cjquik−j+Δx∂u∂xik−j+Δx22!∂2u∂x2ik−j+Δx33!∂3u∂x3ik−j+⋯−2uik−j+uik−j−Δx∂u∂xik−j+Δx22!∂2u∂x2ik−j−Δx33!∂3u∂x3ik−j+⋯−S1∑j=0k−2Cjquik−j−1+Δx∂u∂xik−j−1+Δx22!∂2u∂x2ik−j−1+Δx33!∂3u∂x3ik−j−1+⋯−2uik−j−1+uik−j−1−Δx∂u∂xik−j−1+Δx22!∂2u∂x2ik−j−1−Δx33!∂3u∂x3ik−j−1+⋯,Tik+1=Δt2∂2u∂t2ik+Δt412∂4u∂t4ik+⋯+S1∑j=0kCjqΔx2∂2u∂x2ik−j+1+Δx412∂4u∂x4ik−j+1+⋯−2S1∑j=0k−1CjqΔx2∂2u∂x2ik−j+Δx412∂4u∂x4ik−j+⋯−S1∑j=0k−2CjqΔx2∂2u∂x2ik−j−1+Δx412∂4u∂x4ik−j−1+⋯.The above equation can also be written as(28)Tik+1=Δt2∂2u∂t2ik+Δt212∂4u∂t4ik+⋯+Δx2S1∑j=0kCjq∂2u∂x2ik−j+1−2S1∑j=0k−1Cjq∂2u∂x2ik−j+S1∑j=0k−2Cjq∂2u∂x2ik−j−1.We obtained the following truncation error(29)Tik+1=OΔt2+OΔx2.The IDS of fractional diffusion-wave equation is consistent if space stepΔx and time step Δt approaches zero; then, truncation error approaches zero.Theorem 5.
The consistency and stability analysis are both the necessary and the sufficient for convergence. By Lax-equivalence theorem [33], the proposed scheme is convergent.


## 6. Numerical Results

Here, the FDWE examples reported to find the accuracy and feasibility of the scheme. The numerical examples are coded in Maple 15, i.e., the maximum error, is defined as follows:(30)E∞=max0≤i≤M−1,0≤k≤Nyxi,tk−Yik.TheE2 error is:
(31)E2=∑i=1M−1yxi,tk−Yi,jk2Δx1/2.The convergence order is:(32)C−order=loge1/e2logh1/h2.Example 6.
In equations (1)–(3), take the source term Fx,t=Γ3+q/Γ1+qexx21−x2tq−Γ3+q/Γ1+2qext2q2−8x+x2+6x3+x4 and the closed solution Yx,t=exx21−x2tq+2.Example 7.
In equations (1)–(3) taking the source term Fx,t=exΓ2+q/Γqtq−1−Γ2+q/Γ2qt2q−1 and the closed solution ux,t=ext1+q.


## 7. Discussion

The numerical examples 5.1 and 5.2 are solved by the suggested numerical scheme in equation (12). Tables 1–4 represent the performance of the numerical results. The accuracy level increases as reducing the space and the time step; please refer to Tables 1 and 2 for different and same space and time steps, respectively. In Table 3, fixed Δx=1/100 and for various values of τ and q, and find the convergence order completely demonstrates the analysis theoretically. Similarly, Table 4 also shows the accuracy for the second test example 5.2, which converges to the exact solution to reduce the error by reducing the step size. From the above discussion, it is clear that the proposed scheme proved our theoretical analysis.Table 1
The numerical results of the IDS for Example6 at T=1.0.
τ=Δx2q=0.25q=0.5q=0.751/161.9802×10−21.9053×10−21.8120×10−21/645.0811×10−34.8647×10−34.6537×10−31/1442.2611×10−33.1168×10−32.0756×10−31/2561.1894×10−31.1418×10−31.1693×10−3Table 2
The numerical results of the proposed scheme for Example6 at T=1.0.
τ=Δxq=0.25q=0.5q=0.751/104.9799×10−34.8584×10−34.6045×10−31/309.3870×10−48.3367×10−48.0913×10−41/603.7440×10−43.1316×10−43.1540×10−41/902.0584×10−41.8973×10−41.0265×10−41/1201.7578×10−41.3173×10−49.7360×10−5Table 3
Order of convergence of the IDS for Example6 at q=0.5,Δx=1/100,T=1.0.
τE∞C−orderE2C−order1/51.6672×10−2—1.3142×10−2—1/104.5646×10−31.87913.1168×10−32.07621/201.2792×10−31.85437.7861×10−42.00121/404.0253×10−41.78232.0522×10−41.9239Forq=0.751/51.5879×10−2—1.2531×10−2—1/104.3893×10−31.86522.9858×10−32.06951/201.2362×10−31.82828.5184×10−41.81571/403.4073×10−41.85942.1507×10−41.9926Forq=0.951/51.4832×10−2—1.1722×10−2—1/104.2321×10−31.86522.8672×10−32.06951/201.2239×10−31.82828.4446×10−41.81571/403.4853×10−41.81832.3648×10−41.8427Table 4
The numerical results of the proposed scheme for Example7 at T=1.0.
τ=Δx2q=0.5q=0.751/161.6221×10−23.0081×10−21/644.3498×10−37.8153×10−31/1441.9566×10−33.4826×10−31/3248.6632×10−49.7248×10−4


## 8. Conclusion

A practical and fast numerical scheme has been developed for FDWE. The approximation is based on the discretization of Riemann-Liouville integral in Lemma 1. We have successfully proved the analysis theoretically of stability by mathematical induction, consistency, and convergence. The numerical result confirmed our theoretical analysis and demonstrated that the proposed scheme is fast convergent and more feasible. This approach can also apply to other types of FDWEs.



---

*Source: 6638597-2021-01-06.xml*