Finite Difference and Iteration Methods for Fractional Hyperbolic Partial Differential Equations with the Neumann Condition

Abstract

The numerical and analytic solutions of the mixed problem for multidimensional fractional hyperbolic partial differential equations with the Neumann condition are presented. The stable difference scheme for the numerical solution of the mixed problem for the multidimensional fractional hyperbolic equation with the Neumann condition is presented. Stability estimates for the solution of this difference scheme and for the first- and second-order difference derivatives are obtained. A procedure of modified Gauss elimination method is used for solving this difference scheme in the case of one-dimensional fractional hyperbolic partial differential equations. He's variational iteration method is applied. The comparison of these methods is presented.

1. Introduction

It is known that various problems in fluid mechanics (dynamics, elasticity) and other areas of physics lead to fractional partial differential equations. Methods of solutions of problems for fractional differential equations have been studied extensively by many researchers (see, e.g., [1–15] and the references given therein).

The role played by stability inequalities (well posedness) in the study of boundary-value problems for hyperbolic partial differential equations is well known (see, e.g., [16–29]).

In the present paper, finite difference and He's iteration methods for the approximate solutions of the mixed boundary-value problem for the multidimensional fractional hyperbolic equation 𝜕²𝑢(𝑡,𝑥)𝜕𝑡²−^𝑚_𝑟=1𝑎_𝑟𝑥=(𝑥𝑥)₁𝑢_𝑥𝑟_𝑥𝑟+𝐷^𝑡1/2𝑢(𝑡,𝑥)+𝜎𝑢(𝑡,𝑥)=𝑓(𝑡,𝑥),𝑢(0,𝑥)=0,𝑢_𝑡,…,𝑥_𝑚∈Ω,0<𝑡<1,(0,𝑥)=0,𝑥∈Ω;𝜕𝑢(𝑡,𝑥)𝜕𝑛=0,𝑥∈𝑆,(1.1) are studied. Here Ω is the unit open cube in the 𝑚-dimensional Euclidean space: ℝ^𝑚∶{Ω=𝑥=(𝑥₁,…,𝑥_𝑚)∶0<𝑥_𝑗<1,1≤𝑗≤𝑚} with boundary 𝑆,Ω=Ω∪𝑆; 𝑎_𝑟(𝑥)(𝑥∈Ω) and 𝑓(𝑡,𝑥)(𝑡∈(0,1),𝑥∈Ω) are given smooth functions and 𝑎_𝑟(𝑥)≥𝑎>0.

1.1. Definition

The Caputo fractional derivative of order 𝛼>0 of a continuous function 𝑢(𝑡,𝑥) is defined by 𝐷^𝛼𝑎₊𝑢(𝑡,𝑥)=1Γ(1−𝛼)^𝑡𝑎𝑢^(𝑡,𝑥)(𝑡−𝑠)^𝛼𝑑𝑠,(1.2) where Γ(⋅) is the gamma function.

2. The Finite Difference Method

In this section, we consider the first order of accuracy in 𝑡 and the second-orders of accuracy in space variables’ stable difference scheme for the approximate solution of problem (1.1). The stability estimates for the solution of this difference scheme and its first- and second-order difference derivatives are established. A procedure of modified Gauss elimination method is used for solving this difference scheme in the case of one-dimensional fractional hyperbolic partial differential equations.

2.1. The Difference Scheme: Stability Estimates

