# Solving the Generalized Regularized Long Wave Equation Using a Distributed Approximating Functional Method

Authors: Edson Pindza; Eben Maré
Journal: International Journal of Computational Mathematics
 (2014)

---

## Abstract

The generalized regularized long wave (GRLW) equation is solved numerically by using a distributed approximating functional (DAF) method realized by the regularized Hermite local spectral kernel. Test problems including propagation of single solitons, interaction of two and three solitons, and conservation properties of mass, energy, and momentum of the GRLW equation are discussed to test the efficiency and accuracy of the method. Furthermore, using the Maxwellian initial condition, we show that the number of solitons which are generated can be approximately determined. Comparisons are made between the results of the proposed method, analytical solutions, and numerical methods. It is found that the method under consideration is a viable alternative to existing numerical methods.


---

## Body

## 1. Introduction

Nonlinear partial differential equations (NPDEs) are widely used to describe complex phenomena in various fields of sciences, such as physics, biology, and chemistry. Various methods have been devised to find the exact and approximate solutions of nonlinear partial differential equations [1, 2] in order to provide more information for understanding physical phenomena arising in numerous scientific and engineering fields.In this paper, we consider an important nonlinear wave equation, the generalized regularized long wave (GRLW) equation of the form(1)
u
t
+
u
x
+
ε
u
p
u
x
-
μ
u
x
x
t
=
0
,
(
x
,
t
)
∈
R
×
[
0
T
]
,
where p is a positive integer and ε and μ are positive constants. The boundary conditions are given by u
→
0 as |
x
|
→
∞. This equation has strong connections to the generalized Korteweg-de Vries (GKdV) equation and has solitary solutions [3]. The GRLW equation was proposed many years ago by Peregrine [4, 5] as a model for small-amplitude long waves on the surface of water in a channel. The study of the GRLW equation gives the opportunity of investigating the creation of secondary solitary waves to obtain insight into the corresponding processes of particle physics [6].Forp
=
1, this equation is reduced to the regularized long wave equation as an important equation in physics media, describing phenomena with weak nonlinearity and dispersive waves [7]. The existence and uniqueness of the solution of the RLW equation are given in [8]. Its analytical solutions were obtained under restrictive initial and boundary conditions [8]; therefore, the use of numerical methods has been advocated. Several numerical methods for the solution of the RLW equation have been introduced in the literature. These include finite- difference methods [9], Fourier pseudospectral methods [10], and finite-element methods based on Galerkin and collocation principles [11].Another special case of the GRLW is obtained whenp
=
2. This corresponds to the modified regularized long wave (MRLW) equation. The MRLW was solved numerically by explicit finite-difference methods [12], Galerkin methods using cubic B-spline finite elements [13], Galerkin and Petrov-Galerkin methods using quadratic B-spline finite elements [14], spectral method based on Chebyshev polynomials [15], a least square technique using linear space-time finite elements and Petrov-Galerkin finite-element scheme with shape functions taken as quadratic B-spline functions [16], and the collocation method with quintic B-spline [17]. Recently, Mokhtari and Mohammadi [18] proposed a meshfree technique based on a global Sinc-collocation approximation to solve the GRLW equation.In this paper, we propose a Hermite distributed approximating functional (DAF) method for solving the GRLW equation. The DAF was originally introduced by Hoffman et al. [19] and Hoffman and Kouri [20] as a computational tool for treating a variety of problems in physics and chemistry, with particular focus on both time-dependent and time-independent Schrödinger equations. Later, the DAF was also used to solve the Fokker-Planck equation [21, 22]. The study of the DAFs indicates that these methods deliver spectral accuracy when used to solve the Fokker-Planck equation with nonlinear drift and diffusion coefficients included. DAFs are able to provide an analytical representation of a function and its derivatives in terms of a discrete set of values of this function. This is central to its successful use in various computational applications.The layout of this paper is as follows. We describe the formulation of the DAF formalism and illustrate the choice of an optimal bandwidth in Section2. Numerical results illustrating the merits of the new scheme are given in Section 3. Finally, we present a brief conclusion in Section 4.


## 2. Distributed Approximating Functionals

