Runge-Kutta Integration of the Equal Width Wave Equation Using the Method of Lines

The equal width (EW) equation governs nonlinear wave phenomena like waves in shallow water. Numerical solution of the (EW) equation is obtained by using the method of lines (MOL) based on Runge-Kutta integration. Using von Neumann stability analysis, the scheme is found to be unconditionally stable. Solitary wave motion and interaction of two solitary waves are studied using the proposed method. The three invariants of the motion are evaluated to determine the conservation properties of the generated scheme. Accuracy of the proposed method is discussed by computing the L 2 and L ∞ error norms. The results are found in good agreement with exact solution.


Introduction
Nonlinear equations are widely used to describe complex phenomena in various fields of science, such as fluid mechanics, plasma physics, solid-state physics, elastodynamic problems, nonlinear Schrodinger equations, elasticity problems, free vibration analysis, and optical fibers [1][2][3][4][5][6][7][8][9][10][11][12][13][14][15][16].Benjamin et al. [17] advocated that the partial differential equation (PDE) modeled the same physical phenomena equally well as the (KdV) equation given the same assumptions and approximations that were originally used by Korteweg and de Vries [18].This PDE of Benjamin et al. [17] is now often called the BBM equation, although it is also known as the regularized long wave (RLW) equation.Morrison et al. [19] proposed the onedimensional PDE as an equally valid and accurate model for the same wave phenomena simulated by the KdV and RLW equations.This PDE is called the equal width (EW) equation because the solutions for solitary waves with a permanent form and speed, for a given value of the parameter , are waves with an equal width or wavelength for all wave amplitudes.The equal width (EW) equation, which was introduced by Morrison et al. [19], is an important special kind of a nonlinear dispersive wave equation.It is defined as For a smooth function  = (, ) on a domain Ω × [0, ] with Ω ⊆ R. Except for a single travelling solitary wave solution, no analytic solutions are known, and therefore numerical methods have to be used.The equation has solution of the form That represents a solution for a solitary wave traveling at constant speed  and vanishing at ∓∞.When the EW equation is used to model waves generated in shallow water channel, the variables are normalized so that the distance  and water elevation  are scaled to the water depth ℎ, and time  is scaled to √ℎ/, where  is the acceleration due to gravity.There is experimental evidence to suggest that both descriptions break down if the amplitude of any wave exceeds 0.28 [20].Different numerical solution methods were applied to solve EW equation such as L. R. Gardner and G. A. Gardner [21] using Galerkin's method-based cubic B-spline finite elements, Ramos [22] using finite difference methods, Raslan [23] using quartic B-spline method, Saka [24] using finite element method, Ali et al. [25] using He's Exp-function method, and Cheng and Liew [26] using the improved element-free Galerkin's method.
In this paper, the method of lines (MOL) solution of the EW equation is presented.This method consists of converting the EW equation with auxiliary conditions into a system of ordinary differential equations with corresponding auxiliary conditions.Then the fourth-order Runge-Kutta method is used to solve the system of first-order ordinary differential equations instead of finite difference methods [22] that are accurate and efficient as shown in [27].

The Method of Lines Solution of the EW Equation
The method of lines is a general technique for solving partial differential equations (PDEs) by typically using finite difference relationships for the spatial derivatives and ordinary differential equations for the time derivative.For timedependent partial differential equations, meshless kernelbased methods were similarly based on a fixed spatial interpolation, but now the coefficients are time-dependent, and one obtains a system of ordinary differential equations for these.This is the well-known method of lines, and it turned to be experimentally useful in various cases.
To apply the method of lines for solving the EW equation, firstly we subdivide the solution domain into uniform rectangular meshes by the lines In numerical calculation, Ω = R is usually replaced by Ω = [, ] and with zero boundary condition at both ends.Hence, the solution domain of the EW equation is the rectangle defined as  ≤  ≤ , 0 ≤  ≤ .Then, the partial derivatives depending on spatial variables, / and  2 / 2 , in (1) are replaced by the well-known finite difference approximation at point   : By substituting (4) into (1) and introducing the boundary conditions  0 () =  +1 () = 0, this yields a system of ordinary differential equations which depend on  in the following form: where Thus, we have the system of differential equations of one independent variable .This system can be solved using Runge-Kutta method.The resulting system of ordinary differential equations is integrated with respect to time.Many existing numerical schemes for evolutionary problems in partial differential equations (PDEs) can be viewed as method of lines (MOL) schemes.In 1984, Verwer and Sanz-Serna [28] treated the convergence of one-step MOL schemes.
Their main purpose was to set up a general framework for a convergence analysis applicable to nonlinear problems.

