Some Remarks on the Stability Condition of Numerical Scheme of the KdV-type Equation

This paper has employed a comparative study between the numerical scheme and stability condition. Numerical calculations are carried out based on three different numerical schemes, namely the central finite difference, fourier leap-frog, and fourier spectral RK4 schemes. Stability criteria for different numerical schemes are developed for the KdV equation, and numerical examples are put to test to illustrate the accuracy and stability between the solution profile and numerical scheme. 2000 Mathematics Subject Classification: 35Q53, 65M60, 37K10, 37K05


Introduction
The Korteweg-de Vries (KdV) equation is a nonlinear, dispersive partial differential equation for a function where u (x, t) of two real variables, space x and time t.It is a mathematical model of waves on shallow water surfaces and particularly notable as the prototypical example of an exactly solvable model, i.e., a non-linear partial differential equation whose solutions can be exactly and precisely specified.The solitary solution of the equation was first observed in J. Russell (1837), and the equation itself was later derived by Korteweg and de Vries in D. Korteweg & G. de Vries (1895).Since then, it has been applied in many fields to describe a wide range of physical phenomena such as interaction of nonlinear waves M. J. Ablowitz & D.E. Baldwin (2012), collision-free hydro-magnetic waves in a cold plasma, ion-acoustic waves, interfacial electrohydrodynamics (M.Q. Tran, 1979).
The KdV equation was not studied much until Zabusky and Kruskal (1965) proposed an explicit numerical scheme and discovered numerically that its solutions seemed to decompose at large times into a collection of "solitons", a series of well separated solitary waves.The scheme is described as follows ) , with j = 1, 2, . . . .Here central difference approximations were used for both the first space and first time derivatives to improve the accuracy for given step sizes ∆x and ∆t, respectively.
The study of non-linear waves would not have been so successful had it not done with stable numerical schemes, especially for time-dependent problems, stability guarantees that the numerical method produces a bounded solution whenever the solution of the exact differential equation is bounded.Stability, in general, can be difficult to investigate, especially when the equation under consideration is nonlinear.In this article, we analyze the stability of three numerical schemes on the KdV equation based on von Neumann method and conclude that Fourier RK4 scheme can meet the stability criterion with suitable spatial and time steps.

Finite Difference Scheme
In this section we present a finite difference scheme for the linearized KdV equation in order to proceed with the stability analysis.The equation is described as We note that the KdV equation ( 1) is used to be considered as the weakly nonlinear, weakly dispersive behavior of the long wave case on the free surface.However, equation ( 2) not only can be solved explicitly using Fourier methods but also served as a tool for studying the nonlinear equation (1).
We consider a function u (x i , t n ) with the x-t plane subdivided into a rectangular grid or mesh with each rectangle having sides of length h and k, where x i = ih, t n = nk with i = 1, 2, . . ., m and n = 1, 2, . . ., q, for some integers m and q.The various mesh points may be labeled by a pair of integers.Let the point P n i be denoted by (i, n).The value of u at P n i is approximated by u n i .Expressing the finite-difference approximation in terms of this notation, we have ( ∂u ∂x ) To analyse the numerical scheme for (2), we adopt the central finite difference (FD) scheme with leap-frog time-stepping formula as follows The characteristic equation for this recurrence relation is which we obtain by inserting in (6) the ansatz u n i = g n , and the condition for stability is that both complex roots must lie in the closed unit disk, with only simple roots permitted on the unit circle.As a result, we obtain When x = 2π/3, we have Therefore we have obtained which is the stability criteion for linearized KdV equation ( 2) and can be used to calculate the numerical solution for (1).