Distributed approximating functionals (DAFs) are defined as a mapping that approximates a certain set of continuousL
2 functions onto itself within a specified accuracy [22]. Their main success in various computational practices is derived from their inherent ability to furnish analytical representation of a function and its derivatives at collocation points. These DAFs have different kernels depending on the application of interest. Here, we limit our study to the DAF of Hermite type.DAF methods are a general framework for constructing local spectral methods. They provide discrete approximations to the singular convolution:(2)
f
(
x
)
=
∫
-
∞
∞
‍
δ
(
x
-
x
′
)
f
(
x
′
)
d
x
′
,
where T
(
t
-
x
) is a singular kernel. In this paper, T will designate the delta distribution δ. The basic equation associated with the Dirac delta function δ
(
x
) is
(3)
∫
-
∞
∞
‍
δ
(
x
)
f
(
x
)
d
x
=
f
(
0
)
,
where f is any function that is continuous at x
=
0. The delta function has the following properties:
(4)
∫
-
∞
∞
‍
δ
(
x
-
x
′
)
f
(
x
)
d
x
=
f
(
x
′
)
,
∫
-
∞
∞
‍
δ
′
(
x
)
f
(
x
)
d
x
=
-
f
′
(
0
)
,
δ
(
x
a
)
=
|
a
|
δ
(
x
)
.
As the delta distribution is singular, it is necessary to approximate it so that it can be digitized in a computer. There are many functions that provide accurate approximations of the delta kernel and its derivatives.Some of the approximations, such as the Shannon (Sinc) kernel and the Dirichlet kernel, give rise to classical Fourier spectral methods in essence. In recent work, the so-called delta regularized kernel [23] was proposed to solve some applied mechanics problems. Shannon’s kernel is regularized as
(5)
δ
Δ
,
σ
(
x
-
n
Δ
)
=
sin
⁡
(
π
/
Δ
)
(
x
-
n
Δ
)
⁡
(
π
/
Δ
)
(
x
-
n
Δ
)
exp
⁡
[
-
(
x
-
n
Δ
)
2
2
σ
2
]
.
Dirichlet’s regularized kernel [24] is given by
(6)
δ
Δ
,
σ
(
x
-
n
Δ
)
=
sin
[
(
π
/
Δ
)
(
x
-
n
Δ
)
]
(
2
M
′
+
1
)
sin
[
(
π
/
Δ
)
(
(
x
-
n
Δ
)
/
(
2
M
′
+
1
)
)
]
×
exp
⁡
[
-
(
x
-
n
Δ
)
2
2
σ
2
]
,
where M
′ is a parameter and Δ is the grid spacing. In comparison to Shannon’s kernel, the Dirichlet kernel has one more parameter, namely, M. This parameter can be optimized to achieve a better accuracy. In practice, we choose M sufficiently large to obtain accurate results. Clearly, at the limit M
→
∞, the Dirichlet kernel converts to Shannon’s kernel.Another important kernel, the regularized Lagrange kernel (RLK) [25], defined by
(7)
δ
Δ
,
σ
(
x
-
n
Δ
)
=
∏
i
=
k
-
M
′
,
i
≠
k
i
=
k
+
M
′
‍
x
-
x
i
n
Δ
-
x
i
exp
⁡
[
-
(
x
-
n
Δ
)
2
2
σ
2
]
,
was constructed by regularizing the classic Lagrange polynomial [26].More recently, a new class of DAF approaches, the so-called Hermite DAF (HDAF) [27], was introduced for spectrally accurate numerical solutions of first- and second-order partial differential equations. In this paper, we confine our attention to the HDAF.We first define the regularized Hermite function given by(8)
h
m
(
x
)
=
e
-
x
2
H
m
(
x
)
,
where H
m
(
x
) is the classic Hermite polynomial [28]. The HDAF uses even Hermite polynomials h
2
m to construct δ
(
x
-
x
′
) as
(9)
δ
σ
,
M
(
x
-
x
′
)
=
1
σ
∑
m
=
0
M
/
2
‍
(
-
1
4
)
m
1
2
π
m
!
h
2
m
(
x
-
x
′
2
σ
)
,
where M is the highest polynomial degree involved in the construction of the HDAF and σ is its bandwidth. Note that both parameters control the accuracy of the approximation. For instance, if we fix σ, the HDAF δ
σ
,
M
(
x
-
x
′
) becomes identical to the δ function when the degree of the polynomial M tends to infinity. We write
(10)
lim
⁡
M
→
∞
δ
σ
,
M
(
x
-
x
′
)
=
δ
(
x
-
x
′
)
.
Similarly, if we keep M fixed, the HDAF δ
σ
,
M
(
x
-
x
′
) becomes identical to the δ function when σ tends to zero; that is,
(11)
lim
⁡
σ
→
0
δ
σ
,
M
(
x
-
x
′
)
=
δ
(
x
-
x
′
)
.
The HDAF mapping can be used to sample an arbitrary function at discrete points. If these nodes form an equispaced grid, (2) can be approximated by quadrature as
(12)
f
(
x
)
≈
h
x
∑
j
=
1
M
‍
δ
σ
,
M
(
x
-
x
j
)
f
(
x
j
)
,
where h
x is the grid spacing. Consider the definition of the lth derivative of the f
(
x
) defined in (2) as
(13)
d
l
f
(
x
)
d
x
l
=
f
(
l
)
(
x
)
=
∫
-
∞
∞
‍
δ
(
l
)
(
x
-
x
′
)
f
(
x
′
)
d
x
′
.
Usually, this is a purely formal expression because the derivative of the δ function is not defined and the operation has to be rolled over to the test function f
(
x
) by partial integration. However, the derivative can be evaluated immediately, when using the HDAF approximation. The integral in (13) is computed very accurately using Gauss quadrature rule [29]. Therefore, (13) becomes
(14)
f
(
l
)
(
x
)
≈
h
x
∑
j
=
1
M
‍
δ
σ
,
M
(
l
)
(
x
-
x
j
)
f
(
x
j
)
,
l
=
1,2
,
…
,
where δ
σ
,
M
(
l
)
(
x
-
x
j
), the lth derivative of HDAF, is given by
(15)
δ
σ
,
M
(
l
)
(
x
-
x
′
)
=
(
-
1
)
l
2
l
/
2
σ
l
+
1
∑
m
=
0
M
/
2
‍
(
-
1
4
)
m
1
2
π
m
!
h
2
m
+
l
(
x
-
x
′
2
σ
)
.
We note that the operator matrices (15) can be computed most efficiently, because the differentiating DAF matrices have Toeplitz structure. Equations (2) and (14) are the most useful forms for solving ordinary differential equations and partial differential equations. A significant advantage of the Hermite DAF is that it transforms ordinary and partial differentiations into a matrix vector multiplication that involves highly banded matrix representations of derivatives of local methods such as finite-difference and finite-element methods while preserving exponential accuracy of global methods such as spectral methods. Therefore, the Hermite DAF has sufficient flexibility to handle complicated boundary conditions and geometries, like finite-difference and finite-element methods, but with an accuracy of the same order as spectral methods.Note that (15) assumes that the derivative of δ and f approaches zero as x
→
±
∞; this is convenient for approximating soliton solutions of the RLW equation, as they tend to zero as x
→
±
∞. However, in general this is not the case and the complete numerical approximation must provide a way to handle boundary conditions.In (15), if the kernel, δ
σ
,
M
(
l
), is fixed to be symmetric (or antisymmetric) and invariant by translation, there must be cases where f
(
x
k
) are located outside of the computational domain, [
a
,
b
], and their values are undefined there. In the present algorithm, such f
(
x
k
) are to be obtained by boundary conditions. To handle Dirichlet boundary conditions, such f
(
x
k
) are taken to be f
(
a
) or f
(
b
). For periodic boundary conditions, such f
(
x
k
) are replaced by their corresponding values inside the domain [
a
,
b
]. For Neumann boundary conditions, f
(
x
k
) values are determined by f
(
a
) and f
′
(
a
) or (f
(
b
) and f
′
(
b
)). If f
(
x
) is antisymmetric around the boundary point then values of f
(
x
k
) outside the domain [
a
,
b
] are replaced by their corresponding -
f
(
x
k
) inside the domain [
a
,
b
]. Similarly, if f
(
x
) is symmetric around the boundary point, then values of f
(
x
k
) outside the domain [
a
,
b
] are replaced by their corresponding f
(
x
k
) inside the domain [
a
,
b
].Another significant advantage of the DAF methods over Sinc methods is that their matrices are banded (see Figure1(b)). Therefore, they are more efficient than Sinc methods.Illustration of the matrix structures of Sinc methods (a) and HDAF methods (b).
(a)
(b)


### 2.1. Choosing the Optimal Bandwidth

The limit relations (10) and (11) imply that there exists a connection between the number of expansion terms M and the bandwidth σ of the Hermite DAF. It is natural to think of the expansion order M as a measure of predefined accuracy so one might ask how to choose σ optimally with respect to a given M. One possible approach is to make the approximation as accurate as possible on the discrete grid. Equation (9) is exact, if
(16)
δ
σ
,
M
(
x
i
-
x
j
)
=
Δ
i
j
h
x
,
where Δ
i
j represents the Kronecker-Δ function. Therefore, for i
=
j, using relations (16) and (9), one obtains
(17)
σ
=
h
x
2
π
∑
m
=
0
M
/
2
‍
1
m
!
(
-
1
4
)
m
h
2
m
(
0
)
.
In the remainder of the paper, the bandwidth is chosen according to (17).


### 2.2. Application

In this section, we develop a space discretization of the GRLW equation using the HDAF. This equation was proposed by Peregrine [4, 5] as a model for small-amplitude long waves on the surface of water in a channel and is represented mathematically as
(18)
u
t
+
u
x
+
ε
u
p
u
x
-
μ
u
x
x
t
=
0
,
(
x
,
t
)
∈
R
×
[
0
T
]
.
To apply the HDAF for solving the GRLW equation, firstly, we present the following notations:
(19)
u
(
x
,
t
)
=
h
x
∑
j
=
1
N
‍
u
j
(
t
)
δ
σ
,
M
(
x
-
x
j
)
,
u
x
(
x
,
t
)
=
h
x
∑
j
=
1
N
‍
u
j
(
t
)
δ
σ
,
M
(
1
)
(
x
-
x
j
)
,
u
x
x
(
x
,
t
)
=
h
x
∑
j
=
1
N
‍
u
j
(
t
)
δ
σ
,
M
(
2
)
(
x
-
x
j
)
,
u
t
(
x
,
t
)
=
h
x
∑
j
=
1
N
‍
u
j
′
(
t
)
δ
σ
,
M
(
x
-
x
j
)
.
We substitute (19) into (18) at collocation point x
=
x
i and we obtain
(20)
h
x
∑
j
=
1
N
‍
u
j
′
(
t
)
δ
σ
,
M
(
x
i
-
x
j
)
+
h
x
∑
j
=
1
N
‍
u
j
(
t
)
δ
σ
,
M
(
1
)
(
x
i
-
x
j
)
+
ε
(
h
x
∑
j
=
1
N
‍
u
j
(
t
)
δ
σ
,
M
(
x
i
-
x
j
)
)
p
×
(
h
x
∑
j
=
1
N
‍
u
j
(
t
)
δ
σ
,
M
(
1
)
(
x
i
-
x
j
)
)
-
μ
h
x
∑
j
=
1
N
‍
u
j
′
(
t
)
δ
σ
,
M
(
2
)
(
x
i
-
x
j
)
=
0
.In order to write (20) in matrix form, we introduce the following matrix and vector notations:
(21)
u
=
[
u
1
,
u
2
,
…
,
u
N
]
T
,
D
(
1
)
=
δ
σ
,
M
(
1
)
(
x
i
-
x
j
)
,
D
(
2
)
=
δ
σ
,
M
(
2
)
(
x
i
-
x
j
)
.
Consequently, we have
(22)
M
u
˙
=
f
(
u
,
t
)
,
where M
=
μ
D
2
-
I, with I the identity matrix, f
(
u
,
t
)
=
D
1
(
(
ε
/
(
p
+
1
)
)
u
p
+
1
+
u
).Commonly used time-stepping schemes, either explicit or implicit, can be easily used along with the HDAF algorithm. In the present study, we consider the classical fourth-order Runge-Kutta method computed as(23)
u
m
+
1
=
u
m
+
1
6
(
k
1
+
2
k
2
+
2
k
3
+
k
4
)
,
where
(24)
k
1
=
Δ
t
f
(
u
m
,
t
m
)
,
k
2
=
Δ
t
f
(
u
m
+
1
2
k
1
,
t
m
+
1
2
Δ
t
)
,
k
3
=
Δ
t
f
(
u
m
+
1
2
k
2
,
t
m
+
1
2
Δ
t
)
,
k
4
=
Δ
t
f
(
u
m
+
k
3
,
t
m
+
Δ
t
)
.