Stability Analysis
By solving EW equation using the method of lines, firstly the spatial variables are discretized to obtain an ordinary differential equations system in the time variable that can be discretized by an ordinary differential equation (ODE) solver.The simplest ODE solver is the forward Euler method and it is used widely for analysis of the stability properties of the spatial discretization.However, while forward Euler is ideal for analysis of the stability properties of a given spatial discretization, it is only first-order accurate.In practice, high order time discretization that preserves all the stability properties of forward Euler is needed.In [29], high order strong stability preserving time discretization method for the semidiscrete method of lines approximations of PDE is developed.These methods are derived by assuming that the first-order forward Euler time discretization of the method of lines ODE is strongly stable under a certain norm when the time step Δ is suitably restricted, and then try to find a higher order time discretization (Runge-Kutta or multistep) that maintains strong stability for the same norm perhaps under a different time-step restriction.
The stability analysis constitutes the essential study of the numerical solution of PDEs.In general, this is because such study provides the means by which the step size and the numerical integration scheme for the given differential equation could be selected so as to secure manageable numerical solution.
The stability analysis is based on the Neumann theory in which the growth factor of the error in a typical mode of amplitude ê is where  is a Fourier number and ℎ, the finite difference size, is determined from a linearization of the numerical scheme.
Assuming  in the nonlinear term as a constant  = max    ; this enables one to discuss the stability in the linearized sense.
The numerical method of lines of EW equation gives the system of ordinary differential equations A trial solution is assumed and substituted into (9).The trial solution must take into account the variation of (, ) with both  and  or  and .So In accordance with a method proposed by von Neumann, the function () can be of the following form: Substituting ( 10) and ( 11) into (9) gives Equation (12) shows the growth factor for the error of the following form: In [30], Evans and Raslan proved that || ≤ 1 and the scheme is virtually unconditionally stable, at least for any practical problem, where ℎ is a small quantity and  represents the single speed and will usually be around unity.Kreiss and Scherer in [31] derived the conditions of local stability of Runge-Kutta methods when applied to hyperbolic partial differential equations, as when the time in ( 12) is discretized by using a locally stable Runge-Kutta method, the resulting completely discretized method is stable provided that || ≤  with locally stable Runge-Kutta methods whose stability region contains a half circle.Consider where  denotes the time step and for computational purposes the Runge-Kutta method is only useful if it is stable for sufficient small /ℎ.

Test Problems
In this section, the results of the numerical solution of the EW equation based on the MOL using Runge-Kutta time discretization are presented.The numerical simulation includes the propagation of a solitary wave, the interaction of two solitary waves, the development of an undular bore, and the temporal evaluation of a Maxwellian initial pulse.
In order to confirm the accuracy and efficiency of the method,  2 and  ∞ error norms are used and defined by where  denotes the exact solution and  denotes the numerical solution.According to [32], in many cases, there are the following time invariants: By sufficiently fine spatial resolution, there is no problem to maintain these invariants to reasonable accuracy.Wherever possible, this statement is supported by providing numerical results.The constants of the motion, for a solitary wave of  amplitude 3 and depending on  as given by ( 2), may be evaluated analytically to give 4.1.Single Solitary Wave.We first study the motion of a single solitary wave.This is derived from the initial condition.Consider where  = 1/(2 √ ) and  is a constant.This follows from the analytic solution (2) of the EW equation.We choose  = 1.0,  = 1.0,  0 = 15,  = 1, ℎ = 0.2, and  = 0.2 through the interval [0, 80].In this case the problem reduces to To apply the method, let   = 0.2 ( = 0, 1, 2, . . ., 400) By substituting (4) into (19) and introducing the boundary conditions  0 () =  +1 () = 0, this yields a system of ordinary differential equations depending on  in the following form: where  −1 is the inverse of the matrix  and This system can be solved using Runge-Kutta method.The solution for times from  = 0.0 and 40.0 is shown in Figure 1.
In Table 1, we examine various space/time step combinations and compare with simulation results given in [33].The most accurate simulation is obtained in the case of ℎ =  = 0.04 for which  2 error has a value of 7.5 × 10 −3 .
The simulation is run to time  = 40 and the three invariants  1 ,  2 , and  3 whose analytical values can be obtained as  1 = 12.0,  2 = 28.8, and  3 = 57.6 are listed for the duration of the simulation.We found that the MOL with Runge-Kutta integration is more accurate than Galerkin's method with linear element in studying solitary wave of amplitude 3.0 for ℎ =  = 0.04.
A solitary wave of amplitude 0.3 has also been modeled, and the results of the simulation are given in Table 2.By the time  = 80,  2 and  ∞ error norms remain less than 1.8 × 10 −3 .The profiles of the solitary waves at time  = 0 and  = 80 are compared in Figure 2.
In further simulation of solitary waves of smaller magnitude, (2) is taken as initial condition with  0 = 10.0,  = 0.03, and  = 0.01 that the solitary waves have amplitudes 0.09 and 0.03.The simulations are run to time  = 80 and  2 and  ∞ error norms and the invariants  1 ,  2 , and  3 are recorded throughout the simulation; see Tables 3 and 4.

