On the Modified Laplace Homotopy Perturbation Method for Solving Damped Modified Kawahara Equation and Its Application in a Fluid

Abstract

The manuscript solves a modified Kawahara equation (mKE) within two cases with and without a damping term by applying the Laplace homotopy perturbation method (LHPM). Since the damped mKE is non-integrable (i.e., it does not have analytic integrals) and does not have exact initial conditions, this challenge makes many numerical methods fail to solve non-integrable equations. In this article, we suggested a new modification at LHPM by setting a perturbation parameter and an embedding parameter as the damping parameter and using the initial condition for mKE as the initial condition for non-damped mKE. The results proved that this mathematical approach is an effective method for solving damped mKE. Thus, we believe that the presented method will be helpful for solving many non-integrable equations that describe phenomena in sciences, such as nonlinear symmetrical wave propagation in plasma.

1. Introduction

Nonlinear evolution equations have played an important role in several branches of science, such as engineering physics, geochemistry and fluid mechanics, in recent decades [1,2,3,4,5]. Nonlinear differential equations (NLDEs) have been studied widely by many mathematicians in order to obtain solutions, either numerical, analytical, equivalent, or exact solutions. One of the well-known NLDEs is the Kawahara equation in the following form [6,7]

$$u_t+au{u}_x+b{u}_{3x}-k{u}_{5x}=0,$$

(1)

note that the Korteweg-De Vries KdV equation is a special case of Equation (1) when $k=0$ [8],

$$u_t+au{u}_x+b{u}_{3x}=0.$$

(2)

The modified Kawahara Equation (mKE) takes the form [9]

$$u_t+a{u}^2{u}_x+b{u}_{3x}-k{u}_{5x}=0.$$

(3)

where $a,b$ and k are constants. The mKE depicts shallow water waves with surface pressure or magneto wave propagation in fluid media [10,11]. The mKE has been solved by numerical methods in the literature. For example, it was solved numerically by the Crank–Nicolson discretization method [12], by the Kernel smoothing method [13], and by the septic B-spline collocation method [14].

In this work, we aim to solve mKE by one of the trending methods, which is the Laplace homotopy perturbation method (LHPM). Actually, this method is a combination of the Laplace transform and homotopy perturbation method. Laplace transform is a powerful method that transfers the equation to a simple form and uses the initial condition. The homotopy perturbation method (HPM) was introduced by Ji He. It is a combination of the popular perturbation technique and the homotopy. It is useful because it is applied without any discretization, transformation, small parameter, linearization, or assumption. It only requires a few iterations to reach highly accurate solutions. LHPM is an iteration method and combines the benefits of Laplace transform and HPM.

In the literature, many equations have been solved by LHPM such as the Kawahara equation [15], gas dynamics equations [16], the nonlinear Schrödinger equation with Harmonic Oscillator [17], reaction–diffusion equations [18], nonlinear space-time fractional Fokker–Planck equations [19], Phi-four equation [20], Black–Scholes equation [21], generalized Sylvester matrix equation [22], the foam drainage equation [23], higher order linear and nonlinear boundary value problems (BVPs) [24], Duffing-van der Pol’s cubic-quintic equation [25] and many other problems. Through reviewing the method, we can realize that the algorithm is straightforward and has the ability to solve differential equations in many applications. Thus, we will use the aforementioned method to obtain a very accurate solution. Then, we will add the effect of the damping term into mKE, which converts the Equation (3) to a linear damped non-integrable as follows

$$u_t+a{u}^2{u}_x+b{u}_{3x}-k{u}_{5x}+R{u}_t=0,$$

(4)

where $R<<1$ is a constant [26]. The addition of a damping term results from accounting for viscosity impact in the physical model. It is very important to consider some properties of the problem. Actually, some experimental applications in a complex plasma show that the properties of propagation on collision of the cylindrical dust-acoustic solitons (CDASs), such as velocity, width and density, are different between the bounded nonplanar and unbounded planar geometry [26,27]. In the case of planar geometry, $R=0$, while in the case of nonplanar geometry, such as cylindrical and spherical, $R\ne 0$ [28]. The damping term complicates the mKE and makes Equation (4) a non-integrable equation [29]. Thus, Equation (4) does not have an exact solution or exact initial condition, and this is the reason behind the failure of most numerical methods to reach the approximate solution of Equation (4).

