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

Authors: Allaberen Ashyralyev; Fadime Dal
Journal: Discrete Dynamics in Nature and Society
 (2012)

---

## 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.


---

## Body

## 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
u
(
t
,
x
)
∂
t
2
-
∑
r
=
1
m
‍
(
a
r
(
x
)
u
x
r
)
x
r
+
D
t
1
/
2
u
(
t
,
x
)
+
σ
u
(
t
,
x
)
=
f
(
t
,
x
)
,
x
=
(
x
1
,
…
,
x
m
)
∈
Ω
,
0
<
t
<
1
,
u
(
0
,
x
)
=
0
,
u
t
(
0
,
x
)
=
0
,
x
∈
Ω
¯
;
∂
u
(
t
,
x
)
∂
n
¯
=
0
,
x
∈
S
,
are studied. Here Ω is the unit open cube in the m-dimensional Euclidean space: ℝ
m
:
{
Ω
=
x
=
(
x
1
,
…
,
x
m
)
:
0
<
x
j
<
1
,
1
≤
j
≤
m
} with boundary S
,
Ω
¯
=
Ω
∪
S; a
r
(
x
)
(
x
∈
Ω
) and f
(
t
,
x
)
(
t
∈
(
0,1
)
,
x
∈
Ω
) are given smooth functions and a
r
(
x
)
≥
a
>
0.


### 1.1. Definition

The Caputo fractional derivative of orderα
>
0 of a continuous function u
(
t
,
x
) is defined by
(1.2)
D
a
+
α
u
(
t
,
x
)
=
1
Γ
(
1
-
α
)
∫
a
t
‍
u
′
(
t
,
x
)
(
t
-
s
)
α
d
s
,
where Γ
(
·
) is the gamma function.


## 1.1. Definition

The Caputo fractional derivative of orderα
>
0 of a continuous function u
(
t
,
x
) is defined by
(1.2)
D
a
+
α
u
(
t
,
x
)
=
1
Γ
(
1
-
α
)
∫
a
t
‍
u
′
(
t
,
x
)
(
t
-
s
)
α
d
s
,
where Γ
(
·
) is the gamma function.


## 2. The Finite Difference Method