The discretization of problem (1.1) is carried out in two steps. In the first step, let us define the grid space Ω_ℎ=𝑥=𝑥_𝑟=ℎ₁𝑟₁,…,ℎ_𝑚Ω𝑟_𝑚,𝑟=𝑟₁,…,𝑟_𝑚,0≤𝑟_𝑗≤𝑁_𝑗,ℎ_𝑗𝑁_𝑗=1,𝑗=1,…,𝑚_ℎ,=Ω_ℎ∩Ω,𝑆_ℎ=Ω_ℎ∩𝑆.(2.1) We introduce the Banach space 𝐿_2ℎ=𝐿₂(Ω_ℎ) of the grid functions 𝜑^ℎ(𝑥)={𝜑(ℎ₁𝑟₁,…,ℎ_𝑚𝑟_𝑚)} defined on Ω_ℎ, equipped with the norm ‖‖𝜑^ℎ‖‖_𝐿2(Ω_ℎ)=⎛⎜⎜⎝_𝑥∈Ωℎ||𝜑^ℎ(𝑥)||²ℎ₁⋯ℎ_𝑚⎞⎟⎟⎠^1/2.(2.2) To the differential operator 𝐴^𝑥 generated by problem (1.1), we assign the difference operator 𝐴^𝑥ℎ by the formula 𝐴^𝑥ℎ𝑢^ℎ=−^𝑚_𝑟=1𝑎^𝑟(𝑥)𝑢^ℎ_−𝑥𝑟_{𝑥𝑟,𝑗}+𝜎𝑢^ℎ𝑥(2.3) acting in the space of grid functions 𝑢^ℎ(𝑥), satisfying the conditions 𝐷_ℎ𝑢^ℎ(𝑥)=0 for all 𝑥∈𝑆_ℎ. It is known that 𝐴^𝑥ℎ is a self-adjoint positive definite operator in 𝐿₂(Ω_ℎ). With the help of 𝐴^𝑥ℎ we arrive at the initial boundary value problem 𝑑²𝑣^ℎ(𝑡,𝑥)𝑑𝑡²+𝐷𝑣^𝑡1/2_ℎ𝑣^ℎ(𝑡,𝑥)+𝐴^𝑥ℎ𝑣^ℎ(𝑡,𝑥)=𝑓^ℎ(𝑡,𝑥),0≤𝑡≤1,𝑥∈Ω_ℎ,(0,𝑥)=0,𝑑𝑣^ℎ(0,𝑥)𝑑𝑡=0,𝑥∈Ω,(2.4) for an infinite system of ordinary fractional differential equations.

In the second step, we replace problem (2.4) by the first order of accuracy difference scheme 𝑢^ℎ𝑘+1(𝑥)−2𝑢^ℎ𝑘(𝑥)+𝑢^ℎ𝑘−1(𝑥)𝜏²+1√𝜋^𝑘_𝑚=1Γ(𝑘−𝑚+1/2)(𝑘−𝑚)!𝑢^ℎ𝑚−𝑢^ℎ_m−1𝜏^1/2𝑓_ℎ𝑘+𝐴^𝑥ℎ𝑢^ℎ𝑘+1=𝑓^ℎ𝑘(𝑥),𝑥∈Ω_ℎ,(𝑥)=𝑓𝑢_ℎ1𝑡_𝑘,𝑥,𝑡_𝑘=𝑘𝜏,1≤𝑘≤𝑁−1,𝑁𝜏=1,𝑥∈Ω_ℎ,(𝑥)−𝑢^ℎ0(𝑥)𝜏=0,𝑢^ℎ0(𝑥)=0,𝑥∈Ω_ℎ.(2.5) Here Γ(𝑘−𝑚+1/2)=∫^∞0𝑡^{𝑘−𝑚−1/2}𝑒^−𝑡𝑑𝑡.

