A modified homotopy analysis method for solution of fractional wave equations

Abstract

The aim of this article is to introduce a modified analytical approach to obtain quick and accurate solution of wave-like fractional physical models. This modified analytical approach is an innovative adjustment in Laplace transform algorithm and homotopy analysis method for fractional partial differential equations. The proposed technique solves the problems using Adomian’s polynomials. The homotopy analysis transform method utilizes a simple and powerful method to adjust and control the convergence region of the infinite series solution using an auxiliary parameter. The numerical solutions obtained by this modified proposed method indicate that the approach is easy to implement, highly accurate, and computationally very attractive. A good agreement between the obtained solutions and some well-known results has been obtained.

Keywords

Fractional wave equation fractional derivatives analytical solution Mittag-Leffler function Laplace transform homotopy analysis transform method

Introduction

Fractional differential equations have drawn the interest of many researchers^1–4 due to their important applications in science and engineering. Fractional derivatives provide an excellent instrument for the description of memory and hereditary properties of various materials and processes.

Our concern in this work is to consider the numerical solution of the time-fractional wave equations. Wave-like model can describe many physical problems in different fields of science and engineering. These physical problems describe some nonlinear phenomena such as earthquake stresses⁵ and non-homogeneous elastic waves in soils.⁶ Recently, there has been a growing interest for obtaining the explicit solutions to heat-like and wave-like models by analytic techniques. The different types of heat-like and wave-like physical models in physics have been solved by Wazwaz and Gorguis,⁷ Shou and He,⁸ and Özis and Agırseven.⁹

In this article, the homotopy analysis transform method (HATM) basically illustrates how the Laplace transform can be used to approximate analytical solutions of the wave-like fractional models by manipulating the homotopy analysis method (HAM). The proposed method is coupling of the HAM and Laplace transform. The main advantage of this proposed method is its capability of combining two powerful methods for obtaining rapid convergent series for fractional partial differential equations. HAM was first proposed and applied by Liao^10–13 based on homotopy, a fundamental concept in topology and differential geometry. The HAM has been successfully applied by many researchers for solving linear and nonlinear partial differential equations.^14–22 In recent years, many researchers have paid attention to obtaining solutions to linear and nonlinear differential and integral equations using various methods by combining the Laplace transform method. Among these, we may mention the following: Laplace decomposition methods^23,24 and homotopy perturbation transform method.^25–28 Recently, Khan et al.²⁹ have applied to obtain the solutions of the Blasius flow equation on a semi-infinite domain by coupling homotopy analysis with the Laplace transform method.

The main aim of this article presents approximate analytical solutions of wave-like physical models with time-fractional derivative $α$ ( $1 < α < 2$ ) in the form of a rapidly convergent series with easily computable components using HATM.

Definition 1

The Laplace transform $L [u (x, t)]$ of the Riemann–Liouville fractional integral $I_{t}^{α}$ is defined as³

L [I_{t}^{α} u (x, t)] = s^{- α} L [u (x, t)]

(1)

Definition 2

The Laplace transform $L [u (x, t)]$ of Caputo’s fractional derivative $D_{t}^{α}$ is defined as³

\begin{matrix} L [D_{t}^{α} u (x, t)] = s^{α} L [u (x, t)] \\ - \sum_{k = 0}^{n - 1} s^{(α - k - 1)} u^{(k)} (x, 0), n - 1 < α \leq n \end{matrix}

(2)

Basic idea of HATM

To illustrate the basic idea of the HATM for the fractional partial differential equation, we consider the following fractional partial differential equation as

D_{t}^{n α} u (r, t) + R [r] u (r, t) + N [r] u (r, t) = g (r, t), t > 0, r \in R^{3}, n - 1 < n α \leq n

(3)

where $D_{t}^{n α} = (\partial^{n α} / \partial t^{n α}),$ $R [r]$ is the linear operator in $r \in R^{3},$ $N [r]$ is the general nonlinear operator in $r \in R^{3}$ , and $g (r, t)$ are continuous functions. For simplicity, we ignore all initial and boundary conditions, which can be treated in a similar way. Now, the methodology consists of applying Laplace transform first on both sides of equation (3), and we get

L [D_{t}^{n α} u (r, t)] + L [R [r] u (r, t) + N [r] u (r, t)] = L [g (r, t)]

