Comparison of modified ADM and classical finite difference method for some third-order and fifth-order KdV equations


 The KdV equation, which appears as an asymptotic model in physical systems ranging from water waves to plasma physics, has been studied. In this paper, we are concerned with dispersive nonlinear KdV equations by using two reliable methods: Shehu Adomian decomposition method (STADM) and the classical finite difference method for solving three numerical experiments. STADM is constructed by combining Shehu’s transform and Adomian decomposition method, and the nonlinear terms can be easily handled using Adomian’s polynomials. The Shehu transform is used to accelerate the convergence of the solution series in most cases and to overcome the deficiency that is mainly caused by unsatisfied conditions in other analytical techniques. We compare the approximate and numerical results with the exact solution for the two numerical experiments. The third numerical experiment does not have an exact solution and we compare profiles from the two methods vs the space domain at some values of time. This study provides us with information about which of the two methods are effective based on the numerical experiment chosen. Knowledge acquired will enable us to construct methods for other related partial differential equations such as stochastic Korteweg-de Vries (KdV), KdV-Burgers, and fractional KdV equations.


Introduction
Nonlinear partial differential equations (PDEs) have a significant role in various scientific and engineering fields. Since the discovery of solitons in 1965 by Zabusky and Kruskal [1], numerous nonlinear PDEs have been derived and extensively applied in many branches of sciences; for example, they appear in fluid mechanics, chemical kinetics, plasma physics, nonlinear optics, condensed matter, solid-state physics, theory of turbulence, ocean dynamics, biophysics, and star formation and others.
The well-known Korteweg-de Vries (KdV) equation is a generic evolutionary nonlinear dispersive PDE, which models weakly nonlinear long waves, incorporating a certain balance of leading-order nonlinearity and dispersion. The KdV equation describes solitary water waves (also called solitons) in a shallow water domain and it is given by [2] u uu u 0, t x xxx where u x t , : ( ) × ⟶ + . Equation (1) arises in magnetohydrodynamics, shallow water, stratified flows, internal waves, and plasma physics. The integrability of equation (1) guaranties the existence of an infinite number of invariants of motion and the first three invariants are as follows [3]: There exists a recurrence relation which allows us to construct higher order invariants n , for all n 3 > [4] , , .
The invariant 3 has a special meaning since it is also a Hamiltonian for the KdV equation [4] u u u u x δ δ , , This Hamiltonian structure sets the basis to construct a symplectic discretization for equation (1) [3,4].
There are various techniques for solving KdV equations in the literature; for example for fifth-order KdV equations, we have modified variation iteration technique [5], differential quadrature method [6], and Bspline method [7]. The authors in [8] studied some integral transform-based decomposition methods to solve certain nonlinear PDEs. Jena and Chakraverty [9] used Sumudu transformation method for the Bagley-Torvik equation and the authors in [10] employed homotopy perturbation Elzaki transform method for time-fractional Navier-Stokes equations and also Q-homotopy analysis Aboodh transform method [11]. Sadeway et al. [12] solved three different extended fifth-order Korteweg-de-Vries (efKdV) equations by using some ansatz method. Nuruddeen et al. [13] studied a class of fifth-order KdV equations by formulating suitable novel hyperbolic and exponential ansatz. The author in [14] considered a time-fractional order derivative defined in the conformable fractional derivative sense. The authors in [15] recently investigated classical and multisymplectic schemes for linearized KdV equations using some numerical methods and dispersion analysis was studied. Also, some semi-analytic methods were applied to study dispersive KdV equations by Appadu and Kelil [16]. The paper is organized as follows. The three numerical experiments considered are described in Section 2. Section 3 revisits some preliminaries of Shehu's transform method. Section 4 deals with an outline of Shehu Adomian decomposition method (STADM) and its applications to solve homogeneous and non-homogeneous nonlinear dispersive KdV equations. In Section 4, we present results using STADM. Section 5 is devoted to the application of the classical finite difference method (FDM) to the considered numerical experiments. In this section, we study the stability, consistency, and present numerical results of the considered experiments. The performance of the considered numerical schemes is gauged by computing L L , 1 2 and L ∞ error norms. Besides, we have also solved a nonlinear KdV equation (with no known exact solution) using Zabusky-Kruskal's numerical scheme, and the numerical results obtained for this experiment signify the outbreak of solitons due to perturbation of the parameters involved in the KdV equation. We highlight salient features of the paper in Section 6.
The novelty of this study relies on the comparative study of newly generalized integral transform, i.e., Shehu's integral transform combined with Adomian decomposition method (ADM) together with classical FDM to effectively solve some nonlinear third-order and fifth-order dispersive KdV equations. The importance of this work not only signifies the applicability of STADM to earn realistic series solutions but also to compare it with the classical FDM in order to witness consistent accuracy throughout the domain of the given problems. The newly proposed methods in this work can also be applied to many complicated linear and nonlinear problems since STADM does not require linearization, discretization, or perturbation.