Theorem 2.1. Let 𝜏 and |ℎ|=ℎ²¹+⋯+ℎ^2𝑚 be sufficiently small numbers. Then, the solutions of difference scheme (2.5) satisfy the following stability estimates: max_{1≤𝑘≤𝑁}‖‖𝑢^ℎ𝑘‖‖_𝐿2ℎ+max_{1≤𝑘≤𝑁}‖‖‖‖𝑢^ℎ𝑘−𝑢^ℎ𝑘−1𝜏‖‖‖‖_𝐿2ℎmax_{1≤𝑘≤𝑁−1}≤𝐶₁max_{1≤𝑘≤𝑁−1}‖‖𝑓^ℎ𝑘‖‖_𝐿2ℎ,‖‖𝜏⁻²𝑢^ℎ𝑘+1−2𝑢^ℎ𝑘+𝑢^ℎ𝑘−1‖‖_𝐿2ℎ+max_{1≤𝑘≤𝑁}‖‖𝑢^ℎ𝑘⁻_{𝑥𝑟𝑥𝑟}‖‖_𝐿2ℎ≤𝐶₂‖‖𝑓^ℎ1‖‖_𝐿2ℎ+max_{2≤𝑘≤𝑁−1}‖‖𝜏⁻¹𝑓^ℎ𝑘−𝑓^ℎ𝑘−1‖‖_𝐿2ℎ.(2.6)
Here 𝐶₁ and 𝐶₂ do not depend on 𝜏, ℎ, and 𝑓^ℎ𝑘,1≤𝑘≤𝑁−1.

The proof of Theorem 2.1 is based on the self-adjointness and positive definitness of operator 𝐴^𝑥ℎ in 𝐿_2ℎ and on the following theorem on the coercivity inequality for the solution of the elliptic difference problem in 𝐿_2ℎ.

Theorem 2.2. For the solutions of the elliptic difference problem 𝐴^𝑥ℎ𝐷𝑢^ℎ_ℎ(𝑥)=𝜔^ℎ(𝑥),𝑥∈Ω_ℎ,𝑢^ℎ(𝑥)=0,𝑥∈𝑆_ℎ,(2.7) the following coercivity inequality holds [30]: ^𝑚_𝑟=1‖‖‖𝑢_{ℎ𝑥𝑟−𝑥𝑟}‖‖‖_𝐿2ℎ≤𝐶‖‖𝜔^ℎ‖‖_𝐿2ℎ.(2.8)

Finally, applying this difference scheme, the numerical methods are proposed in the following section for solving the one-dimensional fractional hyperbolic partial differential equation. The method is illustrated by numerical examples.

2.2. Numerical Results

