Analytical and numerical solutions for the variant Boussinseq equations

The principal objective of this paper is to employ the -expansion and the adaptive moving mesh methods to express the exact travelling wave solutions and the numerical solutions, respectively, for the variant Boussinseq equations. Hyperbolic tangent and cotangent functions are utilized to build the exact solutions. The used numerical method uses the finite differences to discretize the proposed equations. The obtained numerical results are compared with other results obtained using the uniform mesh method. The achieved results show that both solutions match with each other. We also illustrate some 2D and 3D figures to confirm the validity of the numerical approach applied here.


Introduction
Various complex physical phenomena such as chemical kinematics, plasma physics, chemical physics and optics are investigated using nonlinear evolution equations (NLEEs). Therefore, effective and efficient approaches to construct the travelling wave solutions have attracted a diverse group of experts. Some scientists have developed an extensive variety of techniques to obtain the exact solutions of partial differential equations (PDEs). Some of these methods are the trial function method [1], the inverse scattering transform [2], the Weierstrass elliptic function approach [3], the sine-cosine method [4,5], the F-expansion method [6], the tanh-sech approach [7,8], the modified tanh-function technique [9], Hirota's bilinear method [10,11], the extended tanhmethod [12,13], the exp(−f (ζ ))-expansion method [14][15][16][17], the truncated Painleve expansion [18] and many others (see Refs. [19][20][21][22][23][24][25][26][27][28][29]). It is notable to mention that some of these approaches cannot be sometimes applied to some NLEEs. Some numerical schemes for solving NLEEs have also been established recently. Among these are the adaptive moving mesh technique [30], the finite element method, the finite differences and the Parabolic Monge-Ampere method [31]. The central purpose of this article is to use the exp(−f (ζ ))-expansion and the adaptive moving mesh techniques in extracting the exact and numerical solutions, respectively, for the following variant Boussinseq equations [32,33]: where U presents the velocity of wave, G denotes the total depth, x is the spatial derivative, t is the temporal derivative, and α and β > 0 are arbitrary constants. We also aim to verify that the used numerical technique gives reliable and successful results. Some figures are shown to verify that the behaviours of the exact and numerical solutions are almost the same.
The boundary conditions are graphically discovered from the behaviour of the analytical solution as time increases. The analytical solution becomes fixed at the end points of the physical domain. Hence, the boundary conditions of the solution are constants. This implies that the solution does not change at the end points of the domain and then we observe that u x = 0 and u xxx = 0 as x → ±∞. Hence, the boundary conditions are given by The rest of this paper is outlined as follows. Section 2 describes the exp(−f (η))-expansion method which is utilized to express the exact travelling wave solutions. In Section 3, the Hamiltonian system is explained to test the stability of the achieved results. In Section 4, we solve the variant Boussinseq equations and point out some solutions. Section 5.1 is devoted to solve the proposed equation by employing the uniform mesh scheme while Section 5.2 presents the moving mesh approach. The last section concludes this article.

The explanation of the exp(−f(η))expansion method
A comprehensive explanation to the exp(−f (η))expansion approach is clearly given in this section. The description shown here was provided in Refs. [15,16]. This technique is described in few steps illustrated as follows. Assume that where is an unknown function and TT is a polynomial of v.v., is a PDE with two independent variables x and t. In order to reduce Equation (3) to an ordinary differential equation (ODE), we introduce the transformation Expanding the derivatives in Equation (3) using the chain rule leads to the ODE where R is a polynomial in v(η) and According to this technique, the solution of Equation (5) is given in the form where the constants a 0 , a 1 , a 2 , . . . , a n will be calculated later. Furthermore, a n = 0, and f = f (η) fulfils the equation We now turn to present various cases for the solutions of Equation (8).
• If λ 2 − 4μ < 0, μ = 0, then where k is an arbitrary constant. The constants μ, λ, w, and a n are evaluated later. When we insert Equation (7) into Equation (5) and equate the coefficients of the same order of exp(−f (η)), we end up with an algebraic system that can be simply solved (using some software such as Maple or Mathematica) to evaluate the values of μ, λ, w, and a n . Substituting μ, λ, w, and a n into Equations (7) and (8) shows the solution of Equation (3).

