POLYNOMIAL AND RATIONAL WAVE SOLUTIONS OF KUDRYASHOV-SINELSHCHIKOV EQUATION AND NUMERICAL SIMULATIONS FOR ITS DYNAMIC MOTIONS

Polynomial and rational wave solutions of Kudryashov-Sinelshchikov equation and numerical simulations for its dynamic motions are investigated. Conservation flows of the dynamic motion are obtained utilizing multiplier approach. Using the unified method, a collection of exact solitary and soliton solutions of Kudryashov-Sinelshchikov equation is presented. Collocation finite element method based on quintic B-spline functions is implemented to the equation to evidence the accuracy of the proposed method by test problems. Stability analysis of the numerical scheme is studied by employing von Neumann theory. The obtained analytical and numerical results are in good agreement.


Introduction
Recently, the concept of soliton has become a widely-studied topic because of its applicability in the modeling of many natural phenomena arising in different fields such as optics, physics, dynamics, fluid and biology [1,2,5,7,13,23]. The theory of soliton has been put forward in the study of nonlinear phenomena, which makes the study of the nonlinear evolution equations (NLEEs) become an essential tool in nonlinear science.
In this paper, we investigate and analyze the solutions of the Kudryashov-Sinelshchikov equation [21,28] analytically and numerically by using the unified method (UM) [18,19] and collocation method, respectively.
The Kudryashov-Sinelshchikov equation, denoted by KSE, is given by where α, β, γ and σ are real parameters. The KSE can be used to describe the nonlinear waves in gas-liquid mixtures in the absence of the viscosity [21]. The KSE is the generalization of the KdV when γ = σ = 0. For more details see [10]. Our work aims firstly to search for the solitary, soliton, and soliton rational wave solutions for the KSE by making use of the UM. Secondly, the conservation law is applied to study the invariance and the multiplier approach for the KSE. Finally, numerical simulation of the obtained results are performed on the behavior of the obtained solutions to the KSE using the quintic B-splines and the modulation instability is discussed and introduced in some figures to clarify the physical meaning for explicit and approximate solutions.
This paper is organized as follows. Section 2 constructs different types of wave solution for the KSE using the UM. Section 3 recalls some necessary preliminaries of the conservation flows for the KSE. Section 4 presents collocation finite element method for the KSE by using quintic B-splines and studies stability of the scheme. To demonstrate accuracy of the scheme, some test problems are studied and the obtained results are discussed in Section 5. Finally, conclusions are given in Section 6.

Exact solutions
In this section, we use the UM [18,19] to derive exact wave solutions for Eq. (1.1). The UM classifies the solutions into two categories namely polynomial and rational wave solutions (for more details see [18,19]).
First, we use the transformation u(x, t) = u(ξ), ξ = κx + νt, (2.1) in Eq. (1.1) and integrating both sides with respect to ξ, we get where the constant of integration is set to be equal zero and κ and ν are arbitrary constants.

Polynomial wave solutions
The UM asserts that the polynomial wave solutions are given by and the auxiliary function Γ(ξ) satisfies the auxiliary equation where p j and b j are real constants. By considering the homogeneous balance between u ′′ and u 2 in Eq. (2.2), we get n = 2(k − 1), k ≥ 1. Here, we confine ourselves to find these solutions when k = 2 and p = 1 or p = 2 to obtain solitary and soliton wave solutions respectively. (2.7) Using Eq. (2.7) in Eq. (2.5) and by solving the auxiliary equation given by Eq.

Rational wave solutions
To find the rational wave solutions of Eq. (1.1), we assume the solution of Eq. (2.2) as where p i , q i and b i are arbitrary constants. By considering the homogeneous balance between U ′′ and U 2 in Eq. (2.2), we get n − m = 2 (k − 1) and k = 1, 2, ... . Here, we find those solutions when k = 1, n = m, and p = 1 (soliton rational solutions).