## 2.1. Choosing the Optimal Bandwidth

The limit relations (10) and (11) imply that there exists a connection between the number of expansion terms M and the bandwidth σ of the Hermite DAF. It is natural to think of the expansion order M as a measure of predefined accuracy so one might ask how to choose σ optimally with respect to a given M. One possible approach is to make the approximation as accurate as possible on the discrete grid. Equation (9) is exact, if
(16)
δ
σ
,
M
(
x
i
-
x
j
)
=
Δ
i
j
h
x
,
where Δ
i
j represents the Kronecker-Δ function. Therefore, for i
=
j, using relations (16) and (9), one obtains
(17)
σ
=
h
x
2
π
∑
m
=
0
M
/
2
‍
1
m
!
(
-
1
4
)
m
h
2
m
(
0
)
.
In the remainder of the paper, the bandwidth is chosen according to (17).


## 2.2. Application

In this section, we develop a space discretization of the GRLW equation using the HDAF. This equation was proposed by Peregrine [4, 5] as a model for small-amplitude long waves on the surface of water in a channel and is represented mathematically as
(18)
u
t
+
u
x
+
ε
u
p
u
x
-
μ
u
x
x
t
=
0
,
(
x
,
t
)
∈
R
×
[
0
T
]
.
To apply the HDAF for solving the GRLW equation, firstly, we present the following notations:
(19)
u
(
x
,
t
)
=
h
x
∑
j
=
1
N
‍
u
j
(
t
)
δ
σ
,
M
(
x
-
x
j
)
,
u
x
(
x
,
t
)
=
h
x
∑
j
=
1
N
‍
u
j
(
t
)
δ
σ
,
M
(
1
)
(
x
-
x
j
)
,
u
x
x
(
x
,
t
)
=
h
x
∑
j
=
1
N
‍
u
j
(
t
)
δ
σ
,
M
(
2
)
(
x
-
x
j
)
,
u
t
(
x
,
t
)
=
h
x
∑
j
=
1
N
‍
u
j
′
(
t
)
δ
σ
,
M
(
x
-
x
j
)
.
We substitute (19) into (18) at collocation point x
=
x
i and we obtain
(20)
h
x
∑
j
=
1
N
‍
u
j
′
(
t
)
δ
σ
,
M
(
x
i
-
x
j
)
+
h
x
∑
j
=
1
N
‍
u
j
(
t
)
δ
σ
,
M
(
1
)
(
x
i
-
x
j
)
+
ε
(
h
x
∑
j
=
1
N
‍
u
j
(
t
)
δ
σ
,
M
(
x
i
-
x
j
)
)
p
×
(
h
x
∑
j
=
1
N
‍
u
j
(
t
)
δ
σ
,
M
(
1
)
(
x
i
-
x
j
)
)
-
μ
h
x
∑
j
=
1
N
‍
u
j
′
(
t
)
δ
σ
,
M
(
2
)
(
x
i
-
x
j
)
=
0
.In order to write (20) in matrix form, we introduce the following matrix and vector notations:
(21)
u
=
[
u
1
,
u
2
,
…
,
u
N
]
T
,
D
(
1
)
=
δ
σ
,
M
(
1
)
(
x
i
-
x
j
)
,
D
(
2
)
=
δ
σ
,
M
(
2
)
(
x
i
-
x
j
)
.
Consequently, we have
(22)
M
u
˙
=
f
(
u
,
t
)
,
where M
=
μ
D
2
-
I, with I the identity matrix, f
(
u
,
t
)
=
D
1
(
(
ε
/
(
p
+
1
)
)
u
p
+
1
+
u
).Commonly used time-stepping schemes, either explicit or implicit, can be easily used along with the HDAF algorithm. In the present study, we consider the classical fourth-order Runge-Kutta method computed as(23)
u
m
+
1
=
u
m
+
1
6
(
k
1
+
2
k
2
+
2
k
3
+
k
4
)
,
where
(24)
k
1
=
Δ
t
f
(
u
m
,
t
m
)
,
k
2
=
Δ
t
f
(
u
m
+
1
2
k
1
,
t
m
+
1
2
Δ
t
)
,
k
3
=
Δ
t
f
(
u
m
+
1
2
k
2
,
t
m
+
1
2
Δ
t
)
,
k
4
=
Δ
t
f
(
u
m
+
k
3
,
t
m
+
Δ
t
)
.


## 3. Numerical Results

The analytical one-soliton solution of the GRLW equation is given by(25)
u
(
x
,
t
)
=
D
sech
2
⁡
[
K
(
x
-
(
1
+
ε
d
)
t
-
x
0
)
]
p
,
where
(26)
D
=
(
p
+
1
)
(
p
+
2
)
d
2
,
K
=
p
2
μ
ε
d
1
+
ε
d
.
Equation (25) represents a solitary wave of amplitude D
p, speed (
1
+
ε
d
), and width K and initially centered at x
0. In this section, we examine the proposed algorithm using different test problems related to the propagation of one soliton and interaction of two soliton solutions for the RLW (p
=
1) and the MRLW (p
=
2). The accuracy of the scheme is measured using the following error norms:
(27)
L
2
=
∥
u
-
u
~
∥
2
=
[
h
∑
j
=
1
N
‍
(
u
j
-
u
~
j
)
2
]
1
/
2
,
L
∞
=
∥
u
-
u
~
∥
∞
=
max
⁡
1
⩽
j
⩽
N
⁡
|
u
j
-
u
~
j
|
,
where u and u
~ represent the exact and approximate solutions, respectively, and h is the minimum distance between any two consecutive points for which the errors are evaluated. In our computational work, we use the collocation points
(28)
{
x
1
=
a
,
…
,
x
i
=
a
+
(
i
-
1
)
h
,
…
,
x
N
=
b
}
,
h
=
|
b
-
a
|
N
-
1
.
The RLW equation possesses infinite conservation laws [30]; the first three are given as follows:
(29)
I
1
=
∫
a
b
‍
u
d
x
,
I
2
=
∫
a
b
‍
(
u
2
+
μ
u
x
2
)
d
x
,
I
3
=
∫
a
b
‍
(
u
3
+
3
u
2
)
d
x
,
related to the mass, momentum, and energy. The quantities I
1, I
2, and I
3 are applied to measure the conservation properties of the collocation scheme, calculated by
(30)
I
1
≈
h
∑
j
=
1
N
‍
u
j
n
,
I
2
≈
h
∑
j
=
1
N
‍
[
(
u
j
n
)
2
+
μ
(
(
u
x
)
j
n
)
2
]
,
I
3
≈
h
∑
j
=
1
N
‍
[
(
u
j
n
)
3
+
3
(
u
j
n
)
2
]
.
Similarly, the first three conservation laws of the MRLW equation are given by
(31)
I
1
=
∫
a
b
‍
u
d
x
,
I
2
=
∫
a
b
‍
(
u
2
+
μ
u
x
2
)
d
x
,
I
3
=
∫
a
b
‍
(
u
4
-
6
μ
u
x
2
)
d
x
,
related to the mass, momentum, and energy. The discrete approximation of the quantities I
1, I
2, and I
3 is obtained by quadrature rule
(32)
I
1
≈
h
∑
j
=
1
N
‍
u
j
n
,
I
2
≈
h
∑
j
=
1
N
‍
[
(
u
j
n
)
2
+
μ
(
(
u
x
)
j
n
)
2
]
,
I
3
≈
h
∑
j
=
1
N
‍
[
(
u
j
n
)
3
-
6
μ
(
(
u
x
)
j
n
)
2
]
.