Stability analysis
The stability of the achieved exact solutions is tested using the Hamiltonian system expressed as where ρ i (w) indicates the momentum and U and G are the obtained solutions of the system (1). The sufficient condition for the stability can be illustrated as follows:

The exact solution of the variant Boussinseq equations
This section concerns with determining the travelling wave solutions of the variant Boussinseq equations. The variant Boussinseq equations [32,33] are given by where α and β are arbitrary constants and β > 0. We begin with introducing the transformation to alter system (14) into the ODEs The ODEs in system (16) are now integrated with respect to η once to yield where p 1 and p 2 are the integration constants. Balancing the highest order u ηη and non-linear term u 2 in the first equation of (17) and the terms u ηη and u g in the second equation, we have 2N = N + 2 and 2N = N + M which leads to N = 2 and M = 2. Thus, the solutions are given by And then, The values of the constants a 0 , a 1 , a 2 , b 0 , b 1 and b 2 are shown later. The first equation in system (19) is now inserted into the first equation in system (17) and the coefficients of e −k f(η) , k ≥ 0, are equalized to zero to introduce the following algebraic system: Similarly, when we substitute the second equation in system (19) into the second equation in system (17) and equate the coefficients of e −m f(η) , m ≥ 0, to zero, we end up with Solving the previous algebraic systems (by Maple or Mathematica) gives Putting λ and μ, so that λ 2 − 4μ > 0, gives • From Case I, the exact solutions are given by • From Case II Setting λ and μ, so that λ 2 − 4μ < 0, yields • From Case I, the exact solutions are illustrated as follows: • From Case II, we have When λ 2 − 4μ = 0, μ = 0, λ = 0, we can observe that • From Case I, we obtain When λ 2 − 4μ = 0, μ = 0, λ = 0, we observe that • From Case I, we obtain • From Case II, we have where η = x − w t. To investigate the numerical results of Equation (14), we semi-discretize the spatial derivatives utilizing the centred finite differences (second-order accurate) and the temporal derivative is kept continuous. Therefore, the main equation has been converted to a system of ODEs solved utilizing the method of lines with the initial condition from Equation (24) at t = 0 and the boundary conditions from Equation (2). In the following section, we use the MATLAB ODE solver (ode15i), which is a variable order implicit time-stepping method based on the numerical differentiation formulas (NDFs), to solve the system.

Numerical results
This section is devoted to extract the numerical solutions of Equation (14) using uniform and adaptive moving mesh methods. To employ these approaches, we first introduce a variable V as follows: Hence, Equation (14) is converted to The used boundary conditions are given by Equation (2).

Numerical solutions using a uniform mesh
A uniform mesh approach, on a physical domain [0, L], is now employed to deal with the numerical results of Equation (37). The used domain is split into N subintervals [x j , x j+1 ] with fixed step size x = L/N such that where x presents a uniform width of each subinterval. The discretization of the spatial derivatives is accomplished by applying the finite difference operators. Note that the temporal differentiation is kept continuous. Thus, the discretization of Equation (37) is given by where U n+1/2 = U n+1 + U n 2 , and n = 0, 1, 2, . . . , N. The associated boundary values of Equation (2) are u t,0 = u t,N = 0 and the initial condition is generated by letting t = 0 in Equation (24).