Numerical experiment
We considered three numerical experiments which are described as follows: (i) Solve the homogeneous nonlinear dispersive KdV equation [17] u uu u 6 0 , t x xxx with x t π , 0,2 0,0.10 The initial condition is u x x , 0 ( ) = and the boundary conditions are u t u t t u t 0, 0, 0, 1 1 6 , 0, 0.
The exact solution is u x t , x t 1 6 ( ) = + . (ii) Solve the non-homogeneous nonlinear dispersive KdV equation [18] u uu u x t x t x cos sin sin 2 2 , , for two cases: Periodic boundary conditions are used. We note that there is no known exact solution for numerical experiment 3.

Preliminaries: Shehu transform
Recently, it became known that integral transform modifications are vital to solve numerous types of differential equations, for example, Laplace and Fourier transforms [19]. New generalized forms of Laplace-type transform methods such as Sumudu and Elzaki transforms have been introduced to solve PDEs, for example [20,21]. The authors in [22] introduced a new integral transform, which generalizes Laplace transform, Sumudu transform, Elzaki transform, and Yang transform, which is known as Shehu's transform. Shehu's transform has become an interesting topic of research for numerous researchers because not many research studies have been carried out on this topic and also due to the need to determine an efficient method for solving nonlinear PDEs [22]. Shehu transform coupled with semi-analytic methods such as ADM has been rarely discussed in the literature. Laplace, Sumudu, and Elzaki transform modifications have been successfully introduced in solving some PDEs as well as fractional differential equations (FDEs) such as the fractional Burgers equation [23].
To the best of our knowledge, Shehu transform coupled with ADM has never been applied before for solving nonlinear dispersive KdV-type equations. Therefore, the results in this work are quite novel.
and is denoted by Definition 2. The function v t ( ) in equation (6) is called inverse Shehu transform of s ρ , ( ) and is denoted by , f o r a l l 0 .
The inversion-type formula for Shehu transform is given by [22] v t sρ πi ρ st ρ s ρ s , where s and ρ are variables in the definition of Shehu transform and α is a real constant and the integral in equation (9) is taken along s α = in the complex plane s x iy = + .
and of exponential order α for t β > .
is the nth derivative of the function v t ( ) ∈ A with respect to "t," then its Shehu transform is defined by [22] v t s ρ s ρ s ρ v , 0 .
When n 1, 2, 3 = in equation (10), we obtain the following derivatives with respect to t: For the proofs of equations (10) and (11), we refer the reader to [22]. Linearity property and scaling property of Shehu transform are given, respectively, by For the reader's convenience, we just give the proof for equation (12).