Conservation laws
We present some preliminaries on conservation laws that will be used in the analysis that follow. The invariance and multiplier approach are resorted based on the well known result that the Euler-Lagrange operator annihilates a total divergence to determine conserved densities and fluxes [8]. Firstly, if (T t , T x ) is a conserved vector corresponding to a conservation law, then along the solutions of the differential equation E(x, t, u, u x , u t , . . .) = 0. Moreover, if there exists a non-trivial differential function Q, called a "multiplier", such that then QE is a total divergence, Thus, a knowledge of each multiplier Q leads to a conserved vector determined by, inter alia, a homotopy operator [9]. and Also, the above is trivially extendable to the multi space situation so that, in the (1 + 2) scenario with independent variables (t, x, y), the conserved flow has the form (T t , T x , T y ) and the total divergence Note. In this case, using the the language of forms, the conserved form (as opposed to the conserved flow) would be It is well known that the vector fields that leave the system of differential equations invariant (generators of Lie point symmetries) contain the algebra of variational symmetries, if the latter exists [12,14,15,22].
The multiplier approach, for σ, γ ̸ = 0, leads to multipliers The corresponding conserved flow with Q 1 has components The conserved density with Q 2 is For σ = 0, the multiplier is Q 3 = ln (γu + β) with conserved density

Numerical simulation
In this section, firstly, quintic B-spline interpolation function and its some properties are defined. The collocation finite element method is implemented to Kudryashov-Sinelshchikov equation. Quintic B-spline functions are chosen as interpolation functions. Then, stability of the applied method is analyzed using von Neumann theory.

Quintic B-splines and properties
Here, Kudryashov-Sinelshchikov equation (1.1) is studied with the physical boundary conditions u → 0 as x → ±∞, where α, β, γ and σ are constants and the subscripts t and x denote the temporal and spatial differentiations, respectively. Firstly, the solution domain is limited over an interval a ≤ x ≤ b to apply the numerical method. Space interval [a, b] is separated into uniformly sized finite and the initial condition are chosen as above. ϕ m (x) quintic B-spline functions at the knots x m are defined over the interval [a, b] by the following relationships for m = −2(1)N + 2 [20]: 3) The values of ϕ m (x) and its derivative are tabulated as in Table 1. The splines ϕ m (x) and its four principle derivatives vanish outside the interval where δ j (t) is time dependent parameters to be determined from the boundary and collocation conditions. Using trial function (4.4) and quintic B-splines (4.3), the values of u, u ′ , u ′′ , u ′′′ and u iv at the knots are determined in terms of the element parameters δ m by where the symbols ′ , ′′ , ′′′ and iv denote first, second, third and fourth differentiation with respect to x, respectively.

Implementation of the method
Substituting Eq. (4.5) into Eq. (1.1), general form of the solution approach is obtained as below. and and · indicates derivative with respect to t. Here, the term u in non-linear terms uu x and uu xxx , the term h 2 20 u xx in non-linear term u x u xx are taken as Eqs. (4.7) and (4.8) by assuming that the quantity u and h 2 20 u xx are locally constants for the linearization technique.
If time parameters δ i 's and its time derivativesδ i 's in Eq. (4.6) are discretized by using the Crank-Nicolson formula and usual finite difference approximation, respectively: A recurrence relationship between two time levels n and n + 1 relating two unknown parameters δ n+1 (4.10) where (4.11) The system (4.10) includes (N + 1) linear equations and (N + 5) unknown param- To obtain a unique solution from this system, four additional constraints are required. These are obtained from the boundary conditions by eliminating δ −2 , δ −1 and δ N +1 , δ N +2 in the system (4.10). In this case, a matrix equation is obtained (N + 1) linear equations and N + 1 unknowns d = (δ 0 , δ 1 , . . . , δ N ) T in the following form P d n+1 = Qd n . (4.12) The matrices P and Q are (N + 1) × (N + 1) pentagonal matrices and this matrix equation (4.12) is solved by using the pentagonal algorithm. To cope with the nonlinearity caused by ω m and ϖ m , two inner iterations are applied to the term δ n * = δ n + 1 2 (δ n −δ n−1 ) at each time step. After the initial vector d 0 = (δ 0 , δ 1 , ..., δ N −1 , δ N ) is determined by using the initial condition and the following derivatives at the boundary conditions, the following matrix form of the initial vector d 0 is obtained Once and for all, this matrix system is solved efficiently by using a variant of Thomas algorithm.