Fourier Leap-frog Scheme
In this section, we first model the KdV equation ( 1) by a Fourier spectral method on [−L, L].In practice, it can be transformed to where L is a given number representing the boundary point of the spatial domain.In practice, we need to transform u, u x into Fourier space and discretize the equation.For any integer N > 0, we consider the collation points x j = j∆x = 2π j/N, j = 0, 1, . . ., N − 1, and note that if u (x, t) is the solution of the KdV equation, then we transform it into the discrete Fourier space as From this, using the inversion formula, we get where we denote the discrete Fourier transform and the inverse Fourier transform by F and F −1 respectively.Therefore, we have in particular, we have ) Then ( 12) can be transformed into a semi-discrete form as follows: Taking into the consideration of the collation points, ( 12) can be further discretized into 19) can be written in the vector form as where F defines the right-hand side of ( 19).
For the stability analysis, it requires further information on the condition imposed on the time step ∆t.Therefore, similar to (12), we could start with the following linearized equation for some suitable boundaries L, and approximate the solution by using Fourier Leap-Frog (FLF) transforms where F is called the fourier transform operator.By using the discrete version of ( 13) and ( 14), the scheme is shown as follows For simplicity, let v = k (2π/L), then (23) can be written as To proceed, we look for a solution to (24) of the form and subsitute it into (23) to get κ (t+∆t)/∆t e ikx − κ (t−∆t)/∆t e ikx + 2iv∆tκ t/∆t e ikx − 2iv 3 ∆tκ t/∆t e ikx = 0, ( i.e., where The scheme is conditionally stable if and only if f (v, ∆t) is real and less than one in magnitude.Since the wave number v takes the values v = 0, ±1, ±2, . . ., ±N/2, (30) and the interval is discretized with N equidistant mesh points, that is we want to find the largest value of ∆t such that is true for all v.Here the most severe restriction on ∆t is imposed for the v, which are largest in magnitude, i.e., for v = ±v max with v max = N/2 = π/∆x.Thus, we obtain and the stability condition becomes

Fourier RK4 Scheme
We have implemented the Fourier spectral method with fourth-order Runge-Kutta (RK4) time differencing to solve the fifth-order KdV equation.The RK4 method is known to have a truncation error of O ( ∆t 4 ) and one of the most widely used methods for solving differential equations.Its algorithm is described below: where For the study of the stability, we use the standard Fourier analysis to find the condition imposed on the time step ∆t.For simplicity, we consider the linearized KdV equation as follows: We could approximate this equation by using RK4 scheme in (36) with the ansatz u (x, t) = κ t/∆t e ivx , which will not only extend the domain of the equation to the whole real line but also enable us to examine max |κ| to decide on stability of the numerical scheme.As a restult, we have Therefore (36) gives We want to find the restriction on ∆t such that max |κ| < 1, ( 47) is true with v max = N/2 = π/∆x.Detailed calculation for (47) gives the condition The analysis of the present article should be more broadly in line with the KdV equation ( 1).Briefly, different numbers of α in (48) will result in different stability regions, hence we assume that α = 1 in (48) to get which is the stability criterion for ∆x and ∆t.

Numerical Test
One of the most interesting features of the KdV equation is the existence of infinitely many conservation laws.A conservation law in differential equation form can be written as T t + X x = 0, in which the "density" T and the "flux" X are polynomials in the solution u and its x-derivatives (P.G. Drazin & R. S. Johnson, 1989).If both T and X x are integrable over the domain (−∞, +∞) , then the assumption that X −→ 0 as x −→ ∞ implies that the conservation law can be integrated over all x to yield d dt where C is a constant.The integral of T over the entire spatial domain is therefore invariant with time and usually called an invariant of motion or a constant of motion (Verheest and Hereman 1994; Goktas and Hereman 1999).The KdV equation (1) itself is already in conservation form, i.e., traveling to the right.In addition, the stability requirement in (35) gives the value ∆t/∆x 3 ≈ 0.00676 < 0.0322 to further support the results and performance of numerical calculations.
On the third case, we discuss the problem when the KdV equation is solved numerically by spectral methods, the pattern is usually the same: spectral differentiation in space, RK4 differences in time.Here we consider the same initial condition as described in (58) and compute the solutions with the numerical scheme in ( 23) with 256 mesh points in domain [−π, π], time step ∆t =1e-7 and plot the solutions from time t = 0 to t =6e-3.Solution profile is plotted at time t ∈[0, 6e-3] in Fig. 8.The results illustrate the phenomenon that the wave propagates at a constant speed coherently.Moreover, neither ripples nor dispersive wave components are emanating from the calculations which brings satisfactory results.On the other hand, Fig. 9 shows the results of the exact and numerical solitary wave solutions and differences in between these two waves.It is seen that the numerical error is maintained to be accurate to 10 −3 in the magnitude of 2000 of wave amplitude.Hence under the FRK4 scheme, the numerical solutions have been undergoing for 60,000 integrations and obtained without any significant expense of accuracy with regard to its height and speed.Furthermore, the velocities of the solutions and their shapes are almost unchanged compared to the exact solution.Moreover, accuracy test for the numerical scheme in Fig. 10 shows that the errors have been remarkably kept as small as O ) .This not only shows that the FRK4 scheme is significantly stable but also the accuracy of the scheme is guranteed to produce numerical calculations.Consequently, the criterion equation in ( 49) gives CRI = 0.99 < 1 and supports the fact that the computation through FRK4 becomes the choice of numerical methods.

Concluding Remarks
Fig. 1 illustrates the stability regions from ( 11), ( 35) and ( 49) for different schemes on the KdV equation with ∆t against ∆x.Three numerical schemes including the finite difference, fourier leap-frog and fourier RK4 procedures are presented for the equation.By using the analysis on the linearized equation, we deduce that the marginal curves of stability region are all nonlinear curves.It is known that the fourier spectral method consisting of space space differentiation and integrating procedure in time can be swiftly convergent and spectrally accurate in spatial domain.In particular, the region of stability of the fourth order Runge-Kutta method is so complicated that it can not be expressed in terms of a closed form of algebraic equation but can be characterized through the symbolic computation.The result also shows that the stability region for FLF scheme is stricter that that of the FRK4 scheme.
Stability criteria of numerical schemes on the KdV equation are developed, namely the central finite difference, fourier leap-frog, and fourier spectral RK4 schemes.Each of the schemes gives us an estimated region for allowed values of ∆x and ∆t.The errors measured by the spread of the values of conservation law from its initial value are recorded, which gives an accuracy test for the numerical schemes.This paper also carried out a comparative study of solitary wave solution for different schemes and numerical settings to test the accuracy and stability condition.It is seen that fourier spectral RK4 scheme strikes the right balance between the finite difference and fourier leap-frog schemes for numerical methods.