STADM
We consider a general nonlinear non-homogeneous PDE of the form where is a general differential operator, and is a boundary operator, f r ( ) is a known analytic function, Θ is the boundary of the domain Ω, and u n ∂ ∂ the is directional derivative in the outwarding normal n to Ω. Equation (15) can be rewritten as Comparison of modified ADM and classical FDM  381 is a linear operator that includes partial derivatives with respect to x, is a nonlinear operator, and g is a nonhomogeneous term, which is u-independent.
STADM consists of applying Shehu transform on both sides of equation (16), obtaining Using the differentiation property of the Shehu transform and given initial conditions, we have In the homogeneous case, g x t , 0 ( ) = . Therefore, we have Now, applying the inverse Shehu transform to both sides of equation (18) gives The nonlinear term can be decomposed, using Adomian polynomials [24,25], as and the series solution of u x t , ( ) via STADM takes the form [26][27][28] , a n d , , , .
x xxx Taking inverse Shehu transform on both sides of equation (25), we obtain By applying the aforesaid decomposition method, we have , , , .
On comparing both sides of equation (27), we get The first few components of Adomian polynomials u u u , , , (cf. [24,25]) are given by Using the iteration formulae (28) and Adomian polynomials in (29), we obtain    The solution in series form (31) is a convergent series, which gives the exact solution This closed form solution for equation (3) agrees with the one obtained by the Elzaki transform method given in [21].
Plots of exact and numerical solutions vs x are displayed in Figure 1. We obtain plots of absolute error vs x at four different values of time in Figure 2. We also compare the absolute and relative errors at some values of x at four different times in Table 2. The absolute errors and relative errors at a given spatial node x j and given time t n are obtained as follows:   Table 2 shows that STADM is very effective at small and medium propagation time but becomes less effective at longer propagation time. Some real-life applications where solving KdV equation over short propagation involves earth quake modeling and simulation of optical laser pulses along fibers.

Solution of numerical experiment 2 using STADM
We note that the authors in [18] applied a numerical approach based on the homotopy perturbation transform method (HPTM) to derive the exact and approximate solutions of nonlinear fifth-order KdV equations for the study of magneto-acoustic waves in plasma. We would like to point out that equation (4) is considered in the paper [18,Example 3.4], where there are some typos in the paper and we have considered a modified and a corrected version of the nonhomogeneous equation, just for clarity to the readers.
By considering the nonhomogeneous equation in equation (4) via STADM and applying Shehu transform on both sides of equation (4), we get Taking inverse Shehu transform on both sides of equation (33), we obtain By applying the aforesaid decomposition method, we have , , , .
On comparing both sides of equation (35), we obtain 0  cos  sin  2  sin 2  , , , For some non-homogeneous PDEs, Wazwaz [29,30] suggested that the construction of zeroth component of the series can be defined in a slightly systematic way by supposing that if the zeroth component is , then the function G is possibly divided into two parts such as G x t , 0 ( ) and G x t , 1 ( ). Then one can construct the recursive algorithm u x t , 0 ( ) in an elegant way.
Comparison of modified ADM and classical FDM  387 By applying the same idea to equation (36), we obtain the modified recursive scheme as are obtained using the formulae (cf. [24,25]) The polynomials u u u , , The first-order approximation u x t , 1 ( ) is given by Thus, one can see from equations (37) and (41) that u x t , 0 n 1 ( ) = + for n 1 ≥ . The sum of first two approximations to build an approximate-analytical solution for u x t , ( ) of equation (4) is given by which immediately coincides with exact solution t x cos( ). Hence, it is important to note from above that this kind of modification of the recursive scheme not only avoids unnecessary computational investment such as computation of the Adomian polynomials but also it reduces large volume of calculations considerably.

Solution of numerical experiment 3 using STADM scheme
Consider the nonlinear KdV equation where x t , 0,2 0,2 , Equation (43) is singularly perturbed KdV-type equation if δ is infinitesimal (see also [31] and for singularly perturbed PDEs [32]). The case when α 1 = and δ 0.022 2 = was considered in [1]. We solve the experiment using STADM. By applying Shehu transform on both sides of equation (43), we obtain Taking inverse Shehu transform on both sides of equation (44) gives cos .
By employing the aforesaid decomposition method, we have where uu u u u , , , are Adomian polynomials that represent the nonlinear terms, and by comparing both sides of equation (46), we obtain the first few components: where δ 0.022 2 = . Thus, the first few sum of terms for approximate-analytical solution for u x t , ( ) is given by Remark 4.2. Figure 3 shows that the numerical results are much affected by the number of terms used.
Firstly, it appears that the first-order sum u x t u x t , , gives more realistic results. KdV equation is a conservative equation and we expect values of solution to be bounded as we progress in time. Secondly, it appears that scheme is unstable as we increase time to