Stability of the scheme
The stability of the linearized numerical scheme is analyzed based on the von Neumann theory. To investigate stability, we use Fourier mode where ξ is growth factor of the error in a typical mode of amplitude ξ n , h is element size and k is mode number. Substituting the Fourier mode (4.16) into the system (4.10) gives the following equality

Results and discussion
In this section, numerical solutions of the Kudryashov-Sinelshchikov equation are considered for three problems which are including the motion of single solitary wave, interaction of two solitary waves, evolution of solitary waves with Gaussian and undular bore initial conditions. To calculate the difference between analytical and numerical solutions and to show accuracy of the applied numerical scheme, the error norm L 2 (5.1) and the error norm L ∞ are used. The Kudryashov-Sinelshchikov equation (1.1) possesses two conserved quantities as follows: 3) The conserved quantities C 1 and C 2 are calculated to check the conservation of some quantities during the wave motion.  Table 2. It can be seen from the Table 2 that L 2 and L ∞ error norms are found to be satisfactorily small and the conserved quantities are pretty much unchanged as the time progresses. The percentage of relative changes of C 1 and C 2 are found to be 5.920 × 10 −5 and 6.488 × 10 −6 for h = ∆t = 0.1; 4.689 × 10 −5 and 4.720 × 10 −6 for h = ∆t = 0.025, respectively. It is observed from the Table 2 that numerical results are better and more accurate for small values of h and ∆t. If we observe Figure 4(a), we find bell shaped solitary wave solutions produced. Also, contour line for motion of the single solitary wave can be seen in Figure 4(b). As it can be clearly seen from the figures, the solitary wave moves to the right at a constant speed. As expected, its amplitude and shape are preserved as time passes. On the other hand, error graphs are plotted for different values of h and ∆t in Figure 5.

Interaction of two solitary waves
Secondly, the interaction of two solitary waves advancing in the same direction is considered by using the initial condition given by the linear sum of two well separated solitary waves having different amplitudes  Table 3 for different values of h and ∆t. It is seen that the obtained values of the conserved quantities remain more constant sensibly for small values of h and ∆t. While interacting, the travelling wave profiles can be seen in Figure 6(a). Contour line for interaction of two solitary waves is depicted in Figure 6

Evolution of solitary waves
As a last test problem, Gaussian and undular bore initial conditions are studied to show birth of solitary waves.

Gaussian initial condition
Evolution of a train of solitary waves is studied on Kudryashov-Sinelshchikov equation using the Gaussian initial condition  Table 4. Here, especially the conserved quantity C 2 is preserved better for small values of h and ∆t. Also, Figure 7(a) illustrates the development of the Gaussian initial condition into solitary waves. As it is seen from the Figure 7(a), a solitary wave plus and oscillating tail are drawn. As seen from the figures, two solitary waves moving is observed. Contour line for evolution of solitary waves with Gaussian initial condition is depicted in Figure 7(b).

Undular bore initial condition
Production of a train of solitary waves is studied on Kudryashov-Sinelshchikov equation using the undular bore initial condition    Table 5. Especially, the conserved quantity C 1 is preserved better for small values of h and ∆t. Undulation bore keeps steady-state during running the which can be observed at specified times in Figure 8(a) and contour line for evolution of solitary waves with undular bore initial condition is shown in Figure 8(b). As it is seen from these figures, the initial perturbation evolves into a good developed train of solitary waves. As the time progresses, five solitary waves moving is observed.

Conclusion
Polynomial and rational wave solutions of Kudryashov-Sinelshchikov equation and numerical simulations for its dynamic motions have been studied. Conservation quantities of the dynamic motion have been calculated using multiplier approach. A collection of exact solitary and soliton solutions of Kudryashov-Sinelshchikov equation has been derived. Collocation finite element method based on quintic Bspline functions has been applied to the equation. The numerical scheme has been (a) Generated waves (b) Contour for generated waves Figure 8. Generated waves and its contour for undular bore initial condition shown to be unconditionally stable. To demonstrate the accuracy of the proposed method, we have considered three test problems including the motion of single solitary wave, interaction of two solitary waves and evolution of waves with Gaussian and undular bore initial conditions. Then, analytical and numerical results have been compared by aid of error norms. The obtained analytical and numerical results are in good agreement. The results of this paper come with a lot of encouragement for further future investigations.