### 3.1. Propagation of Single Solitons

In this experiment, we consider the propagation of single solitons of the RLW and the MRLW equations. The initial and boundary conditions are extracted from (25) at initial profile, that is, at t
=
0. The values of the parameters used in the numerical experiments are μ
=
1 and ε
=
1, the number of collocation points N
=
257, the time step Δ
t
=
0.01, and the HDAF bandwidth M
=
32. The soliton whose initial condition is extracted from (25) moves to the right across the space interval x
∈
[
-
100,100
] and x
∈
[
-
70,70
] for d
=
0.1 and d
=
0.5, respectively, when the time interval is t
∈
[
0,20
]. We perform various numerical tests and compare our results with those obtained using Sinc numerical methods.


#### 3.1.1. RLW Equation

To illustrate the validity of our scheme, we first consider the case of a single soliton solution for the RLW equation. The analytical values of conservation quantities can be found as(33)
I
1
=
6
d
K
,
I
2
=
12
d
2
K
+
48
K
d
2
μ
5
,
I
3
=
36
d
2
K
+
144
d
3
5
K
,
where K is defined in (26). For d
=
0.1, the analytical values of the invariants are I
1
=
3.979949, I
2
=
0.810462, and I
3
=
2.579007, while for d
=
0.5 and μ
=
1, the analytical values of the invariants become I
1
=
10.392304, I
2
=
11.085125, and I
3
=
43.647680.In Figure2, we plot numerical solutions of the RLW equation at different time levels (a) and the absolute error (b) at time t
=
20. The propagation of one-soliton solution of the RLW equation is in very good agreement with the analytical solution.Propagation of soliton solutions (a) and absolute error (b) att
=
20 of the RLW equation for d
=
0.1, μ
=
1, ε
=
1, and x
∈
[
-
100,100
].
(a)
(b)In the next experiment, we investigate the error normsL
∞ and L
2 together with numerical conservation quantities I
1, I
2, and I
3 obtained with the HDAF method and compare the results with those obtained using the Sinc numerical method (see Tables 1 and 2).Table 1
Invariants and errors for a single soliton of the RLW equation ford
=
0.1
,
μ
=
1
,
ε
=
1 and x
∈
[
-
100,100
].
Time
L
∞
L
2
I
1
I
2
I
3
CPU
HDAF
4
1.1366
E
-
12
2.3426
E
-
12
3.979949
0.810462
2.579007
0.297 s
12
1.2127
E
-
12
3.2739
E
-
12
3.979949
0.810462
2.579007
0.789 s
20
2.3382
E
-
12
5.5084
E
-
12
3.979949
0.810462
2.579007
1.488 s
Sinc
4
4.0054
E
-
13
9.1811
E
-
13
3.979949
0.810462
2.579007
0.492 s
12
1.0984
E
-
12
2.5844
E
-
12
3.979949
0.810462
2.579007
1.604 s
20
2.7039
E
-
12
6.4840
E
-
12
3.979949
0.810462
2.579007
3.493 sTable 2
Invariants and errors for a single soliton of the RLW equation ford
=
0.5
,
μ
=
1
,
ε
=
1 and x
∈
[
-
70,70
].
Time
L
∞
L
2
I
1
I
2
I
3
CPU
HDAF
4
1.1268
E
-
09
1.9530
E
-
09
10.392304
11.085125
43.647680
0.287 s
12
1.1372
E
-
09
1.9856
E
-
09
10.392304
11.085125
43.647680
0.851 s
20
1.1480
E
-
09
2.1011
E
-
09
10.392304
11.085125
43.647680
1.503 s
Sinc
4
1.2939
E
-
10
2.3591
E
-
10
10.392304
11.085125
43.647680
0.438 s
12
2.4847
E
-
10
5.1750
E
-
10
10.392304
11.085125
43.647680
1.625 s
20
5.7125
E
-
10
8.5734
E
-
10
10.392304
11.085125
43.647680
3.442 sIn the case of the HDAF method, ford
=
0.1, one can see clearly that the error norms L
∞ and L
2 are of the magnitude 1
0
-
12 throughout the time interval t
∈
(
0,20
]. Numerical invariants I
1, I
2, and I
3, at a given time t, are very satisfactorily conserved (see Table 1). For d
=
0.5, the soliton has a higher amplitude and a faster speed compared to the latter case where d
=
0.1. The error norms L
∞ and L
2 are of magnitude 1
0
-
9 throughout the time interval t
∈
(
0,20
]. Numerical invariants are satisfactorily conserved (see Table 2).In the case of the Sinc method [31], the error is as accurate as using the HDAF method, and we have conserved the three invariants. However, one can see clearly from Table 2 that our scheme is faster than the Sinc method proposed in [31]. The main reason is that the Sinc method is a global method, while the HDAF is local. In the present numerical experiment, the Sinc method solves N
=
257 systems of nonlinear ODEs while the HDAF method solves M
=
32 systems, locally.In the next experiment, we varyM, the highest polynomial degree involved in the construction of the HDAF, and then investigate the convergence of our approach in terms of the number of grid points N at time t
=
20 with μ
=
1, ε
=
1, and d
=
0.1 on the interval x
∈
[
-
100,100
]. To ensure that the error is mostly dominated by the space discretization, we choose the time step as Δ
t
=
0.001.Figure3 shows that as the number of grid points increases, if M is small, for instance, M
=
4,8, the HDAF method shows a poor convergence. For M
=
16, the convergence of the HDAF method is moderate. When M
=
64,128, the convergence of the HDAF and the Sinc numerical methods is very comparable. The error norms L
∞ and L
2 converge exponentially for the propagation of one-soliton solutions for large values of M. In fact, there are two approximations involved in computing the matrix representation of the differential operator. First, the delta function is approximated by an even Hermite-series with leading Gaussian term (15), which is sufficiently smooth everywhere. Second, the integral in (13) is evaluated by Gauss quadrature approximation (14). Therefore, it is important to choose a good combination between N and M to ensure a fast convergence of the HDAF methodology. The accuracy can be improved by increasing the degree of the basis functions M for an increasing number of grid points N as shown in Figure 3.Errors for the propagation of one-soliton solutions of the RLW equation (a)L
2 and (b) L
∞ at t
=
20 of the RLW equation for d
=
0.1, μ
=
1, ε
=
1, and x
∈
[
-
100,100
].
(a)
(b)


#### 3.1.2. MRLW Equation

We model the motion of one-soliton solution of the MRLW equation withd
=
0.1, μ
=
1, and ε
=
1; the number of collocation points N
=
257; the time step Δ
t
=
0.01; and the HDAF bandwidth M
=
32. The soliton whose initial condition is extracted from (25) moves to the right across the space interval x
∈
[
-
100,100
]. The analytical computation of the invariants can be written as
(34)
I
1
=
π
6
d
K
,
I
2
=
12
d
(
K
μ
3
+
1
K
)
,
I
3
=
8
d
(
6
d
-
3
K
2
μ
)
K
.
For d
=
0.1, μ
=
1, and ε
=
1, the analytical values of the invariants are I
1
=
8.070897, I
2
=
4.100554, and I
3
=
0.868352.In Figure4, we plot numerical solutions of the MRLW equation at different time levels (a) and the absolute error (b) at time t
=
20. The propagation of the one-soliton solution of the MRLW equation is in very good agreement with the analytical solution.Propagation of soliton solutions (a) and absolute error (b) att
=
20 of the MRLW equation for d
=
0.1, μ
=
1, ε
=
1, and x
∈
[
-
100,100
].
(a)
(b)The error normsL
∞ and L
2 and numerical conservation quantities I
1, I
2, and I
3 are shown in Table 3. In the case of the HDAF method, it is observed that, throughout the simulation, the error norms L
∞ remain less than 1.05
×
1
0
-
9 and L
2 remains less than 3.08
×
1
0
-
10. Numerical invariants I
1, I
2, and I
3 at a given time t duplicate analytical values. In the case of Sinc method, we obtain the same order of accuracy and one can see that the HDAF method is more efficient than the Sinc method. In the present numerical experiment, the Sinc method solves N
=
257 systems of nonlinear ODEs while the HDAF method solves locally M
=
32 systems. This is more problematic for higher dimension PDEs. For instance, for a two-dimensional RLW equation, we need to solve N
2
=
257
×
257 systems of nonlinear ODEs with the Sinc numerical method while only M
2
=
32
×
32 systems are required to be solved with the HDAF method.Table 3
Invariants and errors for a single soliton of the MRWL equation ford
=
0.1
,
μ
=
1 and ε
=
1.
Time
L
∞
L
2
I
1
I
2
I
3
CPU
HDAF
4
3.4874
E
-
10
9.1053
E
-
10
8.070897
4.100554
0.868352
0.181 s
12
3.5215
E
-
10
1.0913
E
-
09
8.070897
4.100554
0.868352
0.387 s
20
3.0752
E
-
10
1.0499
E
-
09
8.070897
4.100554
0.868352
1.624 s
Sinc
4
1.6085
E
-
10
7.6134
E
-
10
8.070897
4.100554
0.868352
0.922 s
12
1.9117
E
-
10
5.3271
E
-
10
8.070897
4.100554
0.868352
1.765 s
20
2.0221
E
-
10
8.4638
E
-
10
8.070897
4.100554
0.868352
3.722 s