(4)

Now, using the differentiation property of the Laplace transform, we have

\begin{matrix} L [u (r, t)] - \frac{1}{s^{n α}} \sum_{k = 0}^{n - 1} s^{(n α - k - 1)} u^{k} (r, 0) \\ + \frac{1}{s^{n α}} L (R [r] u (r, t) + N [r] u (r, t) - g (r, t)) = 0 \end{matrix}

(5)

We define the nonlinear operator

\begin{matrix} N [ϕ (r, t; q)] = L [ϕ (r, t; q)] - \frac{1}{s^{n α}} \sum_{k = 0}^{n - 1} s^{(n α - k - 1)} u^{k} (r, 0) \\ + \frac{1}{s^{n α}} L (R [r] ϕ (r, t) + N [r] ϕ (r, t) - g (r, t; q)) \end{matrix}

(6)

where $q \in [0, 1]$ is an embedding parameter and $ϕ (r, t; q)$ is the real function of $r, t, and q$ . By means of generalizing the traditional homotopy method, Liao^10–13 constructs the zero-order deformation equation

(1 - q) ℓ [ϕ (r, t; q) - u_{0} (r, t)] = ℏ qH (r, t) N [ϕ (r, t; q)],

(7)

where $ℏ$ is a nonzero auxiliary parameter, $H (r, t) \neq 0$ is an auxiliary function, $u_{0} (r, t)$ is an initial guess of $u (r, t)$ , and $ϕ (r, t; q)$ is an unknown function. It is important that one has great freedom to choose auxiliary thing in HATM. Obviously, when $q = 0$ and $q = 1$ in equation (7), it holds

ϕ (r, t; 0) = u_{0} (r, t), ϕ (r, t; 1) = u (r, t)

(8)

Thus, as $q$ increases from $0$ and $1,$ the solution varies from the initial guess $u_{0} (x, t)$ to the solution $u (r, t) .$ Expanding $ϕ (r, t; q)$ in Taylor’s series with respect to $q,$ we have

ϕ (r, t; q) = u_{0} (r, t) + \sum_{m = 1}^{\infty} q^{m} u_{m} (r, t)

(9)

where

u_{m} (r, t) = \frac{1}{m!} {\frac{\partial^{m} ϕ (r, t; q)}{\partial q^{m}} |}_{q = 0}

(10)

If the auxiliary linear operator $ℓ$ , the initial guess $u_{0} (r, t)$ , the auxiliary parameter $ℏ,$ and the auxiliary function are properly chosen, the series (9) converges at $q = 1,$ and we have

u (r, t) = u_{0} (r, t) + \sum_{m = 1}^{\infty} u_{m} (r, t)

(11)

which must be one of the solutions of original nonlinear equations.

Define the vectors

{\vec{u}}_{n} = {u_{0} (r, t), u_{1} (r, t), u_{2} (r, t), \dots, u_{n} (r, t)}

(12)

Differentiating equation (7) $m$ time with respect to embedding parameter $q$ and then setting $q = 0$ and finally dividing them by $m!,$ we obtain the mth-order deformation equation

L [u_{m} (r, t) - χ_{m} u_{m - 1} (r, t)] = ℏ qH (r, t) R_{m} ({\vec{u}}_{m - 1}, x, t)

(13)

Operating the inverse Laplace transform on both sides, we get

u_{m} (r, t) = χ_{m} u_{m - 1} (r, t) + ℏ q L^{- 1} [H (r, t) R_{m} ({\vec{u}}_{m - 1}, r, t)]

(14)

where

R_{m} ({\vec{u}}_{m - 1}, r, t) = \frac{1}{(m - 1)!} {\frac{\partial^{m - 1} ϕ (r, t; q)}{\partial q^{m - 1}} |}_{q = 0}

(15)

and