Note on the convergence analysis of STADM
Numerous research works on the applications of ADM to the problems arising from various areas of pure and applied sciences are detailed in [33][34][35]. We note here that the theoretical convergence of ADM applied to the general KdV-type equations is proved in [36] and the general framework of providing the convergence of STADM is exactly the same and we omit the details as the convergence proofs are already shown in [36,37]. We just restate a sufficient condition of convergence for the KdV equation as given in [6] and see also [37] for the proofs.
with specified initial and boundary conditions, converges toward a particular solution.
Proof. We refer to [37, Section 3] for the proof (see also [36]). □ In the next section, we make use of a classical FDM to solve the three numerical experiments described in Section 2.

Classical FDM
The authors in [18] used a numerical approach based on the homotopy perturbation method to solve approximate solutions of some fifth-order KdV equation. See also the works of Wazwaz in [38]. In this section, we use the classical FDM to solve third-order and fifth-order nonlinear KdV equations considered in this paper and we also compare results with those of the STADM method.
We use the following difference approximations to approximate derivatives: where h x Δ = and k t Δ = are the spatial and temporal step sizes, respectively, and x ih i = , t nk n = , i 0, 1, = … and n 0, 1, = …, where superscript n denotes a quantity associated with time level t n and subscript i denotes a quantity associated with space mesh point x i .

Solution of numerical experiment 1 using finite difference scheme
The scheme is given by  = from the second expression which gives the maximum value for the inequality.
On substituting this into the inequality, we obtain the region of stability as [41] λ u ω ω ω u u = . In particular, for numerical experiment 1, since u π 2 6.28 , we obtain the stability region as t 0 Δ 0.158 < ≤ and we have used t Δ 0.001 = for this numerical experiment.