In this section, we consider the first order of accuracy int 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
(2.1)
Ω
~
h
=
{
x
=
x
r
=
(
h
1
r
1
,
…
,
h
m
r
m
)
,
r
=
(
r
1
,
…
,
r
m
)
,
0
≤
r
j
≤
N
j
,
h
j
N
j
=
1
,
j
=
1
,
…
,
m
}
,
Ω
h
=
Ω
~
h
∩
Ω
,
S
h
=
Ω
~
h
∩
S
.
We introduce the Banach space L
2
h
=
L
2
(
Ω
~
h
) of the grid functions φ
h
(
x
)
=
{
φ
(
h
1
r
1
,
…
,
h
m
r
m
)
} defined on Ω
~
h, equipped with the norm
(2.2)
‖
φ
h
‖
L
2
(
Ω
~
h
)
=
(
∑
x
∈
Ω
h
¯
‍
|
φ
h
(
x
)
|
2
h
1
⋯
h
m
)
1
/
2
.
To the differential operator A
x generated by problem (1.1), we assign the difference operator A
h
x by the formula
(2.3)
A
h
x
u
h
=
-
∑
r
=
1
m
‍
(
a
r
(
x
)
u
x
-
r
h
)
x
r
,
j
+
σ
u
x
h
acting in the space of grid functions u
h
(
x
), satisfying the conditions D
h
u
h
(
x
)
=
0 for all x
∈
S
h. It is known that A
h
x is a self-adjoint positive definite operator in L
2
(
Ω
~
h
). With the help of A
h
x we arrive at the initial boundary value problem
(2.4)
d
2
v
h
(
t
,
x
)
d
t
2
+
D
t
1
/
2
v
h
(
t
,
x
)
+
A
h
x
v
h
(
t
,
x
)
=
f
h
(
t
,
x
)
,
0
≤
t
≤
1
,
x
∈
Ω
h
,
v
h
(
0
,
x
)
=
0
,
d
v
h
(
0
,
x
)
d
t
=
0
,
x
∈
Ω
~
,
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
(2.5)
u
k
+
1
h
(
x
)
-
2
u
k
h
(
x
)
+
u
k
-
1
h
(
x
)
τ
2
+
1
π
∑
m
=
1
k
‍
Γ
(
k
-
m
+
1
/
2
)
(
k
-
m
)
!
u
m
h
-
u
m
-
1
h
τ
1
/
2
+
A
h
x
u
k
+
1
h
=
f
k
h
(
x
)
,
x
∈
Ω
~
h
,
f
k
h
(
x
)
=
f
(
t
k
,
x
)
,
t
k
=
k
τ
,
1
≤
k
≤
N
-
1
,
N
τ
=
1
,
x
∈
Ω
~
h
,
u
1
h
(
x
)
-
u
0
h
(
x
)
τ
=
0
,
u
0
h
(
x
)
=
0
,
x
∈
Ω
~
h
.
Here Γ
(
k
-
m
+
1
/
2
)
=
∫
0
∞
t
k
-
m
-
1
/
2
e
-
t
d
t.Theorem 2.1.
Letτ and |
h
|
=
h
1
2
+
⋯
+
h
m
2 be sufficiently small numbers. Then, the solutions of difference scheme (2.5) satisfy the following stability estimates:
(2.6)
max
⁡
1
≤
k
≤
N
‖
u
k
h
‖
L
2
h
+
max
⁡
1
≤
k
≤
N
‖
u
k
h
-
u
k
-
1
h
τ
‖
L
2
h
≤
C
1
max
⁡
1
≤
k
≤
N
-
1
‖
f
k
h
‖
L
2
h
,
max
⁡
1
≤
k
≤
N
-
1
‖
τ
-
2
(
u
k
+
1
h
-
2
u
k
h
+
u
k
-
1
h
)
‖
L
2
h
+
max
⁡
1
≤
k
≤
N
‖
(
u
k
h
)
x
-
r
x
r
‖
L
2
h
≤
C
2
[
‖
f
1
h
‖
L
2
h
+
max
⁡
2
≤
k
≤
N
-
1
‖
τ
-
1
(
f
k
h
-
f
k
-
1
h
)
‖
L
2
h
]
.
HereC
1 and C
2 do not depend on τ, h, and f
k
h
,
1
≤
k
≤
N
-
1.The proof of Theorem2.1 is based on the self-adjointness and positive definitness of operator A
h
x in L
2
h and on the following theorem on the coercivity inequality for the solution of the elliptic difference problem in L
2
h.Theorem 2.2.
For the solutions of the elliptic difference problem(2.7)
A
h
x
u
h
(
x
)
=
ω
h
(
x
)
,
x
∈
Ω
h
,
D
h
u
h
(
x
)
=
0
,
x
∈
S
h
,
the following coercivity inequality holds [30]:
(2.8)
∑
r
=
1
m
‍
‖
u
x
r
x
-
r
h
‖
L
2
h
≤
C
‖
ω
h
‖
L
2
h
.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.9)
D
t
2
u
(
t
,
x
)
+
D
t
1
/
2
u
(
t
,
x
)
-
u
x
x
(
t
,
x
)
+
u
(
t
,
x
)
=
f
(
t
,
x
)
,
f
(
t
,
x
)
=
(
2
+
t
2
+
8
t
3
/
2
3
π
+
(
π
t
)
2
)
cos
⁡
⁡
(
π
x
)
,
0
<
t
,
x
<
1
,
u
(
0
,
x
)
=
0
,
u
t
(
0
,
x
)
=
0
,
0
≤
x
≤
1
,
u
x
(
t
,
0
)
=
u
x
(
t
,
1
)
=
0
,
0
≤
t
≤
1
for solving the one-dimensional fractional hyperbolic partial differential equation is considered. Applying difference scheme (2.5), we obtained
(2.10)
u
n
k
+
1
-
2
u
n
k
+
u
n
k
-
1
τ
2
+
1
π
∑
m
=
1
k
‍
Γ
(
k
-
m
+
1
/
2
)
(
k
-
m
)
!
(
u
n
m
-
u
n
m
-
1
τ
1
/
2
)
-
(
u
n
+
1
k
+
1
-
2
u
n
k
+
1
+
u
n
-
1
k
+
1
h
2
)
+
u
n
k
=
φ
n
k
,
φ
n
k
=
f
(
t
k
,
x
n
)
,
1
≤
k
≤
N
-
1
,
1
≤
n
≤
M
-
1
,
u
n
0
=
0
,
τ
-
1
(
u
n
1
-
u
n
0
)
=
0
,
0
≤
n
≤
M
,
u
1
k
-
u
0
k
=
u
M
k
-
u
M
-
1
k
=
0
,
0
≤
k
≤
N
.
We get the system of equations in the matrix form:
(2.11)
A
U
n
+
1
+
B
U
n
+
C
U
n
-
1
=
D
φ
n
,
1
≤
n
≤
M
-
1
,
U
1
=
U
0
,
U
M
=
U
M
-
1
,
where
(2.12)
A
=
[
0
0
0
0
⋯
0
0
0
0
0
0
⋯
0
0
0
0
0
0
⋯
0
0
0
0
a
0
⋯
0
0
0
0
0
a
⋯
0
0
⋯
⋯
⋯
⋯
⋯
⋯
⋯
0
0
0
0
⋯
a
0
0
0
0
0
⋯
0
0
]
(
N
+
1
)
×
(
N
+
1
)
,
B
=
[
b
1,1
0
0
0
⋯
0
0
b
2,1
b
2,2
0
0
⋯
0
0
b
3.1
b
3,2
b
3,3
0
⋯
0
0
b
4,1
b
4,2
b
4,3
b
4,4
⋯
0
0
⋯
⋯
⋯
⋯
⋯
⋯
⋯
b
N
,
1
b
N
,
2
b
N
,
3
b
N
,
4
⋯
b
N
,
N
0
b
N
+
1,1
b
N
+
1,2
b
N
+
1,3
b
N
+
1,4
⋯
b
N
+
1
,
N
b
N
+
1
,
N
+
1
]
(
N
+
1
)
×
(
N
+
1
)
,
C
=
A
,
D
=
[
0
0
0
0
⋯
0
0
0
0
0
0
⋯
0
0
0
0
0
0
⋯
0
0
0
0
1
0
⋯
0
0
0
0
0
1
⋯
0
0
⋯
⋯
⋯
⋯
⋯
⋯
⋯
0
0
0
0
⋯
1
0
0
0
0
0
⋯
0
1
]
(
N
+
1
)
×
(
N
+
1
)
,
U
s
=
[
U
s
0
U
s
1
U
s
2
U
s
3
⋯
U
s
N
-
1
U
s
N
]
(
N
+
1
)
×
(
1
)
,
s
=
n
-
1
,
n
,
n
+
1
.
Here
(2.13)
a
=
-
1
h
2
,
b
1,1
=
1
,
b
2,1
=
-
1
,
b
2,2
=
1
,
b
3,1
=
1
τ
2
-
1
τ
1
/
2
,
b
3,2
=
-
2
τ
2
+
1
τ
1
/
2
,
b
3,3
=
1
+
1
τ
2
+
2
h
2
,
b
k
+
2,1
=
-
1
π
Γ
(
k
-
1
+
1
/
2
)
Γ
(
k
)
τ
1
/
2
,
2
≤
k
≤
N
-
1
,
b
k
+
2
,
k
+
1
=
-
2
τ
2
+
1
τ
1
/
2
,
1
≤
k
≤
N
-
1
,
b
k
+
2
,
k
=
1
τ
2
+
1
π
(
Γ
(
1
+
0.5
)
Γ
(
2
)
-
Γ
(
0.5
)
Γ
(
1
)
)
1
τ
1
/
2
,
2
≤
k
≤
N
-
1
,
b
k
+
2
,
k
+
2
=
1
+
1
τ
2
+
2
h
2
,
1
≤
k
≤
N
-
1
,
b
k
+
2
,
i
+
1
=
1
π
(
Γ
(
k
-
i
+
1
/
2
)
Γ
(
k
-
(
i
-
1
)
)
-
Γ
(
k
-
(
i
+
1
)
+
1
/
2
)
Γ
(
k
-
(
i
-
1
)
-
1
)
)
1
τ
1
/
2
,
3
≤
k
≤
N
-
1
,
1
≤
i
≤
k
-
2
,
φ
n
k
=
(
2
+
(
k
τ
)
2
+
8
(
k
τ
)
3
/
2
3
π
+
(
π
k
τ
)
2
)
cos
⁡
⁡
π
(
n
h
)
,
φ
n
=
[
φ
n
0
φ
n
1
φ
n
2
⋯
φ
n
N
]
(
N
+
1
)
×
1
.
So, we have the second-order difference equation with respect to n matrix coefficients. To solve this difference equation, we have applied a procedure of modified Gauss elimination method for difference equation with respect to k matrix coefficients. Hence, we seek a solution of the matrix equation in the following form:
(2.14)
U
j
=
α
j
+
1
U
j
+
1
+
β
j
+
1
,
n
=
M
-
1
,
…
,
2,1,  α
j
(
j
=
1
,
…
,
M
) are (
N
+
1
)
×
(
N
+
1
) square matrices, and β
j
(
j
=
1
,
…
,
M
) are (
N
+
1
)
×
1 column matrices defined by
(2.15)
α
n
+
1
=
(
B
+
C
α
n
)
-
1
(
-
A
)
,
β
n
+
1
=
(
B
+
C
α
n
)
-
1
(
D
φ
n
-
C
β
n
)
,
n
=
2,3
,
…
,
M
,
where
(2.16)
α
1
=
[
1
0
0
⋯
0
0
1
0
⋯
0
0
0
1
⋯
0
⋯
⋯
⋯
⋯
⋯
0
0
0
⋯
1
]
(
N
+
1
)
×
(
N
+
1
)
,
β
1
=
[
0
0
0
⋯
0
]
(
N
+
1
)
×
1
.Now, we will give the results of the numerical analysis. First, we give an estimate for the constantsC
1 and C
2 figuring in the stability estimates of Theorem 2.1. We have
(2.17)
C
1
=
max
⁡
f
,
u
(
C
t
1
)
,
C
2
=
max
⁡
f
,
u
(
C
t
2
)
,
C
t
1
=
[
max
⁡
1
≤
k
≤
N
‖
u
k
h
‖
L
2
h
+
max
⁡
1
≤
k
≤
N
‖
τ
-
1
(
u
k
h
-
u
k
-
1
h
)
‖
L
2
h
]
×
(
max
⁡
1
≤
k
≤
N
-
1
‖
f
k
h
‖
L
2
h
)
-
1
,
C
t
2
=
[
max
⁡
1
≤
k
≤
N
-
1
‖
τ
-
2
(
u
k
+
1
h
-
2
u
k
h
+
u
k
-
1
h
)
‖
L
2
h
+
max
⁡
1
≤
k
≤
N
∑
r
=
1
n
‍
‖
(
u
k
h
)
x
¯
r
,
x
r
‖
L
2
h
]
×
(
max
⁡
2
≤
k
≤
N
-
1
‖
τ
-
1
(
f
k
h
-
f
k
-
1
h
)
‖
L
2
h
+
‖
f
1
h
‖
L
2
h
)
-
1
.
The constants C
t
1 and C
t
2 in the case of numerical solution of initial-boundary value problem (2.9) are computed. The constants C
t
1 and C
t
2 are given in Table 1 for N
=
20
,
40,80, and M
=
80, respectively.Table 1
Stability estimates for (2.9).
M
=
80
M
=
80
M
=
80
N
=
20
N
=
40
N
=
80
The values ofC
t
1
0.2096
0.2073
0.2061
The values ofC
t
2
0.2075
0.1223
0.0670Second, for the accurate comparison of the difference scheme considered, the errors computed by(2.18)
E
0
=
max
⁡
1
≤
k
≤
N
-
1
(
∑
n
=
1
M
-
1
‍
|
u
(
t
k
,
x
n
)
-
u
n
k
|
2
h
)
1
/
2
,
E
1
=
max
⁡
1
≤
k
≤
N
-
1
(
∑
n
=
1
M
-
1
‍
|
u
t
(
t
k
,
x
n
)
-
(
u
n
k
+
1
-
u
n
k
-
1
)
2
τ
|
2
h
)
1
/
2
,
E
2
=
max
⁡
1
≤
k
≤
N
-
1
(
∑
n
=
1
M
-
1
‍
|
u
t
t
(
t
k
,
x
n
)
-
(
u
n
k
+
1
-
2
u
n
k
+
u
n
k
-
1
)
τ
2
|
2
h
)
1
/
2
of the numerical solutions are recorded for higher values of N
=
M, where u
(
t
k
,
x
n
) represents the exact solution and u
n
k represents the numerical solution at (
t
k
,
x
n
). The errors E
0
,
E
1 and E
2 results are shown in Table 2 for N
=
20,40,60 and M
=
60, respectively.Table 2
Comparison of the errors for the difference scheme.
Method
M = 80
M = 80
M = 60
N = 20
N = 40
N = 60
Comparison of errors (E
0) for approximate solutions
0.0071
0.0037
0.0008
Comparison of errors (E
1) for approximate solutions
0.1030
0.0521
0.0491
Comparison of errors (E
2) for approximate solutions
0.1224
0.0806
0.0882The 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:
(2.19)
u
(
t
,
x
)
=
t
2
cos
⁡
⁡
(
π
x
)
.
The  figure  of  the  exact   solution of  (2.9)  is  shown  by  the  Figure 1.Figure 1
The surface shows the exact solutionu
(
t
,
x
) for (2.9).Figure 2
Difference scheme solution for (2.9).