χ_{m} = {\begin{matrix} 0, m \leq 1 \\ 1, m > 1 \end{matrix}

In this way, it is easy to obtain $u_{m} (r, t)$ for $m \geq 1,$ and at Mth order, we have

u (r, t) = \sum_{m = 0}^{M} u_{m} (r, t)

(16)

when $M \to \infty$ we get an accurate approximation of the original equation (3)

\begin{matrix} u_{m} (r, t) = (χ_{m} + ℏ) u_{, m - 1} - ℏ (1 - χ_{m}) L^{- 1} (\sum_{k = 0}^{j - 1} \frac{u^{k} (x, 0)}{s^{k + 1}}) \\ + ℏ L^{- 1} (\frac{1}{s^{j α}} L (R [t] u (r, t) + N [t] u (r, t) - g (r, t))) \end{matrix}

(17)

In this way, it is easily to obtain $u_{m} (t)$ for $m \geq 1,$ and at Mth order, we have

u_{n} (r, t) = u_{0} (r, t) + lim_{k \to \infty} \sum_{m = 1}^{k} u_{n} (r, t)

(18)

where $u_{0} (r, t) = u (r, 0) + tu' (r, 0) + t^{2} u ″ (r, 0) + t^{3} u ″' (r, 0) + \dots + t^{j - 1} u^{(j - 1)} (r, 0) = \sum_{i = 0}^{j - 1} t^{i} u^{(i - 1)} (r, 0) .$

The nonlinear term $N_{n} [t] u_{n} (t)$ is expanded in terms of Adomian’s polynomials as

\begin{matrix} N [u (r, t)] = N (\sum_{m = 0}^{k} u_{m} (t)) \\ = \sum_{m = 0}^{k} A_{m} (u_{0}, u_{1}, u_{2}, u_{3}, \dots, u_{m}) \end{matrix}

(19)

where $A_{m}$ are Adomian’s polynomials which are calculated by algorithm (20) constructed by Adomian³⁰

= \frac{1}{n!} {[\frac{d^{k}}{d t^{k}} N (\sum_{p = 0}^{k} λ^{p} u_{p})]}_{λ = 0}, k \geq 0

(20)

Combining equations (17) and (20), we get

\begin{matrix} u_{m} (r, t) = (χ_{m} + ℏ) u_{m - 1} - ℏ (1 - χ_{m}) \sum_{i = 0}^{j - 1} t^{i} u^{(i - 1)} (r, 0) \\ + ℏ L^{- 1} (\frac{1}{s^{j α}} L (R [t] u_{m - 1} (t) + \sum_{k = 0}^{m - 1} A_{k} (u_{0}, u_{1}, u_{2}, u_{3}, \dots, u_{m}) - g (r, t))) \cdot m = 1, 2, 3, \dots \end{matrix}

(21)

From equation (21), we calculate the various $u_{m} (t)$ for $m \geq 1$ , and substituting these values in equation (18), we obtain the analytical approximate solution of equation (3). The novelty of our proposed algorithm is that a new correction functional (21) is constructed and the nonlinear term as a series of Adomian’s polynomials in equation (21) is expanded. We combine the new hybrid iterative algorithm for nonlinear fractional partial differential equations arising in science and engineering.

In order to select the optimal value of auxiliary parameter $(ℏ),$ the averaged residual error is introduced as

E_{α} = \frac{1}{(Q + 1) (W + 1)} \sum_{i = 0}^{Q} \sum_{j = 0}^{W} [Re s (\sum_{z = 0}^{α} u (i Δ x, j Δ t))]

(22)

with $Q = W = 1, Δ x = 0.05 / Q, and Δ t = 0.05 / W$ . It should be stated that the selection of optimal value of the auxiliary parameter in the HATM and HAM has an approximately similar procedures (or more details, see Cheng et al.,³¹ Liao^32,33, and Rashidi et al.³⁴). For the current approximation, the optimal value of $ℏ$ is given by the minimum value of $E_{α}$ corresponding to nonlinear algebraic equation

\frac{d E_{α}}{d ℏ} = 0 .

(23)

Numerical experiments

In this section, two examples on time-fractional wave equations are solved to demonstrate the performance and efficiency of the HATM.

Example 1

We consider the one-dimensional wave-like equation with variable coefficients^7–9 as

D_{t}^{2 α} u (x, t) = \frac{1}{2} x^{2} u_{xx}, 0 < x < 1, t > 0, 0 < α \leq 1

(24)

subject to the Neumann boundary conditions

u (0, t) = 0, u (1, t) = 1 + \sin h t

(25)

and the initial conditions

u (x, 0) = x, u_{t} (x, 0) = x^{2}

(26)

We choose linear operator as

£ [ϕ (x, t; q)] = L [ϕ (x, t; q)]

(27)

with property $£ [c] = 0,$ where $c$ is constant.

We now define a nonlinear operator as

N [ϕ (x, t; q)] = L [ϕ (x, t; q)] - \frac{x}{s} - \frac{x^{2}}{s^{2}} - \frac{s^{- 2 α}}{2} L [x^{2} ϕ_{xx}]

(28)

Thus, we obtain the mth_-order deformation equation

£ [u_{m} (x, t) - χ_{m} u_{m - 1} (x, t)] = ℏ R_{m} ({\vec{u}}_{m - 1}, x, t)

(29)

Operating the inverse Laplace transform on both sides in equation (29), we get

u_{m} (x, t) = χ_{m} u_{m - 1} (x, t) + ℏ q L^{- 1} [R_{m} ({\vec{u}}_{m - 1}, x, t)]

(30)

where

\begin{matrix} R_{m} ({\vec{u}}_{m - 1}, x, t) = L [u_{m - 1} (x, t)] \\ - (1 - χ_{m}) (\frac{x}{s} + \frac{x^{2}}{s^{2}}) - \frac{s^{- 2 α}}{2} L [x^{2} {(u_{m - 1})}_{xx}] \end{matrix}

(31)

Now the solution of mth-order deformation equations

\begin{matrix} u_{m} (x, t) = (χ_{m} + ℏ) u_{m - 1} - (1 - χ_{m}) (x + x^{2} t) \\ - \frac{ℏ}{2} L^{- 1} (s^{- 2 α} L [x^{2} {(u_{m - 1})}_{xx}), m \geq 1 \end{matrix}

(32)

Using the initial approximation $u_{0} (x, t) = u (x, 0) + t u_{t} (x, 0) = x + x^{2} t$ and the iterative scheme (32), we obtain the various iterates

u_{1} (x, t) = - \frac{x^{2} ℏ t^{2 α + 1}}{Γ (2 α + 2)}

u_{2} (x, t) = - \frac{ℏ (1 + ℏ) x^{2} t^{2 α + 1}}{Γ (2 α + 2)} + \frac{ℏ^{2} x^{2} t^{4 α + 1}}{Γ (4 α + 2)}

\begin{matrix} u_{3} (x, t) = - \frac{ℏ {(1 + ℏ)}^{2} x^{2} t^{2 α + 1}}{Γ (2 α + 2)} + \frac{{2 ℏ}^{2} (1 + ℏ) x^{2} t^{4 α + 1}}{Γ (4 α + 2)} \\ - \frac{ℏ^{3} x^{2} t^{6 α + 1}}{Γ (6 α + 2)} \end{matrix}

\begin{matrix} u_{4} (x, t) = - \frac{ℏ {(1 + ℏ)}^{3} x^{2} t^{2 α + 1}}{Γ (2 α + 2)} + \frac{{3 ℏ}^{2} {(1 + ℏ)}^{2} x^{2} t^{4 α + 1}}{Γ (4 α + 2)} \\ - \frac{{3 ℏ}^{3} (1 + ℏ) x^{2} t^{6 α + 1}}{Γ (6 α + 2)} + \frac{ℏ^{4} x^{2} t^{8 α + 1}}{Γ (8 α + 2)}, \dots \end{matrix}

Proceeding in this manner, the rest of the components $u_{n} (x, t)$ for $n \geq 5$ can be completely obtained, and the series solutions are thus entirely determined. Consequently, we obtain the solution of equation (24) in a series form

u (x, t) = u_{0} (x, t) + \sum_{m = 0}^{\infty} u_{m} (x, t)

(33)

However, mostly, the results given by the Laplace decomposition method and homotopy perturbation transform method converge to the corresponding numerical solutions in a rather small region. But, different from those two methods, the HATM provides us with a simple way to adjust and control the convergence region of solution series by choosing a proper value for the auxiliary parameter $ℏ .$ If we select $ℏ = - 1,$ then

\begin{matrix} u (x, t) = x + x^{2} \\ (t + \frac{t^{2 α + 1}}{Γ (2 α + 2)} + \frac{t^{4 α + 1}}{Γ (4 α + 2)} + \frac{t^{6 α + 1}}{Γ (6 α + 2)} + \dots) \\ = x + x^{2} \sin h (t^{α}, α) \end{matrix}

(34)

where

\sinh (t^{α}, α) = (t + \frac{t^{2 α + 1}}{Γ (2 α + 2)} + \frac{t^{4 α + 1}}{Γ (4 α + 2)} + \frac{t^{6 α + 1}}{Γ (6 α + 2)} + \dots)

When we choose $α = 1$ , then clearly we can conclude that the obtained solution $\sum_{m = 0}^{\infty} u_{m} (x, t)$ converges to the exact solution $u (x, t) = x + x^{2} \sin ht$ The above result is in complete agreement with Özis and Agırseven.⁹

Figure 1 shows the behavior of obtained solution $u (x, t)$ by the proposed method for different fraction Brownian motions $α = 0.7, 0.8, 0.9$ and for standard motion, that is, at $α = 1 .$ It is seen from Figure 1 that the solution obtained by new fractional homotopy analysis transform method (FHATM) increases very rapidly with increases in $t$ at the value of $x = 0.5, y = 0.5, and ℏ = - 1 .$

Figure 1.

Plot of solutions at different values of $α = 0.7, 0.8, 0.9, and 1$ for Example 1.

The convergence and rate of approximation for the HAM solution strongly depend on the value of auxiliary parameter $ℏ .$ Even if the initial approximation $u_{0} (x, t)$ , the auxiliary linear operator $L,$ and the auxiliary function $H (x, t)$ are given, we still have great freedom to choose the value of the auxiliary parameter $ℏ .$ So, the auxiliary parameter $ℏ$ provides us with an additional way to conveniently adjust and control the convergence region and rate of solution series. By means of the so-called $ℏ - curves$ it is easy to find out the so-called valid regions of $ℏ$ to gain a convergent solution series. When the valid region of $ℏ$ is a horizontal line segment, then the solution is converged.

Figure 2 shows the $ℏ - curve$ obtained from the 20th-order FHATM approximation solution of time-fractional wave equation (24) at $x = 0.5 .$ In our study, it is obvious from Figure 2 that the acceptable range of auxiliary parameter $ℏ$ is $- 2.0 \leq ℏ < 0 .$ We still have freedom to choose the auxiliary parameter according to $ℏ - curve$ . From Figure 2, the valid regions of convergence correspond to the line segments nearly parallel to the horizontal axis.

Figure 2.

Plot of $ℏ - curve$ for different values of $α$ for Example 1.

Example 2

In this example, we consider the following two-dimensional wave-like equation with variable coefficients^7–9 as

D_{t}^{2 α} u (x, y, t) = \frac{1}{12} (x^{2} u_{x x} + y^{2} u_{y y}), 0 < x, y < 1, t > 0, 0 < α \leq 1

(35)

subject to the Neumann conditions

{\begin{matrix} u_{x} (0, y, t) = 0, u_{x} (1, y, t) = 4 \cos h t \\ u_{y} (x, 0, t) = 0, u_{y} (x, 1, t) = 4 \sin h t \end{matrix}

(36)

with initial condition

u (x, y, 0) = x^{4}, u_{t} (x, y, 0) = y^{4}

(37)

Applying the Laplace transform on both sides in equation (35) and after using the differentiation property of Laplace transform, we get

L [u (x, y, t)] - s^{- 1} x^{4} - s^{- 2} y^{4} - \frac{1}{12} s^{- 2 α} L [x^{2} u_{xx} + y^{2} u_{yy}] = 0

(38)

We choose the linear operator as

£ [ϕ (x, y, t; q)] = L [ϕ (x, y, t; q)]

(39)

with property $£ [c] = 0,$ where $c$ is constant.

We now define a nonlinear operator as

\begin{matrix} N [ϕ (x, y, t; q)] = L [ϕ (x, y, t; q)] - (\frac{x^{4}}{s} + \frac{y^{4}}{s^{2}}) \\ - \frac{s^{- 2 α}}{12} L [x^{2} ϕ_{xx} + y^{2} ϕ_{yy}] \end{matrix}

(40)

Using equation (7), with assumption $H (x, y, t) = 1,$ we construct the zero-order deformation equation

(1 - q) £ [ϕ (x, y, t; q) - u_{0} (x, y, t)] = q ℏ N [ϕ (x, y, t; q)]

(41)

Obviously, when $q = 0$ and $q = 1$

ϕ (x, y, t; 0) = u_{0} (x, y, t), ϕ (x, y, t; 1) = u (x, y, t)

(42)

Thus, we obtain the mth-order deformation equation

£ [u_{m} (x, y, t) - χ_{m} u_{m - 1} (x, y, t)] = ℏ R_{m} ({\vec{u}}_{m - 1}, x, y, t)

(43)

Operating the inverse Laplace transform on both sides in equation (43), we get

u_{m} (x, y, t) = χ_{m} u_{m - 1} (x, y, t) + ℏ q L^{- 1} [R_{m} ({\vec{u}}_{m - 1}, x, y, t)]

(44)

where

\begin{matrix} R_{m} ({\vec{u}}_{m - 1}, x, y, t) = L [u_{m - 1} (x, y, t)] - (1 - χ_{m}) \\ (\frac{x^{4}}{s} + \frac{y^{4}}{s^{2}}) - \frac{s^{- 2 α}}{12} L [x^{2} {(u_{m - 1})}_{xx} + y^{2} {(u_{m - 1})}_{yy}] \end{matrix}

(45)

Now, the solution of mth_-order deformation equations

\begin{matrix} u_{m} (x, y, t) = (χ_{m} + ℏ) u_{m - 1} - ℏ (1 - χ_{m}) (x^{4} + y^{4} t) \\ - \frac{ℏ}{12} L^{- 1} (s^{- α} L [x^{2} {(u_{m - 1})}_{xx} + y^{2} {(u_{m - 1})}_{yy}]), m \geq 1 \end{matrix}

(46)

Using the initial approximation $u_{0} (x, y, t) = u (x, y, 0) + t u_{t} (x, y, 0) = x^{4} + t y^{4}$ and the iterative scheme (46), we obtain the various iterates

u_{1} (x, y, t) = - \frac{ℏ x^{4} t^{2 α}}{Γ (2 α + 1)} - \frac{ℏ y^{4} t^{2 α + 1}}{Γ (2 α + 2)}

\begin{matrix} u_{2} (x, y, t) = - \frac{ℏ (1 + ℏ) x^{4} t^{2 α}}{Γ (2 α + 1)} - \frac{ℏ (1 + ℏ) y^{4} t^{2 α + 1}}{Γ (2 α + 2)} \\ + \frac{ℏ^{2} x^{4} t^{4 α}}{Γ (4 α + 1)} + \frac{ℏ^{2} y^{4} t^{4 α + 1}}{Γ (4 α + 2)} \end{matrix}

\begin{matrix} u_{3} (x, y, t) = - \frac{ℏ {(1 + ℏ)}^{2} x^{4} t^{2 α}}{Γ (2 α + 1)} - \frac{ℏ {(1 + ℏ)}^{2} y^{4} t^{2 α + 1}}{Γ (2 α + 2)} \\ + \frac{{2 ℏ}^{2} (1 + ℏ) x^{4} t^{4 α}}{Γ (4 α + 1)} + \frac{{2 ℏ}^{2} (1 + ℏ) y^{4} t^{4 α + 1}}{Γ (4 α + 2)} \\ - \frac{ℏ^{3} x^{4} t^{6 α}}{Γ (6 α + 1)} - \frac{ℏ^{3} y^{4} t^{6 α + 1}}{Γ (6 α + 2)}, \dots \end{matrix}

Proceeding in this manner, the rest of the components $u_{n} (x, y, t)$ for $n \geq 4$ can be completely obtained, and the series solutions are thus entirely determined. Hence, the solution of the given problem is given as

u (x, y, t) = u_{0} (x, y, t) + \sum_{m = 0}^{\infty} u_{m} (x, y, t)

(47)

For $ℏ = - 1,$ we have the following solution

\begin{matrix} u (x, y, t) = x^{4} (1 + \frac{t^{2 α}}{Γ (2 α + 1)} + \frac{t^{4 α}}{Γ (4 α + 1)} + \dots) \\ + y^{4} (t + \frac{t^{2 α + 1}}{Γ (2 α + 2)} + \frac{t^{4 α + 1}}{Γ (4 α + 2)} + \frac{t^{6 α + 1}}{Γ (6 α + 2)} + \dots) \\ = x^{4} \sum_{k = 0}^{\infty} \frac{t^{2 k α}}{Γ (2 k α + 1)} + y^{4} \sum_{k = 0}^{\infty} \frac{t^{(2 k + 1) α}}{Γ ((2 k + 1) α + 1)} \\ = x^{4} \cos h (t^{α}, α) + y^{4} \sin h (t^{α}, α) \end{matrix}

(48)

where the functions $\sin h (z, α)$ and $\cos h (z, α)$ are defined as follows

\begin{matrix} \sinh (z, α) = \frac{E_{α} (z) - E_{α} (- z)}{2} and \cos h (z, α) \\ = \frac{E_{α} (z) + E_{α} (- z)}{2} \end{matrix}

(49)

As $α = 1,$ this series has the closed form $u (x, y, t) = x^{4} \cos h (t) + y^{4} \sin h (t),$ which is an exact solution of the classical heat equation. The above result is in complete agreement with Özis and Agırseven.⁹

Figure 3 shows the behavior of the approximate solution $u (x, y, t)$ for different fraction Brownian motions $α = 0.7, 0.8, 0.9$ and for standard motion, that is, at $α = 1$ for Example 2. It is seen from Figure 3 that the solution obtained by new FHATM increases very rapidly with increases in $t$ at the value of $x = 0.5, y = 0.5, and ℏ = - 1 .$ Figure 4 shows the $ℏ - curve$ obtained from the 20th-order FHATM approximation solution of two-dimensional time-fractional heat equation at $x = 0.5$ and $y = 0.5 .$ In our study, it is obvious from Figure 4 that the acceptable range of auxiliary parameter $ℏ$ is $- 2.0 \leq ℏ < 0 .$ We still have freedom to choose the auxiliary parameter according to $ℏ - curve$ . From Figure 4, the valid regions of convergence correspond to the line segments nearly parallel to the horizontal axis.

Figure 3.

Plot of solutions at different values of α = 0.7, 0.8, 0.9, and 1 for Example 2.

Figure 4.

Plot of $ℏ$ curves for different values of $α$ for Example 2.

Concluding remarks

In this article, we carefully proposed a reliable modification of the HAM which introduces a promising tool for solving wave equations of fractional order. Different from all other analytic methods, it provides us with a simple way to adjust and control the convergence region of solution series by introducing an auxiliary parameter. This work illustrates the validity and great potential of the HATM for fractional differential equations. An excellent agreement is achieved. The solution is very rapidly convergent by utilizing the modified HAM by modification of Laplace operator. It may be concluded that the modified HATM methodology is very powerful and efficient in finding approximate solutions as well as analytical solutions of many fractional physical models.

Footnotes

Acknowledgements

The authors are very grateful to the referees for carefully reading the paper and for their comments and suggestions which have improved the paper. The second author is highly grateful to the Department of Mathematics, National Institute of Technology, Jamshedpur, India, for the provision of some excellent facilities and research environment. This work is also supported by the Natural Science Foundation of Heilongjiang Province, China (Grant No. QC2015069).

Academic Editor: Xiao-Jun Yang

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This work is also supported by the Natural Science Foundation of Heilongjiang Province, China (Grant No. QC2015069)

References

Miller

Ross

. An introduction to the fractional calculus and fractional differential equations. New York: Wiley, 1993.

Oldham

Spanier

. The fractional calculus: integrations and differentiations of arbitrary order. New York: Academic Press, 1974.

Podlubny

. Fractional differential equations. San Diego, CA; New York: Academic Press, 1999.

Diethelm

. An algorithm for the numerical solution of differential equations of fractional order. Electron T Numer Ana 1997; 5: 1–6.

Holliday

Rundle

Tiampo

. Modification of the pattern informatics method for forecasting large earthquake events using complex Eigen factors. Tectonophysics 2006; 413: 87–91.

Manolis

Rangelov

. Non-homogeneous elastic waves in solid: notes on the vector decomposition technique. Soil Dyn Earthq Eng 2006; 26: 952–959.

Wazwaz

Gorguis

. Exact solutions for heat-like and wave-like equations with variable coefficients. Appl Math Comput 2004; 149: 15–29.

Shou

. Beyond adomian method: the variational iteration method for solving heat-like and wave-like equations with variable coefficients. Phys Lett A 2008; 372: 233–237.

Özis

Agırseven

. He’s homotopy perturbation method for solving heat-like and wave-like equations, with variable coefficients. Phys Lett A 2008; 372: 5944–5950.

10.

Liao

. The proposed homotopy analysis technique for the solution of nonlinear problems. PhD Thesis, Shanghai Jiao Tong University, Shanghai, China, 1992.

11.

Liao

Beyond perturbation: introduction to the homotopy analysis method. Boca Raton, FL: CRC Press/Chapman & Hall, 1992.

12.

Liao

. On the homotopy analysis method for nonlinear problems. Appl Math Comput 2004; 147: 499–513.

13.

Liao

. Homotopy analysis method: a new analytical technique for nonlinear problems. Commun Nonlinear Sci 1997; 2: 95–100.

14.

Vishal

Kumar

Das

. Application of homotopy analysis method for fractional swift Hohenberg equation-revisited. Appl Math Model 2012; 36: 3630–3637.

15.

Jafari

Golbabai

Seifi

. Homotopy analysis method for solving multi-term linear and nonlinear diffusion wave equations of fractional order. Comput Math Appl 2010; 59: 1337–1344.

16.

Zhang

Tang

. Homotopy analysis method for higher-order fractional integro-differential equations. Comput Math Appl 2011; 62: 3194–3203.

17.

Khan

Jamil

Ara

. Approximate solution of time fractional Schrodinger equation via homotopy analysis method. ISRN Math Phys 2012; 2012: 197068.

18.

Ghotbi

Bararnia

Domairry

. Investigation of powerful analytical method in to natural convection boundary layer flow. Commun Nonlinear Sci 2009; 15: 2222–2228.

19.

Abbasbandy

Shivanian

Vajravelu

. Mathematical properties of h-curve in the frame work of the homotopy analysis method. Commun Nonlinear Sci 2011; 16: 4268–4275.

20.

Abbasbandy

. Homotopy analysis method for the Kawahara equation. Nonlinear Anal: Real 2010; 11: 307–312.

21.

Abbasbandy

. Approximate solution for the nonlinear model of diffusion and reaction in Porous catalysts by means of the homotopy analysis method. Chem Eng J 2008; 136 144–150.

22.

Domairry

Fazeli

. Homotopy analysis method to determine the fin efficiency of convective straight fins with temperature-dependent thermal conductivity. Commun Nonlinear Sci 2009; 14: 489–499.

23.

Wazwaz

. The combined Laplace transform-Adomian decomposition method for handling nonlinear Volterra integro-differential equations. Appl Math Comp 2010; 216: 1304–1309.

24.

Khan

Gondal

Kumar

. A new analytical procedure for nonlinear integral equation. Math Comput Model 2012; 55: 1892–1897.

25.

Khan

Faraz

Kumar

. A coupling method of homotopy perturbation and Laplace transform for fractional models. UPB Sci Bull A 2012; 74: 57–68.

26.

Kumar

Kocak

Yildirim

. A fractional model of gas dynamics equation and its approximate solution by using Laplace transform. Z Naturforsch 2012; 67a; 389–396.

27.

Kumar

Yildirim

Khan

. A fractional model of the diffusion equation and its analytical solution using Laplace transform. Sci Iran B 2012; 19: 1117–1123.

28.

Singh

Kumar

. New treatment of fractional Fornberg–Whitham equation via Laplace transform. Ain Shams Eng J. Epub ahead of print 16 January 2013. DOI: 10.1016/j.asej.2012.11.009.

29.

Khan

Gondal

Hussain

. A new comparative study between homotopy analysis transform method and homotopy perturbation transform method on a semi-infinite domain. Math Comput Model 2012; 55; 1143–1150.

30.

Adomian

. A review of the decomposition method in applied mathematics. J Math Anal Appl 1988; 135; 501–544.

31.

Cheng

Zhu

Liao

. An explicit series approximation to the optimal exercise boundary of American put options. Commun Nonlinear Sci 2010; 15; 1148–1158.

32.

Liao

. An optimal homotopy-analysis approach for strongly nonlinear differential equations. Commun Nonlinear Sci 2010; 15; 2003–2016.

33.

Liao

. Homotopy analysis method in nonlinear differential equations. Beijing, China: Higher Education Press, Springer, 2012.

34.

Rashidi

Freidoonimehr

Hosseini

. Homotopy simulation of nanofluid dynamics from a non-linearly stretching isothermal permeable sheet with transpiration. Meccanica 2013; 49: 1–14.