1. Introduction
The numerical solution of partial differential equations (PDEs) represents a fundamental component in numerous scientific and engineering applications. In recent years, there has been a growing emphasis on the development of high-order finite difference methods aimed at improving both the accuracy and efficiency of solutions for a wide range of PDEs. A particular area of interest is the Boussinesq equation, which was first introduced by Joseph Boussinesq in 1872 to model nonlinear wave propagation in dispersive media [
1]. As a result, this equation has been extensively applied across multiple domains, including the investigation of ion-sound waves in plasma, solitary waves in vascular systems, tsunami waves in oceanography and coastal sciences, as well as pressure waves in liquid-gas bubble mixtures, among others [
2,
3].
This article investigates the sixth-order Boussinesq equation, taking into account the effects of surface tension, as outlined in [
4]:
with the specified initial conditions
In this context, the parameters
,
,
and
are non-negative, while
k and
c are constants. Furthermore,
is a positive integer, and
and
are two specified smooth functions.
Solitary waves have attracted considerable interest across various fields of physics and engineering. Biswas [
5] successfully derived an exact solitary wave solution for the Korteweg-de Vries equation, which is characterized by power law nonlinearity and time-varying coefficients. The singular solitary wave solution for the Rosenau-KdV equation is presented in [
6], while the specific 1-soliton solution for the Rosenau-KdV-RLW equation is detailed in [
7]. Eq. (
1) is frequently encountered in numerous mathematical models concerning solitary wave propagation, as emphasized by Lu and Helal [
8,
9]. Furthermore, Biswas et al. [
4] provide a variety of solitary wave, shock wave, and singular solitary wave solutions for the Boussinesq equation (
1), which incorporates the effects of surface tension. Under the assumption of solitary wave behavior, it is observed that the solution to the sixth-order Boussinesq equation (
1) and its derivatives approach zero as
, that is [
10,
11]
then Eq. (
1) is linked to the following global conservation law:
This research focuses on the solitary wave solution and seeks to develop a conservative linear finite difference scheme for the sixth-order Boussinesq equation (
1). To effectively implement the numerical method, we analyze a sufficiently large yet finite spatial domain, denoted as
, instead of the unbounded interval
. This strategy ensures that the solution remains sufficiently small at the boundaries of the domain. Consequently, we consider the following boundary conditions:
where
.
The Boussinesq equation and its variants have been thoroughly analyzed through both theoretical and numerical approaches. In the theoretical context, Esfahani and Farah [
12] have conducted an investigation into a nonlinear sixth-order Boussinesq equation
The authors established the local well-posedness of Eq. (
4) within the framework of non-homogeneous Sobolev spaces
, where
and
. Following this, Wang and Esfahani [
13] conducted an investigation into the sixth-order Boussinesq equation as follows:
The authors have demonstrated that the problem described by Eq. (
5) is globally well-posed within the Sobolev space
, where
and
.
In a numerical context, Feng et al. [
14] proposed a symmetric three-level implicit difference scheme that exhibits second-order accuracy for the following sixth-order Boussinesq equation
A linearized stability analysis indicates that the difference scheme, which incorporates a free parameter denoted as
, exhibits stability for values of
that are equal to or greater than 1/4. In 2017, Kolkovska and Vucheva [
15] introduced a nonlinear second-order difference scheme aimed at addressing the sixth-order Boussinesq problem, as detailed below:
However, the scheme presented in [
15] demonstrates conditional stability, which is dependent on a strict limitation on the subsequent ratio
There exists a significant lack of theoretical analysis regarding the solvability and convergence properties of the scheme proposed by Kolkovska [
15]. More recently, Arslan [
16] has introduced a novel methodology that combines the differential transform method with the finite difference technique to obtain approximate solutions for singularly perturbed ill-posed problems and sixth-order Boussinesq equations. In a separate study, Zhang et al. [
17] introduced a meshless numerical technique known as the generalized finite difference method, which has demonstrated efficacy in addressing enhanced Boussinesq-type equations. This method is particularly advantageous for simulating wave propagation over irregular bottom topographies.
In the field of physics, the principle of conservation of mass and energy holds critical importance. A numerical scheme that fails to uphold local conservation may produce non-physical outcomes [
18]. Therefore, the preservation of mass and energy is a fundamental consideration in the development of numerical schemes aimed at solving PDEs. A considerable amount of research has been devoted to the development of conservative finite-difference schemes applicable to various nonlinear wave equations. For instance, Deng et al. [
19] proposed energy-preserving finite difference methods that are applicable to systems of nonlinear wave equations in two dimensions. Furthermore, Bayarassou et al. [
20] investigated a high-order conservative linearized finite difference scheme for the regularized long wave (RLW)-Korteweg-de Vries equation. Additionally, Hou et al. [
21] introduced an energy-preserving, high-order compact finite difference scheme specifically designed for two-dimensional nonlinear wave equations. Nanta et al. [
22] presented a wave model that integrates the classical Camassa-Holm equation with the BBM-KdV equation, incorporating dual-power law nonlinearities while ensuring the preservation of energy conservation properties.
This article aims to present a conservative difference scheme for the resolution of the sixth-order Boussinesq problem as defined in Eqs. (
1)-(
3). Furthermore, the article will examine the solvability, convergence, and stability of the proposed scheme. The problem outlined in Eqs. (
1)-(
3) holds considerable physical significance.
Theorem 1.
Consider u as the solution to the sixth-order Boussinesq problem defined by Eqs. (1)-(3). Let us denote , and assume that the supplementary boundary conditions are given by and . Under these conditions, it can be concluded that for all , where
Proof. Let
. We substitute this expression into Eq. (
1) and subsequently take the antiderivative twice. This process results in the following equation:
Consequently, we obtain further results
Utilizing the boundary conditions delineated in Eq. (
3) along with the additional assumptions, we derive the following results
Since
which yields
Consequently, from Eqs. (
8)-(
12), we deduce that
, which indicates that
for
. This completes the proof. □
It is essential to emphasize that if the parameters
,
,
,
,
k,
c and the initial conditions
and
satisfy
, then
can be defined as the energy at time
t. Moreover, the solution to Eqs. (
1)-(
3) complies with the principle of energy conservation, where
The subsequent sections of this manuscript are organized as follows: In Section (
Section 2), we introduce a compact finite difference scheme designed for the resolution of the sixth-order Boussinesq equation (
1)-(
3). Section (
Section 3) rigorously establishes the discrete conservation properties associated with the mass change rate and energy. A comprehensive theoretical analysis addressing the scheme’s solvability, convergence, and stability is presented in Section (
Section 4). Section (
Section 5) includes numerical experiments that validate the theoretical results. Finally, Section (
Section 6) provides a concise summary of the findings.
2. Compact Difference Scheme
In this section, we present a compact and conservative difference scheme specifically developed to tackle the sixth-order Boussinesq problem as delineated in Eqs. (
1)-(
3). For two positive integers
J and
N, we define the spatial and temporal discretization parameters as
and
, respectively. The uniform grid points are established as
and
, where
and
. We denote
as the approximation of the function
u at the point
and
where
. The Sobolev space is defined as described in [
23]:
We will now proceed to establish the following notations:
In the subsequent section, we will develop a compact difference scheme for Eqs. (
1)-(
3). By setting
Eq. (
1) can be expressed as
. By utilizing the Taylor expansion for the variable
x, we obtain
Furthermore, Eq. (
Section 2) is associated to
Substituting Eq. (
15) into Eq. (
14) yields
Through the application of second-order accuracy in our approximations, we achieve
Consequently, the proposed compact finite difference scheme can be expressed as follows:
It is essential to recognize that since
,
and
hold based on Eq. (18), we may proceed with the assumption that
for the sake of simplicity. Similarly, we can also assume that
, where
and
denote ghost points under the condition
. Consequently, it can be concluded that
.
In this context, we employ the following methodology to derive
:
From Eq. (17), we get
that is
which yields
where
.
5. Numerical Experiments
In this section, we conduct a series of numerical experiments aimed at assessing the efficacy and accuracy of the theoretical analysis presented in the preceding sections. Following this, we will quantify the corresponding errors utilizing the designated error norms [
27,
29,
31]:
Example 1.
which is the following improved Boussinesq equation [2,27]
the solitary wave solution is expressed as follows
where A, and are given constants. The initial conditions may be established by evaluating Eq. (76) and its derivative for t at the point :
First, we choose the following parameters:
We then compare our findings with the numerical results presented in [
2]. The computational accuracy for both temporal and spatial variables is evaluated at
and is detailed in
Table 1 and
Table 2, where
in
Table 1 and
in
Table 2. The data presented in
Table 1 and
Table 2 indicate that the current difference scheme (
16)-(18) demonstrates significantly superior performance in comparison to the finite volume element scheme described in [
2].
The present analysis investigates the solitary wave and its corresponding numerical solutions, as illustrated in
Figure 1. The spatial domain is defined as
and
, with parameters set at
and
at the time instances
.
Figure 1 indicates that the patterns generated by the current scheme, as described by Eqs. (
16)-(18), demonstrate a significant degree of agreement with the solitary wave solutions. Furthermore, the numerical solutions for the auxiliary variable
, as estimated by Eq. (
21), are visually represented in
Figure 2. The graphs displaying the approximate solutions
and
, obtained from the difference scheme (
16)-(18) with the parameters
,
,
, and
at
, are shown in
Figure 3. It is observed that when
, the solitary wave propagates to the right, while for
, the solitary wave propagates to the left. These results are consistent with those reported in [
27].
Finally, we present the discrete mass
, the discrete rate of mass change
and the discrete energy
at various time intervals, as detailed in
Table 3. The parameters are defined with
,
,
and
. Analysis of
Table 3 indicates that the discrete rate of mass change
associated with soliton wave propagation is minimal, indicating that the mass change
is conservative. Additionally, it is noted that the energy
also exhibits conservative properties.
Example 2.
which is the sixth-order Boussinesq equation [15]
the solitary wave solution is given by
In this example, we first calculate the solution over the spatial domain defined by
with a grid size of
and a time step
. The comparison between the solitary wave and the numerical solutions at the time instances
, 20, and 30 is presented in
Figure 4. The results indicate that the errors across all grid points remain below
.
Figure 5 illustrates the numerical solutions at various time intervals, specifically
, 10, 20, 30, and 40. A detailed analysis of
Figure 5 indicates that the wave heights at these specified time points are nearly indistinguishable.
Figure 6 presents the graphs of the approximate solutions
obtained from the current scheme described by Eqs. (
16)-(18) at times
, 10, 15, and 20, with the parameters set to
,
,
and
.
Additionally,
Table 4 presents a summary of the discrete mass
, the discrete rate of mass change
, and the discrete energy
at various time points
, 5, 10, 15, and 20. The findings suggest that both the discrete rate of mass change
and the discrete energy
exhibit conservative behavior in the discrete sense.
Example 3.
Consider the following sixth-order Boussinesq equation with the effects of surface tension
where , , , , . In their research, Biswas et al. [4] examined the solitary wave solution corresponding to Eq. (83), that is
To perform numerical simulations, the initial conditions can be determined by evaluating Eq. (84) and its derivative at the time point , that is
To demonstrate the efficacy of the present difference scheme, we provide numerical results that pertain to the errors and rates of convergence. The parameters are configured with
and
at time
, as illustrated in
Table 5, where
and
. From
Table 5, it is evident that the convergence rate achieved by the current finite difference scheme is approximately 4.0, which is consistent with the theoretical order of convergence specified in Theorem 10.
Figure 7 depicts the numerical solutions and absolute errors at times
and 20 with
,
,
,
and
. Furthermore,
Table 6 presents additional information regarding the discrete mass
, the discrete rate of mass change
and the discrete energy
at various time intervals
,
,
,
,
,
, and
. It is clear that both the rate of mass change
and the energy
exhibit conservative behavior in the discrete sense.
Furthermore, for
, Eq. (
83) exhibits two conserved quantities [
4]:
Furthermore, for any arbitrary
, Eq. (
83) possesses two conserved quantities [
4]:
We will evaluate our numerical method by employing these invariants and selecting the approximate values as follows:
To demonstrate the conservative properties of the current difference scheme represented by Eqs. (
16)-(18), we have compiled the invariants
,
,
,
,
,
and
at various temporal intervals, as presented in
Table 7 and
Table 8. The data illustrated in these tables indicate that the conservative quantities
and
and
remain approximately constant over a time frame of up to 30 units. Consequently, the compact scheme proposed, as outlined by Eqs. (
16)-(18), demonstrates both efficiency and reliability for simulations conducted over long periods.
Figure 1.
Solitary wave solutions of at and numerical solutions at , 10, 15, 20, , (left) and , (right) for Example 1.
Figure 1.
Solitary wave solutions of at and numerical solutions at , 10, 15, 20, , (left) and , (right) for Example 1.
Figure 2.
Numerical solutions of the auxiliary variable V at , 5, 10, 15, 20, , (left) and , (right) for Example 1.
Figure 2.
Numerical solutions of the auxiliary variable V at , 5, 10, 15, 20, , (left) and , (right) for Example 1.
Figure 3.
3D plots of the numerical solutions and at with , for Example 1.
Figure 3.
3D plots of the numerical solutions and at with , for Example 1.
Figure 4.
Numerical and solitary wave solutions and absolute error at , 20 and 30 for Example 2.
Figure 4.
Numerical and solitary wave solutions and absolute error at , 20 and 30 for Example 2.
Figure 5.
Numerical solutions of U and V with and at , 10, 20, 30, 40, , for Example 2.
Figure 5.
Numerical solutions of U and V with and at , 10, 20, 30, 40, , for Example 2.
Figure 6.
3D plots of the numerical solutions at , 10, 15 and 20 with , , , for Example 2.
Figure 6.
3D plots of the numerical solutions at , 10, 15 and 20 with , , , for Example 2.
Figure 7.
Numerical solutions and absolute errors at and 20 with , , , , for Example 3.
Figure 7.
Numerical solutions and absolute errors at and 20 with , , , , for Example 3.
Table 1.
Comparative analysis of spatial errors and convergence rates for Example 1.
Table 1.
Comparative analysis of spatial errors and convergence rates for Example 1.
h |
[2] |
rate |
|
rate |
[2] |
rate |
|
rate |
0.4 |
9.05
|
− |
5.44
|
− |
4.50
|
− |
3.15
|
− |
0.2 |
2.25
|
2.01 |
3.21
|
4.08 |
1.13
|
1.99 |
2.05
|
3.94 |
0.1 |
5.62
|
2.00 |
1.90
|
4.07 |
2.81
|
2.00 |
1.19
|
4.09 |
Table 2.
Comparison of the temporal errors and convergence rates for Example 1.
Table 2.
Comparison of the temporal errors and convergence rates for Example 1.
h |
[2] |
rate |
|
rate |
[2] |
rate |
|
rate |
0.04 |
8.69
|
− |
2.76
|
− |
4.32
|
− |
1.51
|
− |
0.02 |
2.19
|
1.98 |
6.58
|
2.06 |
1.10
|
1.97 |
3.63
|
2.05 |
0.01 |
5.53
|
1.99 |
1.59
|
2.04 |
2.78
|
1.99 |
9.04
|
2.00 |
Table 3.
Discrete mass , discrete rate of mass change and the discrete energy at various time intervals for Example 1.
Table 3.
Discrete mass , discrete rate of mass change and the discrete energy at various time intervals for Example 1.
T |
|
|
|
0 |
3.99999998662202 |
0.00000000770792 |
2.59519634357874 |
5 |
3.99999998684879 |
0.00000000748522 |
2.59519634801190 |
10 |
3.99999998721201 |
0.00000000711875 |
2.59519636602745 |
15 |
3.99999998755727 |
0.00000000676334 |
2.59519639725431 |
20 |
3.99999998788500 |
0.00000000641617 |
2.59519644155837 |
Table 4.
Discrete mass , discrete rate of mass change and the discrete energy with and at various time intervals for Example 2.
Table 4.
Discrete mass , discrete rate of mass change and the discrete energy with and at various time intervals for Example 2.
T |
|
|
|
0 |
8.44807966748734 |
1.88271050488737 |
15.35813895730075 |
5 |
8.44807966747724 |
1.88289006817254 |
15.35813896808486 |
10 |
8.44807966744347 |
1.88285678675310 |
15.35813901287684 |
15 |
8.44807966738941 |
1.88269314441971 |
15.35813909109440 |
20 |
8.44807966731402 |
1.88292585889068 |
15.35813920256866 |
Table 5.
Errors and convergence rates of the current scheme when and at for Example 3.
Table 5.
Errors and convergence rates of the current scheme when and at for Example 3.
p |
|
|
|
|
|
|
1.530067340E-02 |
9.709362247E-04 |
5.777896747E-05 |
|
rate |
− |
3.9780748 |
4.07076021 |
|
|
1.015003674E-02 |
6.556924437E-04 |
4.127486212E-05 |
|
rate |
− |
3.9523219 |
3.9896839 |
|
|
7.862155174E-02 |
4.719448787E-03 |
3.012389794E-04 |
|
rate |
− |
4.0582346 |
3.9696380 |
|
|
7.973541804E-02 |
4.682164328E-03 |
3.018001249E-04 |
|
rate |
− |
4.0899732 |
3.9555103 |
Table 6.
Discrete mass , discrete rate of mass change and the discrete energy with , and at various time intervals for Example 3.
Table 6.
Discrete mass , discrete rate of mass change and the discrete energy with , and at various time intervals for Example 3.
T |
|
|
|
0 |
3.02090976925920 |
1.52260617289068 |
8.70669460321910 |
5 |
3.02090976827031 |
1.52400439056038 |
8.70669481288416 |
10 |
3.02090976455933 |
1.52404485849716 |
8.70669563361528 |
15 |
3.02090975829348 |
1.52179604152168 |
8.70669698162941 |
20 |
3.02090974941932 |
1.52420689642572 |
8.70669874989782 |
25 |
3.02090973801432 |
1.52375610223547 |
8.70670081739405 |
30 |
3.02090972410221 |
1.52240106874032 |
8.70670306135933 |
Table 7.
Conservative quantities , and the discrete quantities , and with , and at various time intervals for Example 3.
Table 7.
Conservative quantities , and the discrete quantities , and with , and at various time intervals for Example 3.
T |
|
|
|
|
|
0 |
0.77574398822 |
0.72568182168 |
2.99999999979 |
1.06047283470 |
3.78558199462 |
5 |
0.77574415451 |
0.72568197184 |
2.99999999869 |
1.06079683168 |
3.78558200635 |
10 |
0.77574474392 |
0.72568250450 |
2.99999999479 |
1.06080624484 |
3.78558205484 |
15 |
0.77574571092 |
0.72568337974 |
2.99999998836 |
1.06028637400 |
3.78558213956 |
20 |
0.77574703291 |
0.72568457913 |
2.99999997943 |
1.06084514881 |
3.78558226063 |
25 |
0.77574867931 |
0.72568607763 |
2.99999996804 |
1.06074072374 |
3.78558241692 |
30 |
0.77575061250 |
0.72568784435 |
2.99999995426 |
1.06042851618 |
3.78558260687 |
Table 8.
Conservative quantities , and the discrete quantities , and with , and at various time intervals for Example 3.
Table 8.
Conservative quantities , and the discrete quantities , and with , and at various time intervals for Example 3.
T |
|
|
|
|
|
0 |
0.00000000879 |
3.02090976945 |
3.02090976927 |
1.52260617358 |
8.70669460355 |
5 |
0.00000001689 |
3.02090976957 |
3.02090976906 |
1.52318643159 |
8.70669464967 |
10 |
0.00000002587 |
3.02090976979 |
3.02090976876 |
1.52365257478 |
8.70669471985 |
15 |
0.00000003535 |
3.02090977012 |
3.02090976836 |
1.52400439292 |
8.70669481421 |
20 |
0.00000004386 |
3.02090977051 |
3.02090976788 |
1.52424172749 |
8.70669493260 |
25 |
0.00000005221 |
3.02090977097 |
3.02090976732 |
1.52436447198 |
8.70669507504 |
30 |
0.00000006009 |
3.02090977148 |
3.02090976668 |
1.52437256971 |
8.70669524020 |