For the numerical result, the mixed problem 𝑓(𝑡𝐷,𝑥^2𝑡)𝑢=(𝑡,𝑥)2++𝐷𝑡^𝑡1/2₂𝑢(𝑡,𝑥)−𝑢_𝑥𝑥(𝑡,𝑥)+𝑢(𝑡,𝑥)=𝑓(𝑡,𝑥),+8𝑡^3/2𝑢(0,𝑥)=03,√𝜋𝑢_𝑡+(𝜋𝑡)²cos(𝜋𝑥),0<𝑡,𝑥<1,𝑢(0,𝑥)=0,0≤𝑥≤1,_𝑥(𝑡,0)=𝑢_𝑥(𝑡,1)=0,0≤𝑡≤1(2.9) for solving the one-dimensional fractional hyperbolic partial differential equation is considered. Applying difference scheme (2.5), we obtained 𝑢^𝑛𝑘+1−2𝑢^𝑘𝑛+𝑢^𝑛𝑘−1𝜏²+1√𝜋^𝑘_𝑚=1Γ(𝑘−𝑚+1/2)(𝑘−𝑚)!𝑢^𝑚𝑛−𝑢^𝑛𝑚−1𝜏^1/2−𝑢^𝑘+1𝑛+1−2𝑢^𝑛𝑘+1+𝑢^{𝑘+1𝑛−1}ℎ²+𝑢^𝑘𝑛𝜑_𝑘𝑛=𝜑^𝑘𝑛,=𝑢𝑓𝑡_𝑘,𝑥_𝑛,1≤𝑘≤𝑁−1,1≤𝑛≤𝑀−1,_0𝑛𝑢=0,𝜏⁻¹𝑢^1𝑛−𝑢^0𝑛=0,0≤𝑛≤𝑀,_𝑘1−𝑢^𝑘0=𝑢^𝑘𝑀−𝑢^𝑘𝑀−1=0,0≤𝑘≤𝑁.(2.10) We get the system of equations in the matrix form: 𝐴𝑈_𝑛+1+𝐵𝑈𝑈_𝑛+𝐶𝑈_𝑛−1=𝐷𝜑_𝑛,1≤𝑛≤𝑀−1,₁=𝑈₀,𝑈_𝑀=𝑈_𝑀−1,(2.11) where 𝐵=⎡⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣𝑏_1,1𝐴=⎡⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⋯0000000⋯0000000⋯𝑎000000⋯𝑎000000⋯⋯⋯⋯⋯⋯⋯⋯⋯𝑎000000⋯0000000⎤⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦_{(𝑁+1)×(𝑁+1)},𝑏000⋯00_2,1𝑏𝑏_2,200⋯00_3.1𝑏𝑏_3,2𝑏_3,30⋯00_4,1𝑏⋯𝑏⋯_4,2𝑏⋯_4,3𝑏⋯_4,4⋯⋯⋯0⋯_𝑁,10𝑏_𝑁+1,1𝑏_𝑁,2𝑏_𝑁,3𝑏_𝑁,4⋯𝑏_𝑁,𝑁0𝑏_𝑁+1,2𝑏_𝑁+1,3𝑏_𝑁+1,4⋯𝑏_𝑁+1,𝑁𝑏_{𝑁+1,𝑁+1}⎤⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦_{(𝑁+1)×(𝑁+1)}𝐷=⎡⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⋯0000000⋯0000000⋯0001000⋯0000100𝐶=⋯⋯⋯⋯⋯⋯⋯⋯𝐴,⋯0000010⋯0000001⎤⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦_{(𝑁+1)×(𝑁+1)},𝑈_𝑠,=⎡⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣𝑈^0𝑠𝑈^1𝑠𝑈^2𝑠𝑈^3𝑠𝑈⋯_𝑠𝑁−1𝑈^𝑁𝑠⎤⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦_{(𝑁+1)×(1)},𝑠=𝑛−1,𝑛,𝑛+1.(2.12) Here 𝑎=−1ℎ²,𝑏_1,1=1,𝑏_2,1=−1,𝑏_2,2=1,𝑏_3,1=1𝜏²−1𝜏^1/2𝑏_3,2,=−2𝜏²+1𝜏^1/2,𝑏_3,3=1+1𝜏²+2ℎ²𝑏_𝑘+2,1,=−1√𝜋Γ(𝑘Γ−(𝑘1)𝜏+1^1/2/2)𝑏,2≤𝑘≤𝑁−1,_{𝑘+2,𝑘+1}=−2𝜏²+1𝜏^1/2𝑏,1≤𝑘≤𝑁−1,_𝑘+2,𝑘=1𝜏²+1√𝜋Γ(1Γ+(20).5)−ΓΓ(0(1.5))𝜏1^1/2𝑏,2≤𝑘≤𝑁−1,_{𝑘+2,𝑘+2}=1+1𝜏²+2ℎ²𝑏,1≤𝑘≤𝑁−1,_{𝑘+2,𝑖+1}=1√𝜋ΓΓ((𝑘𝑘−−𝑖(𝑖+−11/2)))−ΓΓ(𝑘(𝑘−−(𝑖(𝑖+−11))+−11/2))𝜏1^1/2𝜑,3≤𝑘≤𝑁−1,1≤𝑖≤𝑘−2,_𝑘𝑛=2+(𝑘𝜏)²+8(𝑘𝜏)^3/2𝜑3√𝜋+(𝜋𝑘𝜏)²cos𝜋(𝑛ℎ),_𝑛=⎡⎢⎢⎢⎢⎢⎢⎢⎢⎣𝜑^0𝑛𝜑^1𝑛𝜑^2𝑛𝜑⋯_𝑁𝑛⎤⎥⎥⎥⎥⎥⎥⎥⎥⎦_(𝑁+1)×1.(2.13) So, we have the second-order difference equation with respect to 𝑛 matrix coefficients. To solve this difference equation, we have applied a procedure of modified Gauss elimination method for difference equation with respect to 𝑘 matrix coefficients. Hence, we seek a solution of the matrix equation in the following form: 𝑈_𝑗=𝛼_𝑗+1𝑈_𝑗+1+𝛽_𝑗+1,(2.14)𝑛=𝑀−1,…,2,1,  𝛼_𝑗(𝑗=1,…,𝑀) are (𝑁+1)×(𝑁+1) square matrices, and 𝛽_𝑗(𝑗=1,…,𝑀) are (𝑁+1)×1 column matrices defined by 𝛽𝛼_𝑛+1=𝐵+𝐶𝛼_𝑛⁻¹(−𝐴),_𝑛+1=𝐵+𝐶𝛼_𝑛⁻¹𝐷𝜑_𝑛−𝐶𝛽_𝑛,𝑛=2,3,…,𝑀,(2.15) where 𝛼₁=⎡⎢⎢⎢⎢⎢⎢⎢⎢⎣⋯1000𝛽⋯₁0100⋯0010⋯⋯⋯⋯⋯⋯0001⎤⎥⎥⎥⎥⎥⎥⎥⎥⎦_{(𝑁+1)×(𝑁+1)},=⎡⎢⎢⎢⎢⎢⎢⎢⎢⎣⋯0000⎤⎥⎥⎥⎥⎥⎥⎥⎥⎦_(𝑁+1)×1.(2.16)