Two Solitary Waves.
As a second test problem for the EW equation, we chose the interaction of solitary waves, as mentioned in the used literature.The initial function is where and solved the EW equation over the region 0 ≤  ≤ 80 considering  1 = 1.5 and  2 = 0.75,  1 = 10,  2 = 25, ℎ = 0.1, and  = 0.1.The shape of two solitary waves with given parameters is illustrated at the times  = 0.0 and 30 in Figure 3.As seen in Figure 3, two solitary waves at the time  = 0.0 are propagated to the right with velocities dependent upon their magnitudes and reached a stage where the larger solitary wave has passed through the smaller solitary wave and emerged in their original position.The values of the three invariants obtained by the present methods for this numerical experiment are recorded in Table 5.They are found to be very close to values given in [32] and the analytical ones given by  1 = 12( 1 +  2 ) = 27,  2 = 28.8( 2 1 +  2 2 ) = 81, and  3 = 57.6( 3  1 +  3 2 ) = 218.7.Additionally, we have observed that  3 remains almost constant, while  1 and  2 are affected more from the interaction of the two solitary waves during run of the algorithm.

4.3.
The Undular Bore.When a deeper stream of water flows into an area of the still water in a long horizontal channel, a bore is formed.To study development of an undular bore followed earlier by Peregrine [20], the following initial condition is used: where (, 0) denotes the evaluation of the water level above the equilibrium surface at time  = 0.The change in the water level of magnitude  0 is centered on  =  0 and  measures the steepness of the change.We insistently choose the parameters  = 5,  = 0.16666667,  0 = 0.1, and  0 = 0.0 in the region −20 ≤  ≤ 50 to make a comparison with earlier works.The simulation is run until time  = 800 and the values of the quantities  1 ,  2 , and  3 are recorded in Table 6 with the values given in [34].Figure 4 shows the undular bore profiles at times  = 0,  = 200,  = 400,  = 600, and  = 800 for the gentle slope  = 5.For the steep slope  = 2, the results are compared with those given in [34] in Table 7.

Conclusion
The EW equation is numerically solved using MOL with Runge-Kutta integration.The numerical solution leads to an unconditionally stable algorithm with which accurate simulations of the motion of a solitary wave are found over an extended time scale.Single solitary waves of amplitude 3.0 for ℎ =  = 0.04 are found.MOL with Runge-Kutta integration is more accurate than Galerkin's method with linear elements for a single solitary wave; it gives better conservation than that given in [33].Moreover, for a smaller solitary wave amplitude = 0.03, excellent results are obtained.The performance of the method has been examined by studying the propagation of interaction of two solitary waves.Using the proposed method, the magnitude, profile, and position of solitary waves are faithfully retrieved.To make a comparison between results of MOL and some earlier published works, the same step ℎ is used.
The MOL based on Runge-Kutta integration was efficiently applied to an undular bore.The obtained results are good and accurate as compared to the results in [34].In addition, the Maxwellian initial condition is simulated.
The obtained results indicate that the present method is remarkably successful numerical technique for solving the EW equation.The results suggest that proposed method, whose application is easier than many other numerical schemed methods such as finite element methods and Galerkin's method, can be applied to this type of the nonlinear problems with success.The method can be also used efficiently for solving a large number of physically important nonlinear problems.

Table 5 :
Invariants for the interaction of two solitary waves.

Table 8 :
Invariant for the Maxwellian initial condition.