### 3.2. Interaction of Two Solitons

In our second experiment, we study the interaction of two soliton solutions of the RLW and MRLW equations having different amplitudes and traveling in the same direction. The analytical two soliton solutions of the GRLW equation are given by(35)
u
(
x
,
t
)
=
∑
i
2
‍
D
i
sech
2
⁡
[
K
i
(
x
-
(
1
+
ε
d
i
)
t
-
x
i
)
]
p
,
where
(36)
D
i
=
(
p
+
1
)
(
p
+
2
)
d
i
2
,
K
i
=
p
2
μ
ε
d
i
1
+
ε
d
i
.
The analytical computations of conservation quantities of the RLW and MRLW equations can be found as
(37)
I
1
=
∑
i
=
1
2
‍
6
d
i
K
i
,
I
2
=
∑
i
=
1
2
‍
12
d
i
2
K
i
+
48
K
i
d
i
2
μ
5
,
I
3
=
∑
i
=
1
2
‍
36
d
i
2
K
i
+
144
d
i
3
5
K
i
,
(38)
I
1
=
∑
i
=
1
2
‍
π
6
d
i
K
i
,
I
2
=
∑
i
=
1
2
‍
12
d
i
(
K
i
μ
3
+
1
K
i
)
,
I
3
=
∑
i
=
1
2
‍
8
d
i
(
6
d
i
-
3
K
i
2
μ
)
K
i
,
respectively. The initial conditions of the RLW and MRLW are extracted from (35).We choosed
1
=
2, d
2
=
0.75, x
1
=
20, x
2
=
40, μ
=
1, and ε
=
1 with interval [
0,120
]. Figure 5 shows the plot of these solitary waves’ interaction at different time levels. In both cases, we observe that the faster pulse interacts with and emerges ahead of the lower pulse with the shape and velocity of each soliton retained. Analytical values of the invariants for the RLW and MRLW equations are given by I
1
=
39.786181, I
2
=
144.3373671, and I
3
=
960.7366406 and I
1
=
22.753426, I
2
=
47.472854, and I
3
=
209.815585, respectively. In Table 4 we report the three invariants for both RLW and MRLW equations and the two errors. Numerical checks on the conservation mass I
1, momentum I
2, and energy I
3 show that the three quantities remain reasonably constant as the solitary waves evolve in time. We also obtain very accurate errors over longer time periods.Table 4
Invariants and errors for two solitons of RLW and MRLW equation ford
1
=
2
,
d
2
=
0.75
,
x
1
=
20
,
x
2
=
40
,
μ
=
1 and ε
=
1 with interval [
0,120
].
Time
L
∞
L
2
I
1
I
2
I
3
CPU
RLW
4
1.543813
E
-
10
6.353135
E
-
10
39.786172
144.337366
960.736641
0.181 s
20
1.764622
E
-
10
6.646963
E
-
10
39.786196
144.337387
960.736647
0.387 s
30
2.357418
E
-
10
8.155187
E
-
10
39.786192
144.337367
960.736649
1.624 s
MRLW
4
9.651607
E
-
11
4.297970
E
-
10
22.753427
47.472855
209.815556
0.183 s
20
1.004345
E
-
10
5.112434
E
-
10
22.753428
47.472857
209.815538
0.385 s
30
1.585590
E
-
10
6.170251
E
-
10
22.753427
47.472858
209.815566
1.691 sInteraction of two soliton solutions of the RLW (a) and MRLW (b) equations.K
1
=
0.4, K
2
=
0.1, x
1
=
20, x
2
=
40, N
=
257, Δ
t
=
0.001, and x
∈
[
0,120
].
(a)
(b)Propagation of the single solitary wave and two-soliton interaction are simulated well with the proposed algorithms and conservation invariants do not change much during the computer run. Thus, HDAF functions can be used to construct accurate and efficient numerical solutions of the RLW and MRLW equations.


### 3.3. Maxwellian Initial Condition

In final series of numerical experiments, we examined the development of the Maxwellian initial condition(39)
u
(
x
,
0
)
=
exp
⁡
[
-
(
x
-
10
)
2
]
into multiple solitary waves. We apply it to the problem for different cases: μ
=
0.05, 0.01, 0.004.Forμ
=
0.05 the Maxwellian initial condition develops into a single solitary wave plus a small developed oscillating tail as shown in Figure 6(b) at time T
=
12. A bit of radiation is observed between the upstream boundary and the main peak wave. The values of the quantities I
1, I
2, and I
3 are given in Table 5. These three quantities remain relatively constant. The maximum changes are less than 1
0
-
6 for each quantity.Table 5
The invariants of MRLW equation using the Maxwellian initial condition for different values ofμ = 0.05, 0.01, 0.004 and T
=
12.
μ
Time
I
1
I
2
I
3
0.05
0
1.772453
1.253314
4.783269
4
1.772453
1.253314
4.783269
12
1.772453
1.253313
4.783269
20
1.772453
1.253312
4.783268
0.01
0
1.772453
1.253315
4.783269
4
1.772453
1.253318
4.783269
12
1.772453
1.253309
4.783269
20
1.772453
1.253308
4.783269
0.004
0
1.772453
1.253314
4.783269
4
1.772453
1.253307
4.783269
12
1.772453
1.253308
4.783269
20
1.772453
1.253302
4.783269Development of solitary waves using the Maxwellian initial condition for different values ofμ
=
0.05, 0.01, 0.004 and T
=
12.
(a)
(b)
(c)
(d)Forμ
=
0.01 the Maxwellian pulse breaks up into a train of at least three solitary waves as shown in Figure 6(c) at time T
=
12. A bit of radiation is observed between the upstream boundary and the main peak wave. In Table 5, we report the values of the quantities I
1, I
2, and I
3. Clearly, these quantities are satisfactorily constant. The maximum changes are less than 1
0
-
6.Finally, forμ
=
0.004 the Maxwellian initial condition has decayed into at least four stable solitary waves as shown in Figure 6(d). The peaks of the well developed wave lie on a straight line so that their velocities are linearly dependent on their amplitudes and we observe a small oscillating tail appearing behind the last wave as shown in Figure 6. Moreover, the total number of the solitary waves which are generated from the Maxwellian initial condition according to the results obtained from the numerical scheme in test problem follows from a relationship between the number of solitary waves generated and the parameter μ. We notice that the smaller μ becomes, the higher the number of solitons generated is. Raslan [17] showed that this relationship is approximately N
≈
[
1
/
μ
4
] for the RLW equation and N
≈
[
1
/
μ
5
] for the MRLW equation. However, there is no proof available on these relationships.


## 3.1. Propagation of Single Solitons

In this experiment, we consider the propagation of single solitons of the RLW and the MRLW equations. The initial and boundary conditions are extracted from (25) at initial profile, that is, at t
=
0. The values of the parameters used in the numerical experiments are μ
=
1 and ε
=
1, the number of collocation points N
=
257, the time step Δ
t
=
0.01, and the HDAF bandwidth M
=
32. The soliton whose initial condition is extracted from (25) moves to the right across the space interval x
∈
[
-
100,100
] and x
∈
[
-
70,70
] for d
=
0.1 and d
=
0.5, respectively, when the time interval is t
∈
[
0,20
]. We perform various numerical tests and compare our results with those obtained using Sinc numerical methods.