Now, we will give the results of the numerical analysis. First, we give an estimate for the constants 𝐶₁ and 𝐶₂ figuring in the stability estimates of Theorem 2.1. We have 𝐶_𝑡1𝐶₁=max_𝑓,𝑢𝐶_𝑡1,𝐶₂=max_𝑓,𝑢𝐶_𝑡2,=max_{1≤𝑘≤𝑁}‖‖𝑢^ℎ𝑘‖‖_𝐿2ℎ+max_{1≤𝑘≤𝑁}‖‖𝜏⁻¹𝑢^ℎ𝑘−𝑢^ℎ𝑘−1‖‖_𝐿2ℎ×max_{1≤𝑘≤𝑁−1}‖‖𝑓^ℎ𝑘‖‖_𝐿2ℎ⁻¹,𝐶_𝑡2=max_{1≤𝑘≤𝑁−1}‖‖𝜏⁻²𝑢^ℎ𝑘+1−2𝑢^ℎ𝑘+𝑢^ℎ𝑘−1‖‖_𝐿2ℎ+max_{𝑛1≤𝑘≤𝑁}_𝑟=1‖‖𝑢^ℎ𝑘_{𝑥𝑟,𝑥𝑟}‖‖_𝐿2ℎ×max_{2≤𝑘≤𝑁−1}‖‖𝜏⁻¹𝑓^ℎ𝑘−𝑓^ℎ𝑘−1‖‖_𝐿2ℎ+‖‖𝑓^ℎ1‖‖_𝐿2ℎ⁻¹.(2.17) The constants 𝐶_𝑡1 and 𝐶_𝑡2 in the case of numerical solution of initial-boundary value problem (2.9) are computed. The constants 𝐶_𝑡1 and 𝐶_𝑡2 are given in Table 1 for 𝑁=20,40,80, and 𝑀=80, respectively.

Table 1 

Stability estimates for (2.9).