## 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
(2.1)
Ω
~
h
=
{
x
=
x
r
=
(
h
1
r
1
,
…
,
h
m
r
m
)
,
r
=
(
r
1
,
…
,
r
m
)
,
0
≤
r
j
≤
N
j
,
h
j
N
j
=
1
,
j
=
1
,
…
,
m
}
,
Ω
h
=
Ω
~
h
∩
Ω
,
S
h
=
Ω
~
h
∩
S
.
We introduce the Banach space L
2
h
=
L
2
(
Ω
~
h
) of the grid functions φ
h
(
x
)
=
{
φ
(
h
1
r
1
,
…
,
h
m
r
m
)
} defined on Ω
~
h, equipped with the norm
(2.2)
‖
φ
h
‖
L
2
(
Ω
~
h
)
=
(
∑
x
∈
Ω
h
¯
‍
|
φ
h
(
x
)
|
2
h
1
⋯
h
m
)
1
/
2
.
To the differential operator A
x generated by problem (1.1), we assign the difference operator A
h
x by the formula
(2.3)
A
h
x
u
h
=
-
∑
r
=
1
m
‍
(
a
r
(
x
)
u
x
-
r
h
)
x
r
,
j
+
σ
u
x
h
acting in the space of grid functions u
h
(
x
), satisfying the conditions D
h
u
h
(
x
)
=
0 for all x
∈
S
h. It is known that A
h
x is a self-adjoint positive definite operator in L
2
(
Ω
~
h
). With the help of A
h
x we arrive at the initial boundary value problem
(2.4)
d
2
v
h
(
t
,
x
)
d
t
2
+
D
t
1
/
2
v
h
(
t
,
x
)
+
A
h
x
v
h
(
t
,
x
)
=
f
h
(
t
,
x
)
,
0
≤
t
≤
1
,
x
∈
Ω
h
,
v
h
(
0
,
x
)
=
0
,
d
v
h
(
0
,
x
)
d
t
=
0
,
x
∈
Ω
~
,
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
(2.5)
u
k
+
1
h
(
x
)
-
2
u
k
h
(
x
)
+
u
k
-
1
h
(
x
)
τ
2
+
1
π
∑
m
=
1
k
‍
Γ
(
k
-
m
+
1
/
2
)
(
k
-
m
)
!
u
m
h
-
u
m
-
1
h
τ
1
/
2
+
A
h
x
u
k
+
1
h
=
f
k
h
(
x
)
,
x
∈
Ω
~
h
,
f
k
h
(
x
)
=
f
(
t
k
,
x
)
,
t
k
=
k
τ
,
1
≤
k
≤
N
-
1
,
N
τ
=
1
,
x
∈
Ω
~
h
,
u
1
h
(
x
)
-
u
0
h
(
x
)
τ
=
0
,
u
0
h
(
x
)
=
0
,
x
∈
Ω
~
h
.
Here Γ
(
k
-
m
+
1
/
2
)
=
∫
0
∞
t
k
-
m
-
1
/
2
e
-
t
d
t.Theorem 2.1.
Letτ and |
h
|
=
h
1
2
+
⋯
+
h
m
2 be sufficiently small numbers. Then, the solutions of difference scheme (2.5) satisfy the following stability estimates:
(2.6)
max
⁡
1
≤
k
≤
N
‖
u
k
h
‖
L
2
h
+
max
⁡
1
≤
k
≤
N
‖
u
k
h
-
u
k
-
1
h
τ
‖
L
2
h
≤
C
1
max
⁡
1
≤
k
≤
N
-
1
‖
f
k
h
‖
L
2
h
,
max
⁡
1
≤
k
≤
N
-
1
‖
τ
-
2
(
u
k
+
1
h
-
2
u
k
h
+
u
k
-
1
h
)
‖
L
2
h
+
max
⁡
1
≤
k
≤
N
‖
(
u
k
h
)
x
-
r
x
r
‖
L
2
h
≤
C
2
[
‖
f
1
h
‖
L
2
h
+
max
⁡
2
≤
k
≤
N
-
1
‖
τ
-
1
(
f
k
h
-
f
k
-
1
h
)
‖
L
2
h
]
.
HereC
1 and C
2 do not depend on τ, h, and f
k
h
,
1
≤
k
≤
N
-
1.The proof of Theorem2.1 is based on the self-adjointness and positive definitness of operator A
h
x in L
2
h and on the following theorem on the coercivity inequality for the solution of the elliptic difference problem in L
2
h.Theorem 2.2.
For the solutions of the elliptic difference problem(2.7)
A
h
x
u
h
(
x
)
=
ω
h
(
x
)
,
x
∈
Ω
h
,
D
h
u
h
(
x
)
=
0
,
x
∈
S
h
,
the following coercivity inequality holds [30]:
(2.8)
∑
r
=
1
m
‍
‖
u
x
r
x
-
r
h
‖
L
2
h
≤
C
‖
ω
h
‖
L
2
h
.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.9)
D
t
2
u
(
t
,
x
)
+
D
t
1
/
2
u
(
t
,
x
)
-
u
x
x
(
t
,
x
)
+
u
(
t
,
x
)
=
f
(
t
,
x
)
,
f
(
t
,
x
)
=
(
2
+
t
2
+
8
t
3
/
2
3
π
+
(
π
t
)
2
)
cos
⁡
⁡
(
π
x
)
,
0
<
t
,
x
<
1
,
u
(
0
,
x
)
=
0
,
u
t
(
0
,
x
)
=
0
,
0
≤
x
≤
1
,
u
x
(
t
,
0
)
=
u
x
(
t
,
1
)
=
0
,
0
≤
t
≤
1
for solving the one-dimensional fractional hyperbolic partial differential equation is considered. Applying difference scheme (2.5), we obtained
(2.10)
u
n
k
+
1
-
2
u
n
k
+
u
n
k
-
1
τ
2
+
1
π
∑
m
=
1
k
‍
Γ
(
k
-
m
+
1
/
2
)
(
k
-
m
)
!
(
u
n
m
-
u
n
m
-
1
τ
1
/
2
)
-
(
u
n
+
1
k
+
1
-
2
u
n
k
+
1
+
u
n
-
1
k
+
1
h
2
)
+
u
n
k
=
φ
n
k
,
φ
n
k
=
f
(
t
k
,
x
n
)
,
1
≤
k
≤
N
-
1
,
1
≤
n
≤
M
-
1
,
u
n
0
=
0
,
τ
-
1
(
u
n
1
-
u
n
0
)
=
0
,
0
≤
n
≤
M
,
u
1
k
-
u
0
k
=
u
M
k
-
u
M
-
1
k
=
0
,
0
≤
k
≤
N
.
We get the system of equations in the matrix form:
(2.11)
A
U
n
+
1
+
B
U
n
+
C
U
n
-
1
=
D
φ
n
,
1
≤
n
≤
M
-
1
,
U
1
=
U
0
,
U
M
=
U
M
-
1
,
where
(2.12)
A
=
[
0
0
0
0
⋯
0
0
0
0
0
0
⋯
0
0
0
0
0
0
⋯
0
0
0
0
a
0
⋯
0
0
0
0
0
a
⋯
0
0
⋯
⋯
⋯
⋯
⋯
⋯
⋯
0
0
0
0
⋯
a
0
0
0
0
0
⋯
0
0
]
(
N
+
1
)
×
(
N
+
1
)
,
B
=
[
b
1,1
0
0
0
⋯
0
0
b
2,1
b
2,2
0
0
⋯
0
0
b
3.1
b
3,2
b
3,3
0
⋯
0
0
b
4,1
b
4,2
b
4,3
b
4,4
⋯
0
0
⋯
⋯
⋯
⋯
⋯
⋯
⋯
b
N
,
1
b
N
,
2
b
N
,
3
b
N
,
4
⋯
b
N
,
N
0
b
N
+
1,1
b
N
+
1,2
b
N
+
1,3
b
N
+
1,4
⋯
b
N
+
1
,
N
b
N
+
1
,
N
+
1
]
(
N
+
1
)
×
(
N
+
1
)
,
C
=
A
,
D
=
[
0
0
0
0
⋯
0
0
0
0
0
0
⋯
0
0
0
0
0
0
⋯
0
0
0
0
1
0
⋯
0
0
0
0
0
1
⋯
0
0
⋯
⋯
⋯
⋯
⋯
⋯
⋯
0
0
0
0
⋯
1
0
0
0
0
0
⋯
0
1
]
(
N
+
1
)
×
(
N
+
1
)
,
U
s
=
[
U
s
0
U
s
1
U
s
2
U
s
3
⋯
U
s
N
-
1
U
s
N
]
(
N
+
1
)
×
(
1
)
,
s
=
n
-
1
,
n
,
n
+
1
.
Here
(2.13)
a
=
-
1
h
2
,
b
1,1
=
1
,
b
2,1
=
-
1
,
b
2,2
=
1
,
b
3,1
=
1
τ
2
-
1
τ
1
/
2
,
b
3,2
=
-
2
τ
2
+
1
τ
1
/
2
,
b
3,3
=
1
+
1
τ
2
+
2
h
2
,
b
k
+
2,1
=
-
1
π
Γ
(
k
-
1
+
1
/
2
)
Γ
(
k
)
τ
1
/
2
,
2
≤
k
≤
N
-
1
,
b
k
+
2
,
k
+
1
=
-
2
τ
2
+
1
τ
1
/
2
,
1
≤
k
≤
N
-
1
,
b
k
+
2
,
k
=
1
τ
2
+
1
π
(
Γ
(
1
+
0.5
)
Γ
(
2
)
-
Γ
(
0.5
)
Γ
(
1
)
)
1
τ
1
/
2
,
2
≤
k
≤
N
-
1
,
b
k
+
2
,
k
+
2
=
1
+
1
τ
2
+
2
h
2
,
1
≤
k
≤
N
-
1
,
b
k
+
2
,
i
+
1
=
1
π
(
Γ
(
k
-
i
+
1
/
2
)
Γ
(
k
-
(
i
-
1
)
)
-
Γ
(
k
-
(
i
+
1
)
+
1
/
2
)
Γ
(
k
-
(
i
-
1
)
-
1
)
)
1
τ
1
/
2
,
3
≤
k
≤
N
-
1
,
1
≤
i
≤
k
-
2
,
φ
n
k
=
(
2
+
(
k
τ
)
2
+
8
(
k
τ
)
3
/
2
3
π
+
(
π
k
τ
)
2
)
cos
⁡
⁡
π
(
n
h
)
,
φ
n
=
[
φ
n
0
φ
n
1
φ
n
2
⋯
φ
n
N
]
(
N
+
1
)
×
1
.
So, we have the second-order difference equation with respect to n matrix coefficients. To solve this difference equation, we have applied a procedure of modified Gauss elimination method for difference equation with respect to k matrix coefficients. Hence, we seek a solution of the matrix equation in the following form:
(2.14)
U
j
=
α
j
+
1
U
j
+
1
+
β
j
+
1
,
n
=
M
-
1
,
…
,
2,1,  α
j
(
j
=
1
,
…
,
M
) are (
N
+
1
)
×
(
N
+
1
) square matrices, and β
j
(
j
=
1
,
…
,
M
) are (
N
+
1
)
×
1 column matrices defined by
(2.15)
α
n
+
1
=
(
B
+
C
α
n
)
-
1
(
-
A
)
,
β
n
+
1
=
(
B
+
C
α
n
)
-
1
(
D
φ
n
-
C
β
n
)
,
n
=
2,3
,
…
,
M
,
where
(2.16)
α
1
=
[
1
0
0
⋯
0
0
1
0
⋯
0
0
0
1
⋯
0
⋯
⋯
⋯
⋯
⋯
0
0
0
⋯
1
]
(
N
+
1
)
×
(
N
+
1
)
,
β
1
=
[
0
0
0
⋯
0
]
(
N
+
1
)
×
1
.Now, we will give the results of the numerical analysis. First, we give an estimate for the constantsC
1 and C
2 figuring in the stability estimates of Theorem 2.1. We have
(2.17)
C
1
=
max
⁡
f
,
u
(
C
t
1
)
,
C
2
=
max
⁡
f
,
u
(
C
t
2
)
,
C
t
1
=
[
max
⁡
1
≤
k
≤
N
‖
u
k
h
‖
L
2
h
+
max
⁡
1
≤
k
≤
N
‖
τ
-
1
(
u
k
h
-
u
k
-
1
h
)
‖
L
2
h
]
×
(
max
⁡
1
≤
k
≤
N
-
1
‖
f
k
h
‖
L
2
h
)
-
1
,
C
t
2
=
[
max
⁡
1
≤
k
≤
N
-
1
‖
τ
-
2
(
u
k
+
1
h
-
2
u
k
h
+
u
k
-
1
h
)
‖
L
2
h
+
max
⁡
1
≤
k
≤
N
∑
r
=
1
n
‍
‖
(
u
k
h
)
x
¯
r
,
x
r
‖
L
2
h
]
×
(
max
⁡
2
≤
k
≤
N
-
1
‖
τ
-
1
(
f
k
h
-
f
k
-
1
h
)
‖
L
2
h
+
‖
f
1
h
‖
L
2
h
)
-
1
.
The constants C
t
1 and C
t
2 in the case of numerical solution of initial-boundary value problem (2.9) are computed. The constants C
t
1 and C
t
2 are given in Table 1 for N
=
20
,
40,80, and M
=
80, respectively.Table 1
Stability estimates for (2.9).
M
=
80
M
=
80
M
=
80
N
=
20
N
=
40
N
=
80
The values ofC
t
1
0.2096
0.2073
0.2061
The values ofC
t
2
0.2075
0.1223
0.0670Second, for the accurate comparison of the difference scheme considered, the errors computed by(2.18)
E
0
=
max
⁡
1
≤
k
≤
N
-
1
(
∑
n
=
1
M
-
1
‍
|
u
(
t
k
,
x
n
)
-
u
n
k
|
2
h
)
1
/
2
,
E
1
=
max
⁡
1
≤
k
≤
N
-
1
(
∑
n
=
1
M
-
1
‍
|
u
t
(
t
k
,
x
n
)
-
(
u
n
k
+
1
-
u
n
k
-
1
)
2
τ
|
2
h
)
1
/
2
,
E
2
=
max
⁡
1
≤
k
≤
N
-
1
(
∑
n
=
1
M
-
1
‍
|
u
t
t
(
t
k
,
x
n
)
-
(
u
n
k
+
1
-
2
u
n
k
+
u
n
k
-
1
)
τ
2
|
2
h
)
1
/
2
of the numerical solutions are recorded for higher values of N
=
M, where u
(
t
k
,
x
n
) represents the exact solution and u
n
k represents the numerical solution at (
t
k
,
x
n
). The errors E
0
,
E
1 and E
2 results are shown in Table 2 for N
=
20,40,60 and M
=
60, respectively.Table 2
Comparison of the errors for the difference scheme.
Method
M = 80
M = 80
M = 60
N = 20
N = 40
N = 60
Comparison of errors (E
0) for approximate solutions
0.0071
0.0037
0.0008
Comparison of errors (E
1) for approximate solutions
0.1030
0.0521
0.0491
Comparison of errors (E
2) for approximate solutions
0.1224
0.0806
0.0882The 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:
(2.19)
u
(
t
,
x
)
=
t
2
cos
⁡
⁡
(
π
x
)
.
The  figure  of  the  exact   solution of  (2.9)  is  shown  by  the  Figure 1.Figure 1
The surface shows the exact solutionu
(
t
,
x
) 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:
(3.1)
u
n
+
1
(
t
,
x
)
=
u
n
(
t
,
x
)
+
∫
0
t
‍
λ
[
∂
2
u
(
s
,
x
)
∂
s
2
-
∑
r
=
1
m
‍
(
a
r
(
x
)
u
~
x
r
)
x
r
+
D
s
1
/
2
u
~
(
s
,
x
)
+
σ
u
~
(
s
,
x
)
-
f
(
s
,
x
)
]
d
s
,
where λ is a general Lagrangian multiplier (see, e.g., [31]) and u
~ is considered as a restricted variation as a restricted variation (see, e.g., [32]); that is, δ
u
~
=
0, u
0
(
t
,
x
) is its initial approximation. Using the above correction functional stationary and noticing that δ
u
~
=
0, we obtain
(3.2)
δ
u
n
+
1
(
t
,
x
)
=
δ
u
n
(
t
,
x
)
+
∫
0
t
‍
δ
λ
[
∂
u
n
2
(
t
,
x
)
∂
s
2
]
d
s
,
δ
u
n
+
1
(
t
,
x
)
=
δ
u
n
(
t
,
x
)
-
∂
λ
∂
s
δ
u
n
(
s
,
x
)
|
s
=
t
+
λ
∂
∂
s
(
δ
u
n
(
s
,
x
)
)
|
s
=
t
+
∫
0
t
‍
∂
2
λ
(
t
,
s
)
∂
s
2
δ
u
n
(
s
,
x
)
d
s
=
0
.
From the above relation for any δ
u
n, we get the Euler-Lagrange equation:
(3.3)
∂
λ
2
(
t
,
s
)
∂
s
2
=
0
,
with the following natural boundary conditions:
(3.4)
1
-
∂
λ
(
t
,
s
)
∂
s
|
s
=
t
=
0
,
λ
(
t
,
s
)
|
s
=
t
=
0
.
Therefore, the Lagrange multiplier can be identified as follows:
(3.5)
λ
(
t
,
s
)
=
s
-
t
.
Substituting the identified Lagrange multiplier into (3.1), the following variational iteration formula can be obtained:
(3.6)
u
n
+
1
(
t
,
x
)
=
u
n
(
t
,
x
)
+
∫
0
t
‍
(
s
-
t
)
[
∂
2
u
n
(
s
,
x
)
∂
s
2
-
∑
r
=
1
m
‍
(
a
r
(
x
)
u
n
x
r
)
x
r
+
D
s
1
/
2
u
n
(
s
,
x
)
+
u
n
(
s
,
x
)
-
f
(
s
,
x
)
]
d
s
.In this case, let an initial approximationu
0
(
t
,
x
)
=
u
(
0
,
x
)
+
t
u
t
(
0
,
x
). Then approximate solution takes the form u
(
t
,
x
)
=
lim
n
→
∞
u
n
(
t
,
x
).