### 3.1.1. RLW Equation

To illustrate the validity of our scheme, we first consider the case of a single soliton solution for the RLW equation. The analytical values of conservation quantities can be found as(33)
I
1
=
6
d
K
,
I
2
=
12
d
2
K
+
48
K
d
2
μ
5
,
I
3
=
36
d
2
K
+
144
d
3
5
K
,
where K is defined in (26). For d
=
0.1, the analytical values of the invariants are I
1
=
3.979949, I
2
=
0.810462, and I
3
=
2.579007, while for d
=
0.5 and μ
=
1, the analytical values of the invariants become I
1
=
10.392304, I
2
=
11.085125, and I
3
=
43.647680.In Figure2, we plot numerical solutions of the RLW equation at different time levels (a) and the absolute error (b) at time t
=
20. The propagation of one-soliton solution of the RLW equation is in very good agreement with the analytical solution.Propagation of soliton solutions (a) and absolute error (b) att
=
20 of the RLW equation for d
=
0.1, μ
=
1, ε
=
1, and x
∈
[
-
100,100
].
(a)
(b)In the next experiment, we investigate the error normsL
∞ and L
2 together with numerical conservation quantities I
1, I
2, and I
3 obtained with the HDAF method and compare the results with those obtained using the Sinc numerical method (see Tables 1 and 2).Table 1
Invariants and errors for a single soliton of the RLW equation ford
=
0.1
,
μ
=
1
,
ε
=
1 and x
∈
[
-
100,100
].
Time
L
∞
L
2
I
1
I
2
I
3
CPU
HDAF
4
1.1366
E
-
12
2.3426
E
-
12
3.979949
0.810462
2.579007
0.297 s
12
1.2127
E
-
12
3.2739
E
-
12
3.979949
0.810462
2.579007
0.789 s
20
2.3382
E
-
12
5.5084
E
-
12
3.979949
0.810462
2.579007
1.488 s
Sinc
4
4.0054
E
-
13
9.1811
E
-
13
3.979949
0.810462
2.579007
0.492 s
12
1.0984
E
-
12
2.5844
E
-
12
3.979949
0.810462
2.579007
1.604 s
20
2.7039
E
-
12
6.4840
E
-
12
3.979949
0.810462
2.579007
3.493 sTable 2
Invariants and errors for a single soliton of the RLW equation ford
=
0.5
,
μ
=
1
,
ε
=
1 and x
∈
[
-
70,70
].
Time
L
∞
L
2
I
1
I
2
I
3
CPU
HDAF
4
1.1268
E
-
09
1.9530
E
-
09
10.392304
11.085125
43.647680
0.287 s
12
1.1372
E
-
09
1.9856
E
-
09
10.392304
11.085125
43.647680
0.851 s
20
1.1480
E
-
09
2.1011
E
-
09
10.392304
11.085125
43.647680
1.503 s
Sinc
4
1.2939
E
-
10
2.3591
E
-
10
10.392304
11.085125
43.647680
0.438 s
12
2.4847
E
-
10
5.1750
E
-
10
10.392304
11.085125
43.647680
1.625 s
20
5.7125
E
-
10
8.5734
E
-
10
10.392304
11.085125
43.647680
3.442 sIn the case of the HDAF method, ford
=
0.1, one can see clearly that the error norms L
∞ and L
2 are of the magnitude 1
0
-
12 throughout the time interval t
∈
(
0,20
]. Numerical invariants I
1, I
2, and I
3, at a given time t, are very satisfactorily conserved (see Table 1). For d
=
0.5, the soliton has a higher amplitude and a faster speed compared to the latter case where d
=
0.1. The error norms L
∞ and L
2 are of magnitude 1
0
-
9 throughout the time interval t
∈
(
0,20
]. Numerical invariants are satisfactorily conserved (see Table 2).In the case of the Sinc method [31], the error is as accurate as using the HDAF method, and we have conserved the three invariants. However, one can see clearly from Table 2 that our scheme is faster than the Sinc method proposed in [31]. The main reason is that the Sinc method is a global method, while the HDAF is local. In the present numerical experiment, the Sinc method solves N
=
257 systems of nonlinear ODEs while the HDAF method solves M
=
32 systems, locally.In the next experiment, we varyM, the highest polynomial degree involved in the construction of the HDAF, and then investigate the convergence of our approach in terms of the number of grid points N at time t
=
20 with μ
=
1, ε
=
1, and d
=
0.1 on the interval x
∈
[
-
100,100
]. To ensure that the error is mostly dominated by the space discretization, we choose the time step as Δ
t
=
0.001.Figure3 shows that as the number of grid points increases, if M is small, for instance, M
=
4,8, the HDAF method shows a poor convergence. For M
=
16, the convergence of the HDAF method is moderate. When M
=
64,128, the convergence of the HDAF and the Sinc numerical methods is very comparable. The error norms L
∞ and L
2 converge exponentially for the propagation of one-soliton solutions for large values of M. In fact, there are two approximations involved in computing the matrix representation of the differential operator. First, the delta function is approximated by an even Hermite-series with leading Gaussian term (15), which is sufficiently smooth everywhere. Second, the integral in (13) is evaluated by Gauss quadrature approximation (14). Therefore, it is important to choose a good combination between N and M to ensure a fast convergence of the HDAF methodology. The accuracy can be improved by increasing the degree of the basis functions M for an increasing number of grid points N as shown in Figure 3.Errors for the propagation of one-soliton solutions of the RLW equation (a)L
2 and (b) L
∞ at t
=
20 of the RLW equation for d
=
0.1, μ
=
1, ε
=
1, and x
∈
[
-
100,100
].
(a)
(b)


### 3.1.2. MRLW Equation

We model the motion of one-soliton solution of the MRLW equation withd
=
0.1, μ
=
1, and ε
=
1; the number of collocation points N
=
257; the time step Δ
t
=
0.01; and the HDAF bandwidth M
=
32. The soliton whose initial condition is extracted from (25) moves to the right across the space interval x
∈
[
-
100,100
]. The analytical computation of the invariants can be written as
(34)
I
1
=
π
6
d
K
,
I
2
=
12
d
(
K
μ
3
+
1
K
)
,
I
3
=
8
d
(
6
d
-
3
K
2
μ
)
K
.
For d
=
0.1, μ
=
1, and ε
=
1, the analytical values of the invariants are I
1
=
8.070897, I
2
=
4.100554, and I
3
=
0.868352.In Figure4, we plot numerical solutions of the MRLW equation at different time levels (a) and the absolute error (b) at time t
=
20. The propagation of the one-soliton solution of the MRLW equation is in very good agreement with the analytical solution.Propagation of soliton solutions (a) and absolute error (b) att
=
20 of the MRLW equation for d
=
0.1, μ
=
1, ε
=
1, and x
∈
[
-
100,100
].
(a)
(b)The error normsL
∞ and L
2 and numerical conservation quantities I
1, I
2, and I
3 are shown in Table 3. In the case of the HDAF method, it is observed that, throughout the simulation, the error norms L
∞ remain less than 1.05
×
1
0
-
9 and L
2 remains less than 3.08
×
1
0
-
10. Numerical invariants I
1, I
2, and I
3 at a given time t duplicate analytical values. In the case of Sinc method, we obtain the same order of accuracy and one can see that the HDAF method is more efficient than the Sinc method. In the present numerical experiment, the Sinc method solves N
=
257 systems of nonlinear ODEs while the HDAF method solves locally M
=
32 systems. This is more problematic for higher dimension PDEs. For instance, for a two-dimensional RLW equation, we need to solve N
2
=
257
×
257 systems of nonlinear ODEs with the Sinc numerical method while only M
2
=
32
×
32 systems are required to be solved with the HDAF method.Table 3
Invariants and errors for a single soliton of the MRWL equation ford
=
0.1
,
μ
=
1 and ε
=
1.
Time
L
∞
L
2
I
1
I
2
I
3
CPU
HDAF
4
3.4874
E
-
10
9.1053
E
-
10
8.070897
4.100554
0.868352
0.181 s
12
3.5215
E
-
10
1.0913
E
-
09
8.070897
4.100554
0.868352
0.387 s
20
3.0752
E
-
10
1.0499
E
-
09
8.070897
4.100554
0.868352
1.624 s
Sinc
4
1.6085
E
-
10
7.6134
E
-
10
8.070897
4.100554
0.868352
0.922 s
12
1.9117
E
-
10
5.3271
E
-
10
8.070897
4.100554
0.868352
1.765 s
20
2.0221
E
-
10
8.4638
E
-
10
8.070897
4.100554
0.868352
3.722 s