Second, for the accurate comparison of the difference scheme considered, the errors computed by 𝐸₀=max_{1≤𝑘≤𝑁−1}^𝑀−1_𝑛=1𝐸₁||𝑢(𝑡_𝑘,𝑥_𝑛)−𝑢^𝑘𝑛||²ℎ^1/2,=max_{1≤𝑘≤𝑁−1}^𝑀−1_𝑛=1||||𝑢_𝑡(𝑡_𝑘,𝑥_𝑛)−(𝑢^𝑛𝑘+1−𝑢^𝑛𝑘−1)2𝜏||||²𝐸₂ℎ^1/2,=max_{1≤𝑘≤𝑁−1}^𝑀−1_𝑛=1||||𝑢_𝑡𝑡(𝑡_𝑘,𝑥_𝑛)−(𝑢^𝑛𝑘+1−2𝑢^𝑘𝑛+𝑢^𝑛𝑘−1)𝜏²||||²ℎ^1/2(2.18) of the numerical solutions are recorded for higher values of 𝑁=𝑀, where 𝑢(𝑡_𝑘,𝑥_𝑛) represents the exact solution and 𝑢^𝑘𝑛 represents the numerical solution at (𝑡_𝑘,𝑥_𝑛). The errors 𝐸₀,𝐸₁ and 𝐸₂ results are shown in Table 2 for 𝑁=20,40,60 and 𝑀=60, respectively.

Table 2 

Comparison of the errors for the difference scheme.

The figure of the difference scheme solution of (2.9) is given by the Figure 2. The exact solution of (2.9) is given by as follows: 𝑢(𝑡,𝑥)=𝑡²cos(𝜋𝑥).(2.19) The figure of the exact solution of (2.9) is shown by the Figure 1.

Figure 1 

The surface shows the exact solution 𝑢(𝑡,𝑥) for (2.9).

Figure 2 

Difference scheme solution for (2.9).

3. He's Variational Iteration Method

In the present paper, the mixed boundary value problem for the multidimensional fractional hyperbolic equation (1.1) is considered. The correction functional for (1.1) can be approximately expressed as follows: 𝑢_𝑛+1(𝑡,𝑥)=𝑢+_𝑛(𝑡,𝑥)_𝑡0𝜆𝜕²𝑢(𝑠,𝑥)𝜕𝑠²−^𝑚_𝑟=1𝑎_𝑟(𝑥)𝑢̃_𝑥𝑟_𝑥𝑟+D^𝑠1/2𝑢̃(𝑠,𝑥)+𝜎𝑢̃(𝑠,𝑥)−𝑓(𝑠,𝑥)𝑑𝑠,(3.1) where 𝜆 is a general Lagrangian multiplier (see, e.g., [31]) and 𝑢̃ is considered as a restricted variation as a restricted variation (see, e.g., [32]); that is, 𝛿𝑢̃=0, 𝑢₀(𝑡,𝑥) is its initial approximation. Using the above correction functional stationary and noticing that 𝛿𝑢̃=0, we obtain 𝛿𝑢_𝑛+1(𝑡,𝑥)=𝛿𝑢_𝑛(𝑡,𝑥)+^𝑡0𝛿𝜆𝜕𝑢^2𝑛(𝑡,𝑥)𝛿𝑢_𝑛+1𝜕𝑠²𝑑𝑠,(𝑡,𝑥)=𝛿𝑢_𝑛(𝑡,𝑥)−𝜕𝜆𝜕𝑠𝛿𝑢^𝑛(𝑠,𝑥)|||_𝑠=𝑡+𝜆𝜕𝜕𝑠𝛿𝑢^𝑛(𝑠,𝑥)|||_𝑠=𝑡+^𝑡0𝜕²𝜆(𝑡,𝑠)𝜕𝑠²𝛿𝑢_𝑛(𝑠,𝑥)𝑑𝑠=0.(3.2) From the above relation for any 𝛿𝑢_𝑛, we get the Euler-Lagrange equation: 𝜕𝜆²(𝑡,𝑠)𝜕𝑠²=0,(3.3) with the following natural boundary conditions: 1−𝜕𝜆(𝑡,𝑠)𝜕𝑠|||_𝑠=𝑡𝜆(𝑡,𝑠)||_𝑠=𝑡=0,=0.(3.4) Therefore, the Lagrange multiplier can be identified as follows: 𝜆(𝑡,𝑠)=𝑠−𝑡.(3.5) Substituting the identified Lagrange multiplier into (3.1), the following variational iteration formula can be obtained: 𝑢_𝑛+1(𝑡,𝑥)=𝑢_𝑛(𝑡,𝑥)+^𝑡0(𝑠−𝑡)𝜕²𝑢_𝑛(𝑠,𝑥)𝜕𝑠²−^𝑚_𝑟=1𝑎^𝑟(𝑥)𝑢_𝑛𝑥𝑟_𝑥𝑟+𝐷^𝑠1/2𝑢_𝑛(𝑠,𝑥)+𝑢_𝑛(𝑠,𝑥)−𝑓(𝑠,𝑥)𝑑𝑠(3..6)