The novelty of this paper is an approach to modify LHPM in order to solve a non-integrable equation subject to the initial condition of the integrable version of the considered problem. The LHPM is modified by considering a perturbation and an embedding parameter as the damping parameter.

This article is organized as follows: Section 2 contains important definitions to describe the LHPM. Section 3 describes the LHPM in detail and uses LHPM to find the numerical solution of mKE as well as prove the accuracy of the method. Section 4 proposes the modified LHPM to solve linear damping mKE, and the final section is the conclusion of the work.

2. Definitions

This section will present some important definitions in order to explain the LHPM steps.

2.1. Definition of the Laplace Transformation

The Laplace transformation of a function $v\left(t\right)$ is denoted by $\mathcal{L}\left\{v\right(t\left)\right\}$ or $V\left(s\right)$ for $t\ge 0$ and is given by the following integral

$$V\left(s\right)={\int }_0^{\infty }{e}^{-st}v\left(t\right)dt.$$

(5)

An inverse Laplace transformation of a function $V\left(s\right)$ is $v\left(t\right)$. If $v\left(t\right)$ exists for function $V\left(s\right)$ such that $\mathcal{L}\left\{v\right(t\left)\right\}$ = $V\left(s\right)$, then

$${\mathcal{L}}^{-1}\left\{V\left(s\right)\right\}=v\left(t\right).$$

(6)

2.2. He’s Polynomial

The formula to compute He’s polynomials is

$$H_i(u_0,u_1,\cdots ,u_i)=\frac{1}{i!}\frac{{\partial }^i}{\partial p^i}{\left(N\left(\sum _{i=0}^{\infty }{u}_i{p}^i\right)\right)}_{p=0},i=0,1,.$$

(7)

2.3. Definition of Homotopy Perturbation Theory

To characterize the HPM, we use a generalized formula,

$$L\left(u\right)=0,$$

(8)

where L is any integral operator, and $H(u,p)$ is a convex homotopy that is defined as

$$H(u,p)=(1-p)F\left(u\right)+pL\left(u\right),$$

(9)

where $F\left(u\right)$ is a functional operator that has easily acquired known values $v_0$, and p is embedding parameter such that $p\in [0,1]$. It is obvious that

$$H(u,0)=F\left(u\right)andH(u,1)=L\left(u\right).$$

(10)

This shows that $H(u,p)$ follows an implicit curve from a starting point $H(v_0,0)$ to a solution function $H(f,1)$. The solution might be considered as

$$u=\sum _{i=0}^{\infty }{u}_i{p}^i.$$

(11)

If $p\to 1$, then (11) becomes the approximate solution

$$f=\underset{p\to 1}{lim}u=\sum _{i=0}^{\infty }{u}_i.$$

(12)

In the majority of cases, the series (12) is convergent to the solution and the rate of convergence is based on $L\left(u\right)$. Assume the solution for Equation (12) is unique. The solution u can be written as a polynomial of p using the homotopy perturbation method with the homotopy parameter p as follows

$$u\left(x\right)=\sum _{i=0}^{\infty }{u}_i{p}^i.$$

(13)

In the same manner, He’s HPM considers the nonlinear term ($N\left(u\right)$) as follows

$$N\left(u\right)=\sum _{i=0}^{\infty }{H}_i{p}^i,$$

(14)

where $H_i^{\prime }s$ are $H{e}^{’}s$ polynomials and are found by the Formula (7).

3. Description of Laplace Homotopy Perturbation Method

This section describes the LHPM briefly. First, we consider the following PDEs

$$u_t(x,t)+L\left(u(x,t)\right)+N\left(u(x,t)\right)=g(x,t),$$

(15)

where u is an undefined function, L and N are linear and nonlinear operators, respectively, and g is a source term. Equation (15) is subject to the following initial condition

$$u(x,0)=f\left(x\right).$$

(16)

Applying the LHPM algorithm to Equation (15) is summed up in the following steps:

• Take the Laplace transform ($\mathcal{L}$) for Equation (15)

  $$\mathcal{L}\left[u_t(x,t)\right]+\mathcal{L}\left[L\left(u(x,t)\right)\right]+\mathcal{L}\left[N\left(u(x,t)\right)\right]=\mathcal{L}\left[g(x,t)\right].$$

  (17)