## 3.1.1. RLW Equation

To illustrate the validity of our scheme, we first consider the case of a single soliton solution for the RLW equation. The analytical values of conservation quantities can be found as(33)
I
1
=
6
d
K
,
I
2
=
12
d
2
K
+
48
K
d
2
μ
5
,
I
3
=
36
d
2
K
+
144
d
3
5
K
,
where K is defined in (26). For d
=
0.1, the analytical values of the invariants are I
1
=
3.979949, I
2
=
0.810462, and I
3
=
2.579007, while for d
=
0.5 and μ
=
1, the analytical values of the invariants become I
1
=
10.392304, I
2
=
11.085125, and I
3
=
43.647680.In Figure2, we plot numerical solutions of the RLW equation at different time levels (a) and the absolute error (b) at time t
=
20. The propagation of one-soliton solution of the RLW equation is in very good agreement with the analytical solution.Propagation of soliton solutions (a) and absolute error (b) att
=
20 of the RLW equation for d
=
0.1, μ
=
1, ε
=
1, and x
∈
[
-
100,100
].
(a)
(b)In the next experiment, we investigate the error normsL
∞ and L
2 together with numerical conservation quantities I
1, I
2, and I
3 obtained with the HDAF method and compare the results with those obtained using the Sinc numerical method (see Tables 1 and 2).Table 1
Invariants and errors for a single soliton of the RLW equation ford
=
0.1
,
μ
=
1
,
ε
=
1 and x
∈
[
-
100,100
].
Time
L
∞
L
2
I
1
I
2
I
3
CPU
HDAF
4
1.1366
E
-
12
2.3426
E
-
12
3.979949
0.810462
2.579007
0.297 s
12
1.2127
E
-
12
3.2739
E
-
12
3.979949
0.810462
2.579007
0.789 s
20
2.3382
E
-
12
5.5084
E
-
12
3.979949
0.810462
2.579007
1.488 s
Sinc
4
4.0054
E
-
13
9.1811
E
-
13
3.979949
0.810462
2.579007
0.492 s
12
1.0984
E
-
12
2.5844
E
-
12
3.979949
0.810462
2.579007
1.604 s
20
2.7039
E
-
12
6.4840
E
-
12
3.979949
0.810462
2.579007
3.493 sTable 2
Invariants and errors for a single soliton of the RLW equation ford
=
0.5
,
μ
=
1
,
ε
=
1 and x
∈
[
-
70,70
].
Time
L
∞
L
2
I
1
I
2
I
3
CPU
HDAF
4
1.1268
E
-
09
1.9530
E
-
09
10.392304
11.085125
43.647680
0.287 s
12
1.1372
E
-
09
1.9856
E
-
09
10.392304
11.085125
43.647680
0.851 s
20
1.1480
E
-
09
2.1011
E
-
09
10.392304
11.085125
43.647680
1.503 s
Sinc
4
1.2939
E
-
10
2.3591
E
-
10
10.392304
11.085125
43.647680
0.438 s
12
2.4847
E
-
10
5.1750
E
-
10
10.392304
11.085125
43.647680
1.625 s
20
5.7125
E
-
10
8.5734
E
-
10
10.392304
11.085125
43.647680
3.442 sIn the case of the HDAF method, ford
=
0.1, one can see clearly that the error norms L
∞ and L
2 are of the magnitude 1
0
-
12 throughout the time interval t
∈
(
0,20
]. Numerical invariants I
1, I
2, and I
3, at a given time t, are very satisfactorily conserved (see Table 1). For d
=
0.5, the soliton has a higher amplitude and a faster speed compared to the latter case where d
=
0.1. The error norms L
∞ and L
2 are of magnitude 1
0
-
9 throughout the time interval t
∈
(
0,20
]. Numerical invariants are satisfactorily conserved (see Table 2).In the case of the Sinc method [31], the error is as accurate as using the HDAF method, and we have conserved the three invariants. However, one can see clearly from Table 2 that our scheme is faster than the Sinc method proposed in [31]. The main reason is that the Sinc method is a global method, while the HDAF is local. In the present numerical experiment, the Sinc method solves N
=
257 systems of nonlinear ODEs while the HDAF method solves M
=
32 systems, locally.In the next experiment, we varyM, the highest polynomial degree involved in the construction of the HDAF, and then investigate the convergence of our approach in terms of the number of grid points N at time t
=
20 with μ
=
1, ε
=
1, and d
=
0.1 on the interval x
∈
[
-
100,100
]. To ensure that the error is mostly dominated by the space discretization, we choose the time step as Δ
t
=
0.001.Figure3 shows that as the number of grid points increases, if M is small, for instance, M
=
4,8, the HDAF method shows a poor convergence. For M
=
16, the convergence of the HDAF method is moderate. When M
=
64,128, the convergence of the HDAF and the Sinc numerical methods is very comparable. The error norms L
∞ and L
2 converge exponentially for the propagation of one-soliton solutions for large values of M. In fact, there are two approximations involved in computing the matrix representation of the differential operator. First, the delta function is approximated by an even Hermite-series with leading Gaussian term (15), which is sufficiently smooth everywhere. Second, the integral in (13) is evaluated by Gauss quadrature approximation (14). Therefore, it is important to choose a good combination between N and M to ensure a fast convergence of the HDAF methodology. The accuracy can be improved by increasing the degree of the basis functions M for an increasing number of grid points N as shown in Figure 3.Errors for the propagation of one-soliton solutions of the RLW equation (a)L
2 and (b) L
∞ at t
=
20 of the RLW equation for d
=
0.1, μ
=
1, ε
=
1, and x
∈
[
-
100,100
].
(a)
(b)


## 3.1.2. MRLW Equation

We model the motion of one-soliton solution of the MRLW equation withd
=
0.1, μ
=
1, and ε
=
1; the number of collocation points N
=
257; the time step Δ
t
=
0.01; and the HDAF bandwidth M
=
32. The soliton whose initial condition is extracted from (25) moves to the right across the space interval x
∈
[
-
100,100
]. The analytical computation of the invariants can be written as
(34)
I
1
=
π
6
d
K
,
I
2
=
12
d
(
K
μ
3
+
1
K
)
,
I
3
=
8
d
(
6
d
-
3
K
2
μ
)
K
.
For d
=
0.1, μ
=
1, and ε
=
1, the analytical values of the invariants are I
1
=
8.070897, I
2
=
4.100554, and I
3
=
0.868352.In Figure4, we plot numerical solutions of the MRLW equation at different time levels (a) and the absolute error (b) at time t
=
20. The propagation of the one-soliton solution of the MRLW equation is in very good agreement with the analytical solution.Propagation of soliton solutions (a) and absolute error (b) att
=
20 of the MRLW equation for d
=
0.1, μ
=
1, ε
=
1, and x
∈
[
-
100,100
].
(a)
(b)The error normsL
∞ and L
2 and numerical conservation quantities I
1, I
2, and I
3 are shown in Table 3. In the case of the HDAF method, it is observed that, throughout the simulation, the error norms L
∞ remain less than 1.05
×
1
0
-
9 and L
2 remains less than 3.08
×
1
0
-
10. Numerical invariants I
1, I
2, and I
3 at a given time t duplicate analytical values. In the case of Sinc method, we obtain the same order of accuracy and one can see that the HDAF method is more efficient than the Sinc method. In the present numerical experiment, the Sinc method solves N
=
257 systems of nonlinear ODEs while the HDAF method solves locally M
=
32 systems. This is more problematic for higher dimension PDEs. For instance, for a two-dimensional RLW equation, we need to solve N
2
=
257
×
257 systems of nonlinear ODEs with the Sinc numerical method while only M
2
=
32
×
32 systems are required to be solved with the HDAF method.Table 3
Invariants and errors for a single soliton of the MRWL equation ford
=
0.1
,
μ
=
1 and ε
=
1.
Time
L
∞
L
2
I
1
I
2
I
3
CPU
HDAF
4
3.4874
E
-
10
9.1053
E
-
10
8.070897
4.100554
0.868352
0.181 s
12
3.5215
E
-
10
1.0913
E
-
09
8.070897
4.100554
0.868352
0.387 s
20
3.0752
E
-
10
1.0499
E
-
09
8.070897
4.100554
0.868352
1.624 s
Sinc
4
1.6085
E
-
10
7.6134
E
-
10
8.070897
4.100554
0.868352
0.922 s
12
1.9117
E
-
10
5.3271
E
-
10
8.070897
4.100554
0.868352
1.765 s
20
2.0221
E
-
10
8.4638
E
-
10
8.070897
4.100554
0.868352
3.722 s