In this case, let an initial approximation 𝑢₀(𝑡,𝑥)=𝑢(0,𝑥)+𝑡𝑢_𝑡(0,𝑥). Then approximate solution takes the form 𝑢(𝑡,𝑥)=lim_𝑛→∞u_𝑛(𝑡,𝑥).

3.1. Variational Iteration Solution 1

For the numerical result, the mixed problem 𝑓(𝑡𝐷,𝑥^2𝑡)𝑢=(𝑡,𝑥)2++𝐷𝑡^𝑡1/2₂𝑢(𝑡,𝑥)−𝑢_𝑥𝑥(𝑡,𝑥)+𝑢(𝑡,𝑥)=𝑓(𝑡,𝑥),+8𝑡^3/2𝑢(0,𝑥)=03,√𝜋𝑢_𝑡+(𝜋𝑡)²cos(𝜋𝑥),0<𝑡,𝑥<1,𝑢(0,𝑥)=0,0≤𝑥≤1,_𝑥(𝑡,0)=𝑢_𝑥(𝑡,1)=0,0≤𝑡≤1(3.7) for solving the one-dimensional fractional hyperbolic partial differential equation is considered.

According to formula (3.6), the iteration formula for (3.7) is given by 𝑢_𝑛+1(𝑡,𝑥)=𝑢+_𝑛(𝑡,𝑥)_𝑡0(𝑠−𝑡)𝜕𝑢^2𝑛(𝑠,𝑥)𝜕𝑠²+𝐷^𝑠1/2𝑢_𝑛(𝑠,𝑥)−𝜕𝑢^2𝑛(𝑠,𝑥)𝜕𝑥²+𝑢_𝑛(𝑠,𝑥)−𝑓(𝑠,𝑥)𝑑𝑠.(3.8) Now we start with an initial approximation 𝑢₀(𝑡,𝑥)=𝑢(0,𝑥)+𝑡𝑢_𝑡(0,𝑥).(3.9) Using the above iteration formula (3.8), we can obtain the other components as 𝑢𝑢₀(𝑡,𝑥)=0,₁(𝑡,𝑥)=1420√𝜋128𝑡^7/2𝑢+35𝑡⁴√𝜋+35𝑡⁴𝜋^5/2+420𝑡²√𝜋cos(𝜋𝑥),₂(𝑡,𝑥)=1420√𝜋cos(𝜋𝑥)128𝑡^7/2+cos(𝜋𝑥)[−0.90580036+663𝑡5₄𝑡⁴√𝜋+35𝑡⁴𝜋^5/2+420𝑡²√𝜋−0.1510268880𝑡^11/2−0.1719434921𝑡⋮^7/2−0.3281897218𝑡⁶−0.01666666667𝑡⁵](3.10) The figure of (3.10) is given by the Figure 3.

Figure 3 

Variational iteration method for (3.10).

3.2. Variational Iteration Solution 2