• Use the differentiation property of the Laplace transform

  $$s\mathcal{L}\left[u(x,t)\right]-u(x,0)+\mathcal{L}\left[L\left(u\right(x,t\left)\right)\right]+\mathcal{L}\left[N\left(u\right(x,t\left)\right)\right]=\mathcal{L}\left[g(x,t)\right],$$

  (18)

  or in form

  $$\mathcal{L}\left[u(x,t)\right]=\frac{1}{s}u(x,0)-\frac{1}{s}\mathcal{L}\left[L\left(u(x,t)\right)\right]-\frac{1}{s}\mathcal{L}\left[N\left(u(x,t)\right)\right]+\frac{1}{s}\mathcal{L}\left[g(x,t)\right].$$

  (19)
• Operate the inverse Laplace transform (${\mathcal{L}}^{-1}$) in Equation (19)

  $$u(x,t)=G(x,t)-{\mathcal{L}}^{-1}\left[\frac{1}{s}\left(\mathcal{L}\left[Lu(x,t)\right]+\mathcal{L}\left[Nu(x,t)\right]\right)\right],$$

  (20)

  where $G(x,t)=u(x,0)+{\mathcal{L}}^{-1}\left[\frac{1}{s}\mathcal{L}\left[g(x,t)\right]\right]$.
• Use the homotopy approach that was introduced by Liao and assume $u={\sum }_{i=0}^{\infty }{u}_i{p}^i$, where p is the an embedding parameter such that $0<p<1$,

  $$\sum _{i=0}^{\infty }{u}_i{p}^i=G(x,t)-p{\mathcal{L}}^{-1}\left[\frac{1}{s}\left(\mathcal{L}\left[L\left(\sum _{i=0}^{\infty }{u}_i{p}^i\right)\right]+\mathcal{L}\left[N\left(u\right(x,t\left)\right)\right]\right)\right].$$

  (21)
• Substitute $N\left(u\right(x,t\left)\right)$ = ${\sum }_{i=0}^{\infty }{H}_i\left(u\right)p^i$ into Equation (21) as follows

  $$\sum _{i=0}^{\infty }{u}_i{p}^i=G(x,t)-p{\mathcal{L}}^{-1}\left[\frac{1}{s}\left(\mathcal{L}\left[L\left(\sum _{i=0}^{\infty }{u}_i{p}^i\right)\right]+\mathcal{L}\left[N\left(\sum _{i=0}^{\infty }{H}_i\left(u\right)p^i\right)\right]\right)\right],$$

  (22)

  where $H_i$ is called He${}^{’}$s polynomials.