Consistency of the numerical scheme in equation (53)
We consider equation ( The scheme is consistent with the PDE given in equation (3) and is second-order accurate in time and space.

Numerical results for experiment 1
We obtain plots of numerical and exact profiles vs x in Figure 6 and corresponding plots of absolute errors vs x are shown in Figure 7. Table 3 displays absolute and relative errors at the spatial nodes for the four different values of time, which indicate that the scheme is very efficient at short, medium, and long time propagation.

Rate of convergence for numerical experiment 1
The rate of convergence R T as given in [40] is obtained by where e k and e ℓ are the error norms (w.r.t. L ∞ as in Table 4) corresponding to the spatial (or temporal) steps h k and h ℓ , respectively. Table 4 indicates that the numerical rate of convergence in time is two, which is the theoretical rate of convergence.
A single expression for the scheme is Region of stability for equation (4) using the classical scheme can be obtained by considering equation (4) with inhomogeneous term zero as the inhomogeneous term does not depend on u; that is, we consider A condition for stability is obtained by finding a condition for t x Δ , Δ so that for all ω π π , Since the second expression in the bracket for the above inequality dominates the first for small values of h, we obtain ω π 0.7325 rad = from the second expression which gives the maximum value for the inequality. On substituting this into the inequality, the region of stability is given as = . In particular, for numerical experiment 2, since u , we obtain the stability region as t 0 Δ 1.18254 10 3 < ≤ × − and we have used t Δ 0.0001 = for this numerical experiment.

Consistency of the scheme in equation (64)
We apply Taylor series expansion about t x , , Thus, the scheme in equation (66) is consistent with the PDE given in equation (4) and is second-order accurate in time and space as shown in equation (74) ( Table 5).
The following results for numerical experiment 2 show the method is effective for long propagation times when t 4.0 = and t 5.0 = as depicted in Figures 10 and 11. = . The classical finite difference scheme performs quite well for short, medium, and long time propagation.
and considered two cases: The classical finite difference scheme, using α 1 = , is given by    and the solitary waves retain their identity when they meet each other, apart from a phase shift. When two solitons of different amplitude meet, the smaller one is negative shifted, the bigger one positive shifted. These experimental observations are validated by the methods in this work and coincide with the FPU experiment in [1]. From Figure 13, we depict that there is a wiggling behavior around x 1 = for some time and also solitary wave nature appears as time fluctuates.
The recurrence time T R π 30.4 ( ) = as in [1], which is the time that all the solitons "focus" at a common spatial point. Such a recurring observation is not clear from the above plots and this may be due to the accumulative error of the numerical computation. However, later we see that at t  (80), we obtain the stability region t 0 Δ 1.155 10 , 2 < < × − and we have used t Δ 0.0001 = for this numerical experiment. The cause of the minor oscillations in Figure 13 could be due to numerical dispersion, which may lead to some noise-like structure in the numerical solution. We note that, for Case 2, it is difficult to obtain a numerical solution of the dispersive shock wave of the KdV equation [42] due to the weak dispersion limit (δ 10 4 = − ) to counterbalance the nonlinearity. Thus, as future work, we will construct some other methods such as NSFD and symplectic methods to solve numerical experiment 3 and compare the results especially in regard to the nonphysical oscillations. It is worth mentioning that for Case 1 as shown in Figure 12, the obtained numerical investigation is in good agreement with Zabusky-Kruskal's experimental findings [1]. We should point out here that the numerical simulations were done on a dell computer, with a windows 8 operating system (Intel core I5 with processor speed 1.8 GHz) ( Table 6).
Remark 5.2. We note here from equation (75) that (i) when δ 0 → , it behaves like Burger's equation; (ii) when δ 1 = , it is the classical third-order KdV equation; (iii) when δ 1 ≪ , it behaves like singularly perturbed KdV equation, which is studied when α 1 = , δ 0.022 2 = using t Δ 0.000002 = and x Δ 0.01 = (cf. [1]). We here point out that the authors in [43] proposed some composite schemes to reduce the dispersive effect on the numerical solution for some dissipative and dispersive phenomena [44].
We roughly observe from Figures 12 and 13 that as 0 δ α → (i.e., for δ sufficiently small), the numerical simulations using the Zabusky-Kruskal scheme exhibit shock-like structure as expected [1]. In other words, when the coefficient of δ gets smaller, the emergence of more oscillatory profiles become realistic.

Conclusion
In this paper, we constructed two reliable methods, namely, STADM and classical FDM to solve third-order and fifth-order KdV equations described by three numerical experiments. The first numerical experiment consisted of a homogeneous third-order KdV equation. STADM is very effective at small propagation time t 0.01, 0.05 ( ) = giving relative error of order 10 7 − , 10 3 − but less effective at medium propagation time t 0.10, 0.125 ( ) = giving relative error of 10 2 − , 10 1 − . The classical FDM is very efficient at small, medium, and long propagation times (t 0.01, 0.05, 0.10 = , and 1.0), with relative error of order 10 7 − to 10 4 − . The numerical rate of convergence in time for the classical scheme is obtained and it coincides with second-order theoretical rate of convergence in time.
The second numerical experiment consisted of nonhomogeneous nonlinear dispersive KdV equation. For some nonhomogeneous PDEs, it is possible to construct zeroth component in a systematic way for semianalytic methods. This was done and STADM gave exact results. We note that exact results are not always guaranteed as this depends on the initial conditions. The classical finite difference performs well at t 0.05, 0.10, 0.50, 1.0 = with relative error of order 10 5 − to 10 3 − and even at long propagation times of t 4.0 = and 5.0. The third numerical experiment consisted of homogeneous third-order KdV equation but no exact solution is known. STADM gives very different results when using first-order sum and second-order sum. The classical scheme gives realistic results and does not experience blow up at long propagation times. The KdV equation is conservative and therefore the solution values must always be bounded and blow up phenomenon is not expected.
The study is important in many ways. It allows us to conclude that for the first and third numerical experiment that the classical scheme is better one. STADM is very efficient for first numerical experiment at short time propagation only and it gives exact results for the second numerical experiment. The study gives insight in constructing other PDEs such as time-space fractional KdV, KdV-Burgers, stochastic KdV which model many real-life phenomena such as wave dynamics, geophysical flows etc.