For the numerical result, the mixed problem 𝐷𝑓^2𝑡(𝑢𝑡,(𝑥𝑡,,𝑥𝑦,)𝑦=)+𝐷2+^𝑡1/2𝑡𝑢₂(𝑡,𝑥,𝑦)−𝑢_𝑥𝑥(𝑡,𝑥,𝑦)−𝑢_𝑦𝑦(𝑡,𝑥,𝑦)+𝑢(𝑡,𝑥,𝑦)=𝑓(𝑡,𝑥,𝑦),+8𝑡^3/2𝑢(0,𝑥,𝑦)=30√,𝜋𝑢+_𝑡(𝜋𝑡)²cos(𝜋𝑥)cos(𝜋𝑦),0<𝑡,𝑥<1,𝑦<1,𝑢(0,𝑥,𝑦)=0,0≤𝑥≤1,0≤𝑦≤1,_𝑥𝑢(𝑡,0,𝑦)=𝑢_𝑥(𝑡,1,𝑦)=0,0≤𝑡≤1,0≤𝑦≤1,_𝑦(𝑡,𝑥,0)=𝑢_𝑦(𝑡,𝑥,1)=0,0≤𝑡≤1,0≤𝑥≤1(3.11) for solving the two-dimensional fractional hyperbolic partial differential equation is considered.

According to formula (3.6), the iteration formula for (3.11) is given by 𝑢_𝑛+1(𝑡,𝑥,𝑦)=𝑢+_𝑛(𝑡,𝑥,𝑦)_𝑡0(𝑠−𝑡)𝜕²𝑢_𝑛(𝑠,𝑥,𝑦)𝜕𝑠²+𝐷^𝑠1/2𝑢_𝑛(𝑠,𝑥,𝑦)−𝜕𝑢^2𝑛(𝑠,𝑥,𝑦)𝜕𝑥²−𝜕𝑢^2𝑛(𝑠,𝑥,𝑦)𝜕𝑦²+𝑢_𝑛(𝑠,𝑥,𝑦)−𝑓(𝑠,𝑥,𝑦)𝑑𝑠;(3.12) we start with an initial approximation 𝑢₀(𝑡,𝑥,𝑦)=𝑢(0,𝑥,𝑦)+𝑡𝑢_𝑡(0,𝑥,𝑦).(3.13) Using the above iteration formula (3.12), we can obtain the other components as 𝑢𝑢₀(𝑡,𝑥,𝑦)=0,₁(𝑡,𝑥,𝑦)=1420√𝜋128𝑡^7/2𝑢+35𝑡⁴√𝜋+35𝑡⁴𝜋^5/2+420𝑡²√𝜋cos(𝜋𝑥)cos(𝜋𝑦),₂(𝑡,𝑥,𝑦)=cos(𝜋𝑥)cos(𝜋𝑦)321051𝑡^5/2𝑡^7/2+1121𝑡^5/2+1𝑡⁴√𝜋121𝑡^5/2𝑡⁴𝜋^5/2+1𝑡^5/2×−0.2195931203𝑡_11/2𝑡²√𝜋+cos(𝜋𝑥)cos(𝜋𝑦)√𝜋1𝑡^5/2−0.1719434921√𝜋𝑡^7/2×1𝑡^5/2−0.6261860981𝑡⁶√𝜋1𝑡^5/2−1.728267400√𝜋𝑡⁴1𝑡^5/2−0.01666666667√𝜋𝑡^5/2⋮(3.14) The exact solution of (3.11) is given by as follows: 𝑢(𝑡,𝑥,𝑦)=𝑡²cos(𝜋𝑥)cos(𝜋𝑦).(3.15) The figure of the exact solution of (3.11) is shown by the Figure 5.

The figure of (3.14) is given by the Figure 4, and so on; in the same manner the rest of the components of the iteration formula (3.12) can be obtained using the Maple package.

Figure 4 

Variational iteration method for (3.14).

Figure 5 

The surface shows the exact solution 𝑢(𝑡,𝑥,𝑦) for (3.11).