• Equate the coefficients of $p^i$ on the right side to left in Equation (22) as follows:
  $p^0:u_0=G(x,t)$,
  $p^1:u_1=-{\mathcal{L}}^{-1}\left[\frac{1}{s}(\mathcal{L}\left[L{u}_0\right]+\mathcal{L}\left[N{H}_0\right]\right]$,
  ⋮
  $p^{i+1}:u_{i+1}=-{\mathcal{L}}^{-1}\left[\frac{1}{s}(\mathcal{L}\left[L{u}_i\right]+\mathcal{L}\left[N{H}_i\right]\right]$.
  If $p\to 1$, then the form of the approximate solution becomes

  $$u(x,t)=\underset{p\to 1}{lim}u=\sum _{i=0}^{\infty }{u}_i.$$

  (23)
  It may eventually lead to the exact solution of Equation (15)

Numerical Solution for Modified Kawahara Equation by LHPM

In this subsection, we employ the LHPM into mKE,

$$u_t+a{u}^2{u}_x+b{u}_{3x}-k{u}_{5x}=0,$$

(24)

subject to the initial conditions [9]

$$u(x,0)=u_0=-\frac{3b}{\sqrt{10ak}}sec{h}^2\left[\frac{1}{2}\sqrt{\frac{b}{5k}}\left(x\right)\right].$$

(25)

The linear part is $L\left(u\right)=-b{u}_{3x}+k{u}_{5x}$, and the nonlinear part is $N\left(u\right)=-a{u}^2{u}_x$. Applying LHPM into Equation (24) and following the steps in a previous section with MATLAB help leads to the following value of $u_i$:

$$\begin{array}{cc}\hfill u_0=& -\frac{3b}{\sqrt{10ak}}sec{h}^2\left[\frac{1}{2}\sqrt{\frac{b}{5k}}\left(x\right)\right],\hfill \\ \hfill u_1=& -{\mathcal{L}}^{-1}\left[\frac{1}{s}\left[-\mathcal{L}\left[a\left(H_0(x,t)\right)\right]-\mathcal{L}\left[b\left(u_0{(x,t)}_{3x}\right)\right]+\mathcal{L}\left[k\left(u_0{(x,t)}_{5x}\right)\right]\right]\right],\hfill \\ \hfill u_2=& -{\mathcal{L}}^{-1}\left[\frac{1}{s}\left[-\mathcal{L}\left[a\left(H_1(x,t)\right)\right]-\mathcal{L}\left[b\left(u_1{(x,t)}_{3x}\right)\right]+\mathcal{L}\left[k\left(u_1{(x,t)}_{5x}\right)\right]\right]\right],\hfill \\ \hfill ⋮\end{array}$$

where $H_k(x,t)$ is defined for nonlinear terms as

$$\begin{array}{cc}\hfill H_0\left(u_0\right)& =u_0^2{u}_0^{\prime },\hfill \\ \hfill H_1(u_0,u_1)& =2u_0{u}_1{u}_0^{\prime }+u_0^2{u}_1^{\prime },\hfill \\ \hfill H_2(u_0,u_1,u_2)& =(u_1^2+2u_0{u}_2)u_0^{\prime }+2u_0{u}_1{u}_1^{\prime }+u_0{u}_2^{\prime },\hfill \\ \hfill ⋮\\ \hfill H_k(u_0,u_1,\cdots ,u_k)=& \frac{1}{k!}\frac{{\partial }^k}{\partial p^k}\left[a\left({\left(\sum _{i=0}^{\infty }{u}_k{p}^k\right)}^2\left(\sum _{i=0}^{\infty }{u}_k{p}^k\right)\right)\right].\hfill \end{array}$$

The parts of the series can be gained in this manner. The sequence solution formulation is expressed as follows

$$u_t=\sum _{i=0}^9{u}_i.$$

(26)

We used ten iteration terms to obtain and plot the solution in Figure 1 and Figure 2.

The absolute error is defined as

$$Error=\left|U\right(x,t)-u(x,t\left)\right|,$$

where U is the exact solution.

Table 1. The absolute difference between the numerical and exact solution at $t=2,4,6,8,10$ and $1\le x\le 10$.

x/t	2	4	6	8	10
1	$1.366\times {10}^{-14}$	$2.892\times {10}^{-11}$	$2.522\times {10}^{-9}$	$5.991\times {10}^{-8}$	$6.964\times {10}^{-7}$
2	$3.997\times {10}^{-15}$	$7.391\times {10}^{-12}$	$5.639\times {10}^{-10}$	$1.145\times {10}^{-8}$	$1.102\times {10}^{-7}$
3	$4.774\times {10}^{-15}$	$1.038\times {10}^{-11}$	$9.290\times {10}^{-10}$	$2.268\times {10}^{-8}$	$2.710\times {10}^{-7}$
4	$1.221\times {10}^{-15}$	$2.001\times {10}^{-12}$	$1.573\times {10}^{-10}$	$3.310\times {10}^{-9}$	$3.322\times {10}^{-8}$
5	$9.437\times {10}^{-16}$	$2.397\times {10}^{-12}$	$2.127\times {10}^{-10}$	$5.153\times {10}^{-9}$	$6.123\times {10}^{-8}$
6	$5.551\times {10}^{-17}$	$2.873\times {10}^{-13}$	$2.812\times {10}^{-11}$	$7.483\times {10}^{-10}$	$9.732\times {10}^{-9}$
7	$1.388\times {10}^{-16}$	$2.565\times {10}^{-13}$	$2.199\times {10}^{-11}$	$5.139\times {10}^{-10}$	$5.878\times {10}^{-9}$
8	$9.714\times {10}^{-17}$	$1.563\times {10}^{-13}$	$1.386\times {10}^{-11}$	$3.361\times {10}^{-10}$	$4.005\times {10}^{-9}$
9	$1.388\times {10}^{-17}$	$4.523\times {10}^{-14}$	$4.089\times {10}^{-11}$	$1.012\times {10}^{-10}$	$1.233\times {10}^{-9}$
10	$2.776\times {10}^{-17}$	$3.705\times {10}^{-15}$	$3.658\times {10}^{-13}$	$9.927\times {10}^{-12}$	$1.316\times {10}^{-10}$

shows the error in different values of t and x.

The error in Table 1 is very small and is based on the number of iteration. This confirmed the results in reference [15], which solved KE and proved that the results are more accurate than the optimal homotopy asymptotic method (OHAM) homotopy perturbation and variational iteration method (VHPM), and homotopy perturbation method (HPM).

Figure 1 presents the surface of the exact solution with the surface of the numerical solution, while Figure 2 plots u and U in fixed t and proves the high accuracy of LHPM. Herein, we found that the LHPM is more effective and accurate in solving integrable mKE.

4. Description of Improved Laplace Homotopy Perturbation Method

In sum, the non-integrable equation does not have an exact solution; thus, we do not have an initial condition. The suggesting technique assumes the solution is a finite series in terms of the damping parameter and the initial condition $u_0$ is for well-known integrable equations. Consider the following form of nonlinear equation

$$u_t+Lu+Nu=0,$$

(27)

where u is an undefined function, and $Lu$ and $Nu$ are linear and nonlinear parts, respectively. Assume the non-integrable equation as follows

$$u_t+Lu+Nu+wu=0,$$

(28)

where $wu$ is a damping parameter that causes the non-integrability for Equation (27). Then, we apply the LHPM as in the previous section and take an embedding parameter as the damping parameter.

Numerical Simulation for Linear Damped Modified Kawahara Equation

The equation of the damped mKE is given as

$$u_t+a{u}^2{u}_x+b{u}_{3x}-k{u}_{5x}+wu=0,$$

(29)

subject to [9]

$$u(x,0)=-\frac{3b}{\sqrt{10ak}}sec{h}^2\left[\frac{1}{2}\sqrt{\frac{b}{5k}}\left(x\right)\right].$$

(30)

In order to solve non-integrable equations, we consider the embedding parameter as the damping parameter, which is usually $0<w<1$. Thus, the solution is defined as $u={\sum }_{i=0}^N{w}^i{u}^i$. Applying LHPM to Equation (29), following the steps in the previous subsection and using MATLAB yields the following scheme:

$$u^{\left(0\right)}=-\frac{3b}{\sqrt{10ak}}sec{h}^2\left[\frac{1}{2}\sqrt{\frac{b}{5k}}\left(x\right)\right],$$

(31)

$$u^{\left(1\right)}(x,t)=-{\mathcal{L}}^{-1}\left[\frac{1}{s}\mathcal{L}\left[a{H}^{\left(k\right)}+b{u}_{3x}^{\left(k\right)}-k{u}_{5x}^{\left(k\right)}\right]\right],$$

(32)

$$u^{(k+1)}(x,t)=-{\mathcal{L}}^{-1}\left[\frac{1}{s}\mathcal{L}\left[a{H}^{\left(k\right)}+b{u}_{3x}^{\left(k\right)}-k{u}_{5x}^{\left(k\right)}+w{u}^{(k-1)}\right]\right],i=1,2,3...,N-1.$$

(33)

Figure 3 demonstrates the dependence of the modified Kawahara solution on the damping coefficient w. We obtain the rarefaction wave or negative wave. In Figure 4, it is obvious that when the damping term w is increased, the amplitude of the wave decreases. Equation (4) has been solved by the ansatz method, and this method used a suitable hypothesis based on the exact solutions for non-damped mKE [30]. However, the ansatz method is an analytical method, which assumes the solution in ansatz and requires some specific computations to find unknown functions in the ansatz. The comparison between the ansatz method and LHPM reveals that LHPM is a numerical method and only requires initial conditions.

5. Conclusions

In this manuscript, the LHPM has been used to find an accurate numerical solution for the modified Kawahara problem. The accuracy of the numerical solutions was tested by determining the maximum absolute error throughout the whole space-time domain. We realized that LHPM is a very powerful method for solving the integrable equation. By considering the damping term in mKE, we obtain the damped mKE, which is a non-integrable equation. The non-integrable equation does not have an exact initial solution, which makes many numerical methods fail to approximate the solution. Thus, the LHPM numerical methods, similar to other numerical methods in the literature, cannot find the approximate solution for non-integrable equations. This article modified LHPM to be able to solve non-integrable damped mKE by considering the perturbation and embedding parameter as the damping parameter in the LHPM algorithm while using the initial conditions for the integrable mKE. In addition, the obtained result can help many researchers investigate numerous phenomena in plasma physics, notably plasma oscillations [31].

In future work, we aim to improve other numerical methods to solve non-integrable equations. Furthermore, many non-integrable equations in the literature can be solved by the improved LHPM.