### 3.1. Variational Iteration Solution 1

For the numerical result, the mixed problem(3.7)
D
t
2
u
(
t
,
x
)
+
D
t
1
/
2
u
(
t
,
x
)
-
u
x
x
(
t
,
x
)
+
u
(
t
,
x
)
=
f
(
t
,
x
)
,
f
(
t
,
x
)
=
(
2
+
t
2
+
8
t
3
/
2
3
π
+
(
π
t
)
2
)
cos
⁡
⁡
(
π
x
)
,
0
<
t
,
x
<
1
,
u
(
0
,
x
)
=
0
,
u
t
(
0
,
x
)
=
0
,
0
≤
x
≤
1
,
u
x
(
t
,
0
)
=
u
x
(
t
,
1
)
=
0
,
0
≤
t
≤
1
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
(3.8)
u
n
+
1
(
t
,
x
)
=
u
n
(
t
,
x
)
+
∫
0
t
‍
(
s
-
t
)
[
∂
u
n
2
(
s
,
x
)
∂
s
2
+
D
s
1
/
2
u
n
(
s
,
x
)
-
∂
u
n
2
(
s
,
x
)
∂
x
2
+
u
n
(
s
,
x
)
-
f
(
s
,
x
)
]
d
s
.
Now we start with an initial approximation
(3.9)
u
0
(
t
,
x
)
=
u
(
0
,
x
)
+
t
u
t
(
0
,
x
)
.
Using the above iteration formula (3.8), we can obtain the other components as
(3.10)
u
0
(
t
,
x
)
=
0
,
u
1
(
t
,
x
)
=
1
420
π
(
128
t
7
/
2
+
35
t
4
π
+
35
t
4
π
5
/
2
+
420
t
2
π
)
cos
⁡
⁡
(
π
x
)
,
u
2
(
t
,
x
)
=
1
420
π
cos
⁡
(
π
x
)
(
128
t
7
/
2
+
35
t
4
π
+
35
t
4
π
5
/
2
+
420
t
2
π
)
+
cos
⁡
⁡
(
π
x
)
[
-
0.9058003666
t
4
-
0.1510268880
t
11
/
2
-
0.1719434921
t
7
/
2
-
0.3281897218
t
6
-
0.01666666667
t
5
]
⋮
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(3.11)
D
t
2
u
(
t
,
x
,
y
)
+
D
t
1
/
2
u
(
t
,
x
,
y
)
-
u
x
x
(
t
,
x
,
y
)
-
u
y
y
(
t
,
x
,
y
)
+
u
(
t
,
x
,
y
)
=
f
(
t
,
x
,
y
)
,
f
(
t
,
x
,
y
)
=
(
2
+
t
2
+
8
t
3
/
2
3
π
+
(
π
t
)
2
)
cos
⁡
⁡
(
π
x
)
cos
⁡
⁡
(
π
y
)
,
0
<
t
,
x
<
1
,
y
<
1
,
u
(
0
,
x
,
y
)
=
0
,
u
t
(
0
,
x
,
y
)
=
0
,
0
≤
x
≤
1
,
0
≤
y
≤
1
,
u
x
(
t
,
0
,
y
)
=
u
x
(
t
,
1
,
y
)
=
0
,
0
≤
t
≤
1
,
0
≤
y
≤
1
,
u
y
(
t
,
x
,
0
)
=
u
y
(
t
,
x
,
1
)
=
0
,
0
≤
t
≤
1
,
0
≤
x
≤
1
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
(3.12)
u
n
+
1
(
t
,
x
,
y
)
=
u
n
(
t
,
x
,
y
)
+
∫
0
t
‍
(
s
-
t
)
[
∂
2
u
n
(
s
,
x
,
y
)
∂
s
2
+
D
s
1
/
2
u
n
(
s
,
x
,
y
)
-
∂
u
n
2
(
s
,
x
,
y
)
∂
x
2
-
∂
u
n
2
(
s
,
x
,
y
)
∂
y
2
+
u
n
(
s
,
x
,
y
)
-
f
(
s
,
x
,
y
)
]
d
s
;
we start with an initial approximation
(3.13)
u
0
(
t
,
x
,
y
)
=
u
(
0
,
x
,
y
)
+
t
u
t
(
0
,
x
,
y
)
.
Using the above iteration formula (3.12), we can obtain the other components as
(3.14)
u
0
(
t
,
x
,
y
)
=
0
,
u
1
(
t
,
x
,
y
)
=
1
420
π
(
128
t
7
/
2
+
35
t
4
π
+
35
t
4
π
5
/
2
+
420
t
2
π
)
cos
⁡
⁡
(
π
x
)
cos
⁡
⁡
(
π
y
)
,
u
2
(
t
,
x
,
y
)
=
cos
⁡
⁡
(
π
x
)
cos
⁡
⁡
(
π
y
)
[
32
105
(
1
t
)
5
/
2
t
7
/
2
+
(
1
12
)
(
1
t
)
5
/
2
t
4
π
+
1
12
(
1
t
)
5
/
2
t
4
π
5
/
2
+
(
1
t
)
5
/
2
t
2
π
(
1
t
)
5
/
2
]
+
cos
⁡
⁡
(
π
x
)
cos
⁡
⁡
(
π
y
)
×
[
-
0.2195931203
t
11
/
2
π
(
1
t
)
5
/
2
-
0.1719434921
π
t
7
/
2
×
(
1
t
)
5
/
2
-
0.6261860981
t
6
π
(
1
t
)
5
/
2
-
1.728267400
π
t
4
(
1
t
)
5
/
2
-
0.01666666667
π
t
5
/
2
]
⋮
The exact solution of (3.11)  is   given  by  as  follows:
(3.15)
u
(
t
,
x
,
y
)
=
t
2
cos
⁡
⁡
(
π
x
)
cos
⁡
⁡
(
π
y
)
.
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 solutionu
(
t
,
x
,
y
) for (3.11).