## 3.2. Interaction of Two Solitons

In our second experiment, we study the interaction of two soliton solutions of the RLW and MRLW equations having different amplitudes and traveling in the same direction. The analytical two soliton solutions of the GRLW equation are given by(35)
u
(
x
,
t
)
=
∑
i
2
‍
D
i
sech
2
⁡
[
K
i
(
x
-
(
1
+
ε
d
i
)
t
-
x
i
)
]
p
,
where
(36)
D
i
=
(
p
+
1
)
(
p
+
2
)
d
i
2
,
K
i
=
p
2
μ
ε
d
i
1
+
ε
d
i
.
The analytical computations of conservation quantities of the RLW and MRLW equations can be found as
(37)
I
1
=
∑
i
=
1
2
‍
6
d
i
K
i
,
I
2
=
∑
i
=
1
2
‍
12
d
i
2
K
i
+
48
K
i
d
i
2
μ
5
,
I
3
=
∑
i
=
1
2
‍
36
d
i
2
K
i
+
144
d
i
3
5
K
i
,
(38)
I
1
=
∑
i
=
1
2
‍
π
6
d
i
K
i
,
I
2
=
∑
i
=
1
2
‍
12
d
i
(
K
i
μ
3
+
1
K
i
)
,
I
3
=
∑
i
=
1
2
‍
8
d
i
(
6
d
i
-
3
K
i
2
μ
)
K
i
,
respectively. The initial conditions of the RLW and MRLW are extracted from (35).We choosed
1
=
2, d
2
=
0.75, x
1
=
20, x
2
=
40, μ
=
1, and ε
=
1 with interval [
0,120
]. Figure 5 shows the plot of these solitary waves’ interaction at different time levels. In both cases, we observe that the faster pulse interacts with and emerges ahead of the lower pulse with the shape and velocity of each soliton retained. Analytical values of the invariants for the RLW and MRLW equations are given by I
1
=
39.786181, I
2
=
144.3373671, and I
3
=
960.7366406 and I
1
=
22.753426, I
2
=
47.472854, and I
3
=
209.815585, respectively. In Table 4 we report the three invariants for both RLW and MRLW equations and the two errors. Numerical checks on the conservation mass I
1, momentum I
2, and energy I
3 show that the three quantities remain reasonably constant as the solitary waves evolve in time. We also obtain very accurate errors over longer time periods.Table 4
Invariants and errors for two solitons of RLW and MRLW equation ford
1
=
2
,
d
2
=
0.75
,
x
1
=
20
,
x
2
=
40
,
μ
=
1 and ε
=
1 with interval [
0,120
].
Time
L
∞
L
2
I
1
I
2
I
3
CPU
RLW
4
1.543813
E
-
10
6.353135
E
-
10
39.786172
144.337366
960.736641
0.181 s
20
1.764622
E
-
10
6.646963
E
-
10
39.786196
144.337387
960.736647
0.387 s
30
2.357418
E
-
10
8.155187
E
-
10
39.786192
144.337367
960.736649
1.624 s
MRLW
4
9.651607
E
-
11
4.297970
E
-
10
22.753427
47.472855
209.815556
0.183 s
20
1.004345
E
-
10
5.112434
E
-
10
22.753428
47.472857
209.815538
0.385 s
30
1.585590
E
-
10
6.170251
E
-
10
22.753427
47.472858
209.815566
1.691 sInteraction of two soliton solutions of the RLW (a) and MRLW (b) equations.K
1
=
0.4, K
2
=
0.1, x
1
=
20, x
2
=
40, N
=
257, Δ
t
=
0.001, and x
∈
[
0,120
].
(a)
(b)Propagation of the single solitary wave and two-soliton interaction are simulated well with the proposed algorithms and conservation invariants do not change much during the computer run. Thus, HDAF functions can be used to construct accurate and efficient numerical solutions of the RLW and MRLW equations.


## 3.3. Maxwellian Initial Condition

In final series of numerical experiments, we examined the development of the Maxwellian initial condition(39)
u
(
x
,
0
)
=
exp
⁡
[
-
(
x
-
10
)
2
]
into multiple solitary waves. We apply it to the problem for different cases: μ
=
0.05, 0.01, 0.004.Forμ
=
0.05 the Maxwellian initial condition develops into a single solitary wave plus a small developed oscillating tail as shown in Figure 6(b) at time T
=
12. A bit of radiation is observed between the upstream boundary and the main peak wave. The values of the quantities I
1, I
2, and I
3 are given in Table 5. These three quantities remain relatively constant. The maximum changes are less than 1
0
-
6 for each quantity.Table 5
The invariants of MRLW equation using the Maxwellian initial condition for different values ofμ = 0.05, 0.01, 0.004 and T
=
12.
μ
Time
I
1
I
2
I
3
0.05
0
1.772453
1.253314
4.783269
4
1.772453
1.253314
4.783269
12
1.772453
1.253313
4.783269
20
1.772453
1.253312
4.783268
0.01
0
1.772453
1.253315
4.783269
4
1.772453
1.253318
4.783269
12
1.772453
1.253309
4.783269
20
1.772453
1.253308
4.783269
0.004
0
1.772453
1.253314
4.783269
4
1.772453
1.253307
4.783269
12
1.772453
1.253308
4.783269
20
1.772453
1.253302
4.783269Development of solitary waves using the Maxwellian initial condition for different values ofμ
=
0.05, 0.01, 0.004 and T
=
12.
(a)
(b)
(c)
(d)Forμ
=
0.01 the Maxwellian pulse breaks up into a train of at least three solitary waves as shown in Figure 6(c) at time T
=
12. A bit of radiation is observed between the upstream boundary and the main peak wave. In Table 5, we report the values of the quantities I
1, I
2, and I
3. Clearly, these quantities are satisfactorily constant. The maximum changes are less than 1
0
-
6.Finally, forμ
=
0.004 the Maxwellian initial condition has decayed into at least four stable solitary waves as shown in Figure 6(d). The peaks of the well developed wave lie on a straight line so that their velocities are linearly dependent on their amplitudes and we observe a small oscillating tail appearing behind the last wave as shown in Figure 6. Moreover, the total number of the solitary waves which are generated from the Maxwellian initial condition according to the results obtained from the numerical scheme in test problem follows from a relationship between the number of solitary waves generated and the parameter μ. We notice that the smaller μ becomes, the higher the number of solitons generated is. Raslan [17] showed that this relationship is approximately N
≈
[
1
/
μ
4
] for the RLW equation and N
≈
[
1
/
μ
5
] for the MRLW equation. However, there is no proof available on these relationships.


## 4. Conclusion

A numerical technique based on the Hermite distributed approximating functional method has been presented for the numerical solution of the GRLW equation. The efficiency of the method is tested on the problems of propagation of the single soliton and the interaction of two solitons for the RLW and MRLW equations. The accuracy is examined in terms of theL
∞ and L
2 error norms. Numerical results have illustrated the high accuracy, efficiency, and exponential convergence of the HDAF method. We also investigate the conservation quantities I
1, I
2, and I
3 which were satisfactorily maintained. Furthermore, a Maxwellian initial condition has been used and a relation between μ and the number of waves was explored; results show the interaction of three and four solitons, with conservation laws adhered to satisfactorily.



---

*Source: 178024-2014-08-12.xml*