Numerical solutions using an adaptive mesh
We now turn to employ the adaptive mesh method to obtain the numerical results of Equation (37) We now divide the physical domain to equal subintervals as follows: where x n (ζ ) = x(ζ n , t), n = 0, 1, . . . , N, with x 0 = 0, x N = L. Additionally, the nodes are presented by ζ n = n(1/ (N − 1)), n = 0, 1, . . . , N. Applying the chain rule for Equation (37) yieldṡ It is worth pointing out that the error and the convergence of the adaptive moving scheme for different MMPDEs (and their parameters) and various monitor functions are deeply considered to obtain the new mesh. We found that all these mesh equations give the same results. Consequently, we here utilize the most commonly used one which is MMPDE5 [30,[34][35][36] given by where the smoothing operator (I − γ ∂ ζ ζ ) is important to make the computation easier and faster [37,38], τ ∈ (0, 1) is a relaxation parameter, γ = √ max(ϕ) is selected by Budd et al. [31] and Walsh [37] and ϕ(U, G, x) is called a monitor function. The monitor function depends on the solutions' changes and controls the movement of the mesh so that it increases the number of the mesh points where the solution has significant variations and fewer points elsewhere. Here, we use an exceptional monitor function given by where r 1 and r 2 are constants. In this section, we fix the parameter values by r 1 = 10, r 2 = 10, k = −6, λ = 2.1, μ = 10 −1 , α = 10 −1 , β = 2 × 10 −1 , w = 1 and τ = 10 −3 . The semi-discretizations of MMPDE5 (Equation (42)) and system equation (41) are shown as We use here the average value of the monitor function so that ϕ i+1/2 = (ϕ i+1 + ϕ i )/2 and ϕ i−1/2 = (ϕ i + ϕ i−1 )/2.
The initial condition is taken by The boundary conditions of the equations given in systems (44) and (45) are converted tȯ (47) Some fictitious points are required to compute the boundaries of U xx , G x and U x The sufficient conditions for the stability (Equations (12) and (13)) are employed, and then we discover that all of the achieved solutions are stable in the interval [0,20]. Figure 1(a,b) presents the evolution time of the exact and numerical results and Figure 1(c) presents the mesh behaviour taken at t = 0 → 10 with fixed number of points N = 1000. Figure 2 shows the exact and numerical results for both U and G at t = 5 and Figure 2. The comparison between the exact solutions evaluated by (24) and the numerical results obtained by solving the system (45) for U and G at t = 5 and N = 1000. N = 1000. All of the parameter values are fixed, as mentioned above. Figure 3(a,b) illustrates the behaviour of the travelling waves for both the exact and numerical solutions. We note that, from all of the above figures, the exact and the numerical results are almost identical. The insets in Figure 2 show the steep front regions which take more points than elsewhere. Therefore, the results appear almost equal in these regions for both U(x, t) and G(x, t). Table 1 illustrates L 2 error and CPU time taken to reach t = 5 for the adaptive moving and uniform mesh (using MMPDE5 (42) and the modified monitor function (43)) approaches. The numerical results are obtained at t = 5. Figure 4 summaries the error columns in Table 1 (solid blue line for the adaptive mesh scheme and solid red line for the uniform mesh scheme). It can be easily observed that the error for the adaptive moving mesh is much smaller and is accomplished utilizing less number of nodes compared to the    uniform mesh method. However, regarding the CPU time consumed, the uniform mesh scheme consumes less time to reach t = 5 compared to the adaptive moving mesh technique for the same number of points. This can be attributed to the additional adaptive functions which are required to be simultaneously resolved along with the PDE. As the number of points increases, the error for the adaptive moving scheme sharply decreases with a slight increase in the CPU time. Therefore, it can be definitely concluded that the adaptive moving mesh technique is more computationally effective than the uniform mesh method.

Conclusions
The discussion of this article concentrates on constructing the travelling wave solutions and the numerical solutions of the variant Boussinseq equations by applying the exp(−f (η))-expansion and the adaptive mesh approaches, respectively. The given 2D and 3D figures show that the solutions agree and coincide together. As can be seen in Table 1, L 2 error for the adaptive mesh scheme vanishes for a very small x. The error for the uniform mesh scheme is found larger than the error for the adaptive mesh scheme. The achieved results have been compared to each other and found that the performance of the adaptive scheme is effective and appropriate to be utilized in high-order PDEs.

Disclosure statement
No potential conflict of interest was reported by the author(s).