## 3.1. Variational Iteration Solution 1

For the numerical result, the mixed problem(3.7)
D
t
2
u
(
t
,
x
)
+
D
t
1
/
2
u
(
t
,
x
)
-
u
x
x
(
t
,
x
)
+
u
(
t
,
x
)
=
f
(
t
,
x
)
,
f
(
t
,
x
)
=
(
2
+
t
2
+
8
t
3
/
2
3
π
+
(
π
t
)
2
)
cos
⁡
⁡
(
π
x
)
,
0
<
t
,
x
<
1
,
u
(
0
,
x
)
=
0
,
u
t
(
0
,
x
)
=
0
,
0
≤
x
≤
1
,
u
x
(
t
,
0
)
=
u
x
(
t
,
1
)
=
0
,
0
≤
t
≤
1
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
(3.8)
u
n
+
1
(
t
,
x
)
=
u
n
(
t
,
x
)
+
∫
0
t
‍
(
s
-
t
)
[
∂
u
n
2
(
s
,
x
)
∂
s
2
+
D
s
1
/
2
u
n
(
s
,
x
)
-
∂
u
n
2
(
s
,
x
)
∂
x
2
+
u
n
(
s
,
x
)
-
f
(
s
,
x
)
]
d
s
.
Now we start with an initial approximation
(3.9)
u
0
(
t
,
x
)
=
u
(
0
,
x
)
+
t
u
t
(
0
,
x
)
.
Using the above iteration formula (3.8), we can obtain the other components as
(3.10)
u
0
(
t
,
x
)
=
0
,
u
1
(
t
,
x
)
=
1
420
π
(
128
t
7
/
2
+
35
t
4
π
+
35
t
4
π
5
/
2
+
420
t
2
π
)
cos
⁡
⁡
(
π
x
)
,
u
2
(
t
,
x
)
=
1
420
π
cos
⁡
(
π
x
)
(
128
t
7
/
2
+
35
t
4
π
+
35
t
4
π
5
/
2
+
420
t
2
π
)
+
cos
⁡
⁡
(
π
x
)
[
-
0.9058003666
t
4
-
0.1510268880
t
11
/
2
-
0.1719434921
t
7
/
2
-
0.3281897218
t
6
-
0.01666666667
t
5
]
⋮
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(3.11)
D
t
2
u
(
t
,
x
,
y
)
+
D
t
1
/
2
u
(
t
,
x
,
y
)
-
u
x
x
(
t
,
x
,
y
)
-
u
y
y
(
t
,
x
,
y
)
+
u
(
t
,
x
,
y
)
=
f
(
t
,
x
,
y
)
,
f
(
t
,
x
,
y
)
=
(
2
+
t
2
+
8
t
3
/
2
3
π
+
(
π
t
)
2
)
cos
⁡
⁡
(
π
x
)
cos
⁡
⁡
(
π
y
)
,
0
<
t
,
x
<
1
,
y
<
1
,
u
(
0
,
x
,
y
)
=
0
,
u
t
(
0
,
x
,
y
)
=
0
,
0
≤
x
≤
1
,
0
≤
y
≤
1
,
u
x
(
t
,
0
,
y
)
=
u
x
(
t
,
1
,
y
)
=
0
,
0
≤
t
≤
1
,
0
≤
y
≤
1
,
u
y
(
t
,
x
,
0
)
=
u
y
(
t
,
x
,
1
)
=
0
,
0
≤
t
≤
1
,
0
≤
x
≤
1
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
(3.12)
u
n
+
1
(
t
,
x
,
y
)
=
u
n
(
t
,
x
,
y
)
+
∫
0
t
‍
(
s
-
t
)
[
∂
2
u
n
(
s
,
x
,
y
)
∂
s
2
+
D
s
1
/
2
u
n
(
s
,
x
,
y
)
-
∂
u
n
2
(
s
,
x
,
y
)
∂
x
2
-
∂
u
n
2
(
s
,
x
,
y
)
∂
y
2
+
u
n
(
s
,
x
,
y
)
-
f
(
s
,
x
,
y
)
]
d
s
;
we start with an initial approximation
(3.13)
u
0
(
t
,
x
,
y
)
=
u
(
0
,
x
,
y
)
+
t
u
t
(
0
,
x
,
y
)
.
Using the above iteration formula (3.12), we can obtain the other components as
(3.14)
u
0
(
t
,
x
,
y
)
=
0
,
u
1
(
t
,
x
,
y
)
=
1
420
π
(
128
t
7
/
2
+
35
t
4
π
+
35
t
4
π
5
/
2
+
420
t
2
π
)
cos
⁡
⁡
(
π
x
)
cos
⁡
⁡
(
π
y
)
,
u
2
(
t
,
x
,
y
)
=
cos
⁡
⁡
(
π
x
)
cos
⁡
⁡
(
π
y
)
[
32
105
(
1
t
)
5
/
2
t
7
/
2
+
(
1
12
)
(
1
t
)
5
/
2
t
4
π
+
1
12
(
1
t
)
5
/
2
t
4
π
5
/
2
+
(
1
t
)
5
/
2
t
2
π
(
1
t
)
5
/
2
]
+
cos
⁡
⁡
(
π
x
)
cos
⁡
⁡
(
π
y
)
×
[
-
0.2195931203
t
11
/
2
π
(
1
t
)
5
/
2
-
0.1719434921
π
t
7
/
2
×
(
1
t
)
5
/
2
-
0.6261860981
t
6
π
(
1
t
)
5
/
2
-
1.728267400
π
t
4
(
1
t
)
5
/
2
-
0.01666666667
π
t
5
/
2
]
⋮
The exact solution of (3.11)  is   given  by  as  follows:
(3.15)
u
(
t
,
x
,
y
)
=
t
2
cos
⁡
⁡
(
π
x
)
cos
⁡
⁡
(
π
y
)
.
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 solutionu
(
t
,
x
,
y
) for (3.11).



---

*Source: 434976-2012-06-28.xml*