# Petrov-Galerkin Method for the Coupled Schrödinger-KdV Equation

Authors: M. S. Ismail; Farida M. Mosally; Khadeejah M. Alamoudi
Journal: Abstract and Applied Analysis
 (2014)

---

## Abstract

Petrov-Galerkin method is used to derive a numerical scheme for the coupled Schrödinger-KdV (SKdV) equations, where we have used the cubic B-splines as a test functions and a linear B-splines as a trial functions. Product approximation technique is used to deal with the nonlinear terms. An implicit midpoint rule and the Runge-Kutta method of fourth-order (RK4) are used to discretize in time. A block nonlinear pentadiagonal system is obtained. We solve this system by the fixed point method. The resulting scheme has a fourth-order accuracy in space direction and second-order in time direction in case of the implicit midpoint rule and it is unconditionally stable by von Neumann method. Using the RK4 method the scheme will be linear and fourth-order in time and space directions, and it is also conditionally stable. The exact soliton solution and the conserved quantities are used to assess the accuracy and to show the robustness and the efficiency of the proposed schemes.


---

## Body

## 1. Introduction

The coupled nonlinear Schrödinger-KdV equation has attracted extensive interest in physics and mathematics. The nonlinear Schrödinger-KdV equation describing the nonlinear dynamics.Describing the nonlinear dynamics of one-dimensional Langmuir and ion-acoustic waves in a system of coordinates moving at the ion-acoustic speed has the following form [1, 2]:
(1)
i
ϵ
u
t
+
3
2
u
x
x
-
1
2
u
v
=
0
,
(
x
,
t
)
∈
(
a
,
b
)
×
(
0
,
T
)
v
t
+
1
2
v
x
x
x
+
1
2
(
|
u
|
2
+
v
2
)
x
=
0
and the initial condition
(2)
u
(
x
,
0
)
=
-
6
5
3
α
tanh
⁡
(
α
10
x
)
sech
⁡
(
α
10
x
)
exp
⁡
{
i
α
[
-
ϵ
x
3
]
}
v
(
x
,
0
)
=
-
9
5
α
sec
h
2
(
α
10
x
)
,
where α and ϵ are constants, together with the homogenous Dirichlet boundary conditions
(3)
u
(
x
l
,
t
)
=
u
(
x
r
,
t
)
=
0
v
(
x
l
,
t
)
=
v
(
x
r
,
t
)
=
0
.Analytical and numerical methods are very important tools to understand the feature and the behavior of the nonlinear wave equations. Many numerical methods have been used to solve numerically the single nonlinear Schrödinger and the single KdV equation using finite element and finite difference methods [3–6]. Ismail et al. [7–12] solved numerically the coupled nonlinear Schrodinger equation and the coupled KdV equation using the finite difference and finite element methods. Recently, Bhatt and Khaliq [13] solve the coupled nonlinear Schrodinger using high order exponential differencing method.In this work we are going to study numerical the SKdV equation. Analytical solutions for this system using different methods are given in [14–16]. For numerical methods for the SKdV equation the literature is not rich. Appert and Vaclavik [17] solved the SKdV equations using a Crank-Nicolson scheme. A very good numerical works by Bai and Zhang [1, 2], they have solved this system using a split-step quadratic B-spline finite element method.The exact solution of coupled Schrödinger-KdV equations (1) is
(4)
u
(
x
,
t
)
=
-
6
5
3
α
tanh
ξ
cosh
⁡
ξ
exp
⁡
{
i
α
[
(
3
20
ϵ
-
ϵ
α
6
)
t
-
ϵ
x
3
]
}
v
(
x
,
t
)
=
-
9
5
α
1
cosh
⁡
2
ξ
,
ξ
=
α
10
(
x
+
α
t
)
,
where α is a free positive parameter. The SKdV equations has the following conserved quantities [1, 2, 18]:
(5)
I
1
=
∫
Ω
|
u
|
2
d
x
I
2
=
∫
Ω
u
3
d
x
I
3
=
∫
Ω
[
3
|
u
x
|
2
+
u
3
|
u
|
2
+
1
3
u
3
3
-
1
2
(
u
3
x
)
2
]
d
x
.To avoid complex computation, we assume(6)
u
(
x
,
t
)
=
u
1
(
x
,
t
)
+
i
u
2
(
x
,
t
)
,
i
2
=
-
1
v
(
x
,
t
)
=
u
3
(
x
,
t
)
,
where u
1
(
x
,
t
)
,
u
2
(
x
,
t
), and u
3
(
x
,
t
) are real functions. This will reduce Schrödinger-KdV equations to the coupled system
(7)
ϵ
(
u
1
)
t
+
3
2
(
u
2
)
x
x
-
1
2
u
2
u
3
=
0
ϵ
(
u
2
)
t
-
3
2
(
u
1
)
x
x
+
1
2
u
1
u
3
=
0
(
u
3
)
t
+
1
2
(
u
3
)
x
x
x
+
1
2
(
u
1
2
+
u
2
2
+
u
3
2
)
x
=
0
.


## 2. Numerical Method

Application of the numerical method requires truncation of the infinite interval to a finite interval[
a
,
b
]. We assume that the solution of the coupled Schrödinger-KdV equation is negligible outside this interval. Also we assume all space derivatives at the boundaries approaches to zero in the region (
x
,
t
)
∈
[
a
,
b
]
×
(
0
,
T
].A standard weak formulation [3–5] of this problem is obtained by multiplying (7) by a twice differentiable test function ψ
(
x
) and integrating by parts to obtain
(8)
ϵ
(
u
1
,
t
,
ψ
)
+
3
2
(
u
2
,
ψ
x
x
)
-
1
2
(
u
2
u
3
,
ψ
)
=
0
ϵ
(
u
2
,
t
,
ψ
)
-
3
2
(
u
1
,
ψ
x
x
)
+
1
2
(
u
1
u
3
,
ψ
)
=
0
(
u
3
,
t
,
ψ
)
+
1
2
(
u
3
,
x
,
ψ
x
x
)
+
1
2
(
(
u
1
2
+
u
2
2
+
u
3
2
)
x
,
ψ
)
=
0
,
where (
,
) denotes the usual L
2 inner product
(9)
(
f
,
g
)
=
∫
a
b
‍
f
(
x
)
g
(
x
)
d
x
.The space interval[
a
,
b
] is discretized by uniform (
N
+
1
) grid points
(10)
x
m
=
a
+
m
h
,
m
=
0,1
,
…
,
N
,
where the grid spacing h is given by h
=
(
b
-
a
)
L
/
N. We introduce finite elements in space in (7) and approximate the exact solution of the SKdV equation by
(11)
U
(
x
,
t
)
=
∑
m
=
0
N
‍
U
m
(
t
)
ϕ
m
(
x
)
,
and the product approximation technique is used for treating the nonlinear terms in the following manner:
(12)
U
1
(
x
,
t
)
U
3
(
x
,
t
)
=
∑
m
=
0
N
‍
U
1
,
m
(
t
)
U
3
,
m
(
t
)
ϕ
m
(
x
)
,
U
2
(
x
,
t
)
U
3
(
x
,
t
)
=
∑
m
=
0
N
‍
U
2
,
m
(
t
)
U
3
,
m
(
t
)
ϕ
m
(
x
)
,
[
U
1
2
(
x
,
t
)
+
U
2
2
(
x
,
t
)
+
U
3
2
(
x
,
t
)
]
=
∑
m
=
0
N
‍
[
U
1
,
m
2
(
t
)
+
U
2
,
m
2
(
t
)
+
U
3
,
m
2
(
t
)
]
ϕ
m
(
x
)
.The trial functionsϕ
m
(
x
)
=
ϕ
(
(
x
-
x
m
)
/
h
)
,
m
=
0,1
,
…
,
N, are chosen to be the piecewise linear functions
(13)
ϕ
(
x
)
=
{
1
+
x
if
-
1
<
x
≤
0
1
-
x
if
0
<
x
≤
1
0
otherwise
.The unknown functions[
U
1
,
m
(
t
)
,
U
2
,
m
(
t
)
,
U
3
,
m
(
t
)
]
,
m
=
0,1
,
2
,
…
,
N are determined from the numerical solution of the first-order ordinary differential system
(14)
ϵ
(
U
˙
1
,
ψ
j
)
+
3
2
(
U
2
,
x
x
,
ψ
j
)
-
1
2
(
U
2
U
3
,
ψ
j
)
=
0
ϵ
(
U
˙
2
,
ψ
j
)
-
3
2
(
U
1
,
x
x
,
ψ
j
)
+
1
2
(
U
1
U
3
,
ψ
j
)
=
0
(
U
˙
3
,
ψ
j
)
+
1
2
(
U
3
,
x
,
(
ψ
j
)
x
x
)
+
1
2
(
(
U
1
2
+
U
2
2
+
U
3
2
)
x
,
ψ
j
)
=
0
,
where U
˙
l
=
(
d
/
d
t
)
(
U
l
)
,
l
=
1,2
,
3.We choose the test functionsψ
j
(
x
)
=
ψ
(
(
x
-
x
i
)
/
h
),  j
=
1,2
,
…
,
N, to be the cubic B-spline with compact support
(15)
ψ
(
x
)
=
1
4
{
(
x
+
2
)
3
if
-
2
≤
x
≤
-
1
[
(
2
+
x
)
3
-
4
(
1
+
x
)
3
]
if
-
1
<
x
≤
0
[
(
2
-
x
)
3
-
4
(
x
-
1
)
3
]
if
0
<
x
≤
1
(
2
-
x
)
3
if
1
<
x
≤
2
0
otherwise
.Direct calculations of the inner products in (14), using (11)–(13) and (15), will produce the first-order ordinary differential system
(16)
ϵ
80
[
U
˙
1
,
m
-
2
+
26
U
˙
1
,
m
-
1
+
66
U
˙
1
,
m
+
26
U
˙
1
,
m
+
1
+
U
˙
1
,
m
+
2
]
+
3
8
h
2
[
U
2
,
m
+
2
+
2
U
2
,
m
+
1
-
6
U
2
,
m
+
2
U
2
,
m
-
1
+
U
2
,
m
-
2
]
-
1
160
[
(
U
2
U
3
)
m
-
2
+
26
(
U
2
U
3
)
m
-
1
+
66
(
U
2
U
3
)
m
+
26
(
U
2
U
3
)
m
+
1
+
(
U
2
U
3
)
m
+
2
]
=
0
,
ϵ
80
[
U
˙
2
,
m
-
2
+
26
U
˙
2
,
m
-
1
+
66
U
˙
2
,
m
+
26
U
˙
2
,
m
+
1
+
U
˙
2
,
m
+
2
]
-
3
8
h
2
[
U
2
,
m
+
2
+
2
U
2
,
m
+
1
-
6
U
2
,
m
+
2
U
2
,
m
-
1
+
U
2
,
m
-
2
]
+
1
160
[
(
U
1
U
3
)
m
-
2
+
26
(
U
1
U
3
)
m
-
1
+
66
(
U
1
U
3
)
m
+
26
(
U
1
U
3
)
m
+
1
+
(
U
1
U
3
)
m
+
2
]
=
0
,
1
80
[
U
˙
3
,
m
-
2
+
26
U
˙
3
,
m
-
1
+
66
U
˙
3
,
m
+
26
U
˙
3
,
m
+
1
+
U
˙
3
,
m
+
2
]
+
3
8
h
3
[
U
3
,
m
+
2
-
2
U
3
,
m
+
1
+
2
U
3
,
m
-
1
-
U
3
,
m
-
2
]
+
1
32
h
[
G
(
U
1
,
m
+
2
,
U
2
,
m
+
2
,
U
3
,
m
+
2
)
+
10
G
(
U
1
,
m
+
1
,
U
2
,
m
+
1
,
U
3
,
m
+
1
)
-
10
G
(
U
1
,
m
-
1
,
U
2
,
m
-
1
,
U
3
,
m
-
1
)
-
G
(
U
1
,
m
-
2
,
U
2
,
m
-
2
,
U
3
,
m
-
2
)
]
=
0
,
where
(17)
G
(
u
1
,
u
2
,
u
3
)
=
u
1
2
+
u
2
2
+
u
2
3
.The first order ordinary differential system (16) can be solved by many ordinary differential equation solver. In this work, we will choose two methods, the implicit midpoint rule and the explicit Runge-Kutta method of fourth-order (RK4). These methods will be discussed in the following subsections.


### 2.1. The Implicit Midpoint Rule

By making use of the following substitution for the implicit midpoint rule(18)
U
˙
=
(
U
j
,
m
n
+
1
-
U
j
,
m
n
)
k
,
U
j
,
m
=
U
j
,
m
*
=
(
U
j
,
m
n
+
1
+
U
j
,
m
n
)
2
,
h
h
h
h
h
j
=
1,2
,
3
,
where k is the time step size, into the system (16), this will lead us to the nonlinear block pentadiagonal system
(19)
ϵ
80
k
{
[
U
1
,
m
-
2
n
+
1
+
26
U
1
,
m
-
1
n
+
1
+
66
U
1
,
m
n
+
1
+
26
U
1
,
m
+
1
n
+
1
+
U
1
,
m
+
2
n
+
1
]
-
[
U
1
,
m
-
2
n
+
26
U
1
,
m
-
1
n
+
66
U
1
,
m
n
+
26
U
1
,
m
+
1
n
+
U
1
,
m
+
2
n
]
}
+
3
8
h
2
[
U
2
,
m
+
2
*
+
2
U
2
,
m
+
1
*
-
6
U
2
,
m
*
+
2
U
2
,
m
-
1
*
+
U
2
,
m
-
2
*
]
-
1
160
[
(
U
2
*
U
3
*
)
m
-
2
+
26
(
U
2
*
U
3
*
)
m
-
1
+
66
(
U
2
*
U
3
*
)
m
+
26
(
U
2
*
U
3
*
)
m
+
1
+
(
U
2
*
U
3
*
)
m
+
2
]
=
0
ϵ
80
k
{
[
U
2
,
m
-
2
n
+
1
+
26
U
2
,
m
-
1
n
+
1
+
66
U
2
,
m
n
+
1
+
26
U
2
,
m
+
1
n
+
1
+
U
2
,
m
+
2
n
+
1
]
-
[
U
2
,
m
-
2
n
+
26
U
2
,
m
-
1
n
+
66
U
2
,
m
n
+
26
U
2
,
m
+
1
n
+
U
2
,
m
+
2
n
]
}
-
3
8
h
2
[
U
1
,
m
+
2
*
+
2
U
1
,
m
+
1
*
-
6
U
1
,
m
*
+
2
U
1
,
m
-
1
*
+
U
1
,
m
-
2
*
]
+
1
160
[
(
U
1
*
U
3
*
)
m
-
2
+
26
(
U
1
*
U
3
*
)
m
-
1
+
66
(
U
1
*
U
3
*
)
m
+
26
(
U
1
*
U
3
*
)
m
+
1
+
(
U
1
*
U
3
*
)
m
+
2
]
=
0
1
80
k
{
[
U
3
,
m
-
2
n
+
1
+
26
U
3
,
m
-
1
n
+
1
+
66
U
3
,
m
n
+
1
+
26
U
3
,
m
+
1
n
+
1
+
U
3
,
m
+
2
n
+
1
]
-
[
U
3
,
m
-
2
n
+
26
U
3
,
m
-
1
n
+
66
U
3
,
m
n
+
26
U
3
,
m
+
1
n
+
U
3
,
m
+
2
n
]
}
+
3
8
h
3
[
U
3
,
m
+
2
*
-
2
U
3
,
m
+
1
*
+
2
U
3
,
m
-
1
*
-
U
3
,
m
-
2
*
]
+
1
32
h
[
G
(
U
1
,
m
+
2
*
,
U
2
,
m
+
2
*
,
U
3
,
m
+
2
*
)
+
10
G
(
U
1
,
m
+
1
*
,
U
2
,
m
+
1
*
,
U
3
,
m
+
1
*
)
-
10
G
(
U
1
,
m
-
1
*
,
U
2
,
m
-
1
*
,
U
3
,
m
-
1
*
)
-
G
(
U
1
,
m
-
2
*
,
U
2
,
m
-
2
*
,
U
3
,
m
-
2
*
)
]
=
0
.The numerical solution of the nonlinear system (19) can be obtained by many iterative methods like Newton’s method, fixed point method, and in this work we will adopt the fixed point method and this will be given next.


### 2.2. Fixed Point Method

The fixed point method for solving the nonlinear system (19) can be given as follows:
(20)
ϵ
80
k
{
[
U
1
,
m
-
2
n
+
1
,
s
+
1
+
26
U
1
,
m
-
1
n
+
1
,
s
+
1
+
66
U
1
,
m
n
+
1
,
s
+
1
+
26
U
1
,
m
+
1
n
+
1
,
s
+
1
+
U
1
,
m
+
2
n
+
1
,
s
+
1
]
}
+
3
8
h
2
[
U
2
,
m
+
2
n
+
1
,
s
+
1
+
2
U
2
,
m
+
1
n
+
1
,
s
+
1
-
6
U
2
,
m
n
+
1
,
s
+
1
+
2
U
2
,
m
-
1
n
+
1
,
s
+
1
+
U
2
,
m
-
2
n
+
1
,
s
+
1
]
-
1
160
[
(
U
2
n
+
1
,
s
+
1
U
~
3
)
m
-
2
+
26
(
U
2
n
+
1
,
s
+
1
U
~
3
)
m
-
1
+
66
(
U
2
n
+
1
,
s
+
1
U
~
3
)
m
+
26
(
U
2
n
+
1
,
s
+
1
U
~
3
)
m
+
1
+
(
U
2
n
+
1
,
s
+
1
U
~
3
)
m
+
2
]
=
ϵ
80
k
{
[
U
1
,
m
-
2
n
+
26
U
1
,
m
-
1
n
+
66
U
1
,
m
n
+
26
U
1
,
m
+
1
n
+
U
1
,
m
+
2
n
]
}
+
3
8
h
2
[
U
2
,
m
+
2
n
+
2
U
2
,
m
+
1
n
-
6
U
2
,
m
n
+
2
U
2
,
m
-
1
n
+
U
2
,
m
-
2
n
]
ϵ
80
k
{
[
U
2
,
m
-
2
n
+
1
,
s
+
1
+
26
U
2
,
m
-
1
n
+
1
,
s
+
1
+
66
U
2
,
m
n
+
1
,
s
+
1
+
26
U
2
,
m
+
1
n
+
1
,
s
+
1
+
U
2
,
m
+
2
n
+
1
,
s
+
1
]
}
+
3
8
h
2
[
U
1
,
m
+
2
n
+
1
,
s
+
1
+
2
U
1
,
m
+
1
n
+
1
,
s
+
1
-
6
U
1
,
m
n
+
1
,
s
+
1
+
2
U
1
,
m
-
1
n
+
1
,
s
+
1
+
U
1
,
m
-
2
n
+
1
,
s
+
1
]
-
1
160
[
(
U
1
n
+
1
,
s
+
1
U
~
3
)
m
-
2
+
26
(
U
1
n
+
1
,
s
+
1
U
~
3
)
m
-
1
+
66
(
U
1
n
+
1
,
s
+
1
U
~
3
)
m
+
26
(
U
1
n
+
1
,
s
+
1
U
~
3
)
m
+
1
+
(
U
1
n
+
1
,
s
+
1
U
~
3
)
m
+
2
]
=
ϵ
80
k
{
[
U
2
,
m
-
2
n
+
26
U
2
,
m
-
1
n
+
66
U
2
,
m
n
+
26
U
2
,
m
+
1
n
+
U
2
,
m
+
2
n
]
}
-
3
8
h
2
[
U
1
,
m
+
2
n
+
2
U
1
,
m
+
1
n
-
6
U
1
,
m
n
+
2
U
1
,
m
-
1
n
+
U
1
,
m
-
2
n
]
,
1
80
k
{
[
U
3
,
m
-
2
n
+
1
,
s
+
1
+
26
U
3
,
m
-
1
n
+
1
,
s
+
1
+
66
U
3
,
m
n
+
1
,
s
+
1
+
26
U
3
,
m
+
1
n
+
1
,
s
+
1
+
U
3
,
m
+
2
n
+
1
,
s
+
1
]
}
+
3
8
h
3
[
U
3
,
m
+
2
n
+
1
,
s
+
1
-
2
U
3
,
m
+
1
n
+
1
,
s
+
1
+
2
U
3
,
m
-
1
n
+
1
,
s
+
1
-
U
3
,
m
-
2
n
+
1
,
s
+
1
]
=
1
80
k
{
[
U
3
,
m
-
2
n
+
26
U
3
,
m
-
1
n
+
66
U
3
,
m
n
+
26
U
3
,
m
+
1
n
+
U
3
,
m
+
2
n
]
}
+
3
8
h
3
[
U
3
,
m
+
2
n
-
2
U
3
,
m
+
1
n
+
2
U
3
,
m
-
1
n
-
U
3
,
m
-
2
n
]
+
1
32
h
[
G
(
U
~
1
,
m
+
2
,
U
~
2
,
m
+
2
,
U
~
3
,
m
+
2
)
+
10
G
(
U
~
1
,
m
+
1
,
U
~
2
,
m
+
1
,
U
~
3
,
m
+
1
)
-
10
G
(
U
~
1
,
m
-
1
,
U
~
2
,
m
-
1
,
U
~
3
,
m
-
1
)
-
G
(
U
~
1
,
m
-
2
,
U
~
2
,
m
-
2
,
U
~
3
,
m
-
2
)
]
=
0
,
where
(21)
U
~
j
,
m
=
(
U
j
,
m
n
+
1
,
s
+
U
j
,
m
n
)
2
,
j
=
1,2
,
3
.
U
j
,
m
n
+
1,0
=
U
j
,
m
n
,
j
=
1,2
,
3
.We apply the iterative scheme (20) until the condition
(22)
∥
U
n
+
1
,
s
+
1
-
U
n
+
1
,
s
∥
∞
≤
1
0
-
8
,
s
=
0,1
,
…
is satisfied. Concerning the accuracy of the scheme, the scheme is of second-order in time and fourth-order in space.


### 2.3. Stability of the Scheme

To study the stability of the resulting scheme (19), we use von Neumann stability method. This method can only be applied for linear schemes, so we consider the linear version of the proposed scheme (19) by freezing all terms which make the scheme nonlinear. The linear version of the proposed method can be displayed as follows:
(23)
[
U
1
,
m
-
2
n
+
1
+
26
U
1
,
m
-
1
n
+
1
+
66
U
1
,
m
n
+
1
+
26
U
1
,
m
+
1
n
+
1
+
U
1
,
m
+
2
n
+
1
]
-
[
U
1
,
m
-
2
n
+
26
U
1
,
m
-
1
n
+
66
U
1
,
m
n
+
26
U
1
,
m
+
1
n
+
U
1
,
m
+
2
n
]
+
p
1
[
U
2
,
m
+
2
*
+
2
U
2
,
m
+
1
*
-
6
U
2
,
m
*
+
2
U
2
,
m
-
1
*
+
U
2
,
m
-
2
*
]
-
p
2
[
U
2
,
m
-
2
*
+
26
U
2
,
m
-
1
*
+
66
U
2
,
m
*
+
26
U
2
,
m
+
1
*
+
U
2
,
m
+
2
*
]
=
0
[
U
2
,
m
-
2
n
+
1
+
26
U
2
,
m
-
1
n
+
1
+
66
U
2
,
m
n
+
1
+
26
U
2
,
m
+
1
n
+
1
+
U
2
,
m
+
2
n
+
1
]
-
[
U
2
,
m
-
2
n
+
26
U
2
,
m
-
1
n
+
66
U
2
,
m
n
+
26
U
2
,
m
+
1
n
+
U
2
,
m
+
2
n
]
-
p
1
[
U
1
,
m
+
2
*
+
2
U
1
,
m
+
1
*
-
6
U
1
,
m
*
+
2
U
1
,
m
-
1
*
+
U
1
,
m
-
2
*
]
+
p
2
[
U
1
,
m
-
2
*
+
26
U
1
,
m
-
1
*
+
66
U
1
,
m
*
+
26
U
1
,
m
+
1
*
+
U
1
,
m
+
2
*
]
=
0
1
80
k
{
[
U
3
,
m
-
2
n
+
1
+
26
U
3
,
m
-
1
n
+
1
+
66
U
3
,
m
n
+
1
+
26
U
3
,
m
+
1
n
+
1
+
U
3
,
m
+
2
n
+
1
]
-
[
U
3
,
m
-
2
n
+
26
U
3
,
m
-
1
n
+
66
U
3
,
m
n
+
26
U
3
,
m
+
1
n
+
U
3
,
m
+
2
n
]
}
+
3
8
h
3
[
U
3
,
m
+
2
*
-
2
U
3
,
m
+
1
*
+
2
U
3
,
m
-
1
*
-
U
3
,
m
-
2
*
]
=
0
,
where p
1
=
15
k
/
ϵ
h
2,  p
2
=
k
U
-
/
4
ϵ, and U
- is assumed to be constant on the whole range.We consider the first two equations and the third equation can be done in the similar way. We assume that the solution of (23) to be of the form
(24)
U
1
,
m
n
=
W
1
n
e
i
β
m
h
,
U
2
,
m
n
=
W
2
n
e
i
β
m
h
.By substituting this into (23), this will lead us after some manipulations to the following system:
(25)
γ
1
W
1
n
+
1
+
(
p
1
γ
2
-
p
2
γ
1
)
W
2
n
+
1
=
γ
1
W
1
n
-
(
p
1
γ
2
-
p
2
γ
1
)
W
2
n
-
(
p
1
γ
2
-
p
2
γ
1
)
W
1
n
+
1
+
γ
1
W
2
n
+
1
=
(
p
1
γ
2
-
p
2
γ
1
)
W
1
n
+
γ
1
W
1
n
,
where
(26)
γ
1
=
66
+
2
cos
⁡
(
2
β
h
)
+
52
cos
⁡
(
β
h
)
,
γ
2
=
2
cos
⁡
(
2
β
h
)
+
4
cos
⁡
(
β
h
)
-
6
.
The system (25) can be written in a matrix-vector form as
(27)
Ψ
n
+
1
=
B
Ψ
n
,
where Ψ
n
=
[
W
1
n
,
W
2
n
]
t and B is the (
2
×
2
) matrix
(28)
B
=
[
c
d
-
d
c
]
-
1
[
c
-
d
d
c
]
,
where
(29)
c
=
γ
1
,
d
=
p
1
γ
2
-
p
2
γ
1
.The von Neumann stability condition for the system (27) is the maximum modulus of the eigenvalues of the matrix B are to be less than or equal to one. The eigenvalues of the matrix B are
(30)
λ
1
=
c
+
i
d
c
+
i
d
,
λ
2
=
c
-
i
d
c
+
i
d
,
with modulus equal to one. This means the scheme which we have derived is unconditionally stable according to von Neumann stability analysis.The same methodology of stability can be applied for the third equation.


### 2.4. Runge-Kutta Method of Fourth Order (RK4)

Another approach for solving the semidiscrete system (16) is to use the Runge-Kutta method of fourth-order. The semidiscrete system (16) can be written in a matrix vector form as
(31)
M
d
U
1
d
t
=
F
1
(
U
1
,
U
2
,
U
3
)
,
M
d
U
2
d
t
=
F
2
(
U
1
,
U
2
,
U
3
)
,
M
d
U
3
d
t
=
F
3
(
U
1
,
U
2
,
U
3
)
,
where M is a pentadiagonal matrix of the form
(32)
M
=
1
80
(
66
26
1
0
⋯
⋯
⋯
0
26
66
26
1
0
⋯
⋯
0
1
26
66
26
1
0
⋯
0
0
⋱
⋱
⋱
⋱
⋱
⋱
⋮
⋮
⋱
⋱
⋱
⋱
⋱
⋱
0
⋮
⋱
⋱
⋱
⋱
⋱
⋱
1
⋮
⋱
⋱
⋱
1
26
66
26
0
⋯
⋯
⋯
0
1
26
66
)
,
which can be in a compact form as
(33)
W
˙
=
M
-
1
G
,
where
(34)
W
=
[
U
1
,
U
2
,
U
3
]
t
,
G
=
[
F
1
(
W
)
,
F
2
(
W
)
,
F
3
(
W
)
]
t
.
The RK4 method for the block system (33) can be given as follows:
(35)
W
n
+
1
=
W
n
+
1
6
(
K
1
+
2
K
2
+
2
K
3
+
K
4
)
,
where
(36)
K
1
=
k
M
-
1
G
(
W
n
)
K
2
=
k
M
-
1
G
(
W
n
+
1
2
K
1
)
K
3
=
k
M
-
1
G
(
W
n
+
1
2
K
2
)
K
4
=
k
M
-
1
G
(
W
n
+
K
3
)
.
So in order to calculate K
1
,
K
2
,
K
3, and K
4, we need to solve a linear pentadiagonal system, and this can be done by Crout’s method as follows.At the beginning of the calculations, we factor the matrixM into L
U, where L and U are lower and upper triangular matrices one time only, and then we left only with a lower and an upper triangular systems which can be solved easily by forward and backward substitution process, respectively. This can be easily accomplished with minimum cost. The RK4 method is conditionally stable. In order to get stable results we choose the time step k as (
k
≈
h
3
). Regarding the accuracy in this case the scheme is of fourth-order in time and space directions.


## 2.1. The Implicit Midpoint Rule

By making use of the following substitution for the implicit midpoint rule(18)
U
˙
=
(
U
j
,
m
n
+
1
-
U
j
,
m
n
)
k
,
U
j
,
m
=
U
j
,
m
*
=
(
U
j
,
m
n
+
1
+
U
j
,
m
n
)
2
,
h
h
h
h
h
j
=
1,2
,
3
,
where k is the time step size, into the system (16), this will lead us to the nonlinear block pentadiagonal system
(19)
ϵ
80
k
{
[
U
1
,
m
-
2
n
+
1
+
26
U
1
,
m
-
1
n
+
1
+
66
U
1
,
m
n
+
1
+
26
U
1
,
m
+
1
n
+
1
+
U
1
,
m
+
2
n
+
1
]
-
[
U
1
,
m
-
2
n
+
26
U
1
,
m
-
1
n
+
66
U
1
,
m
n
+
26
U
1
,
m
+
1
n
+
U
1
,
m
+
2
n
]
}
+
3
8
h
2
[
U
2
,
m
+
2
*
+
2
U
2
,
m
+
1
*
-
6
U
2
,
m
*
+
2
U
2
,
m
-
1
*
+
U
2
,
m
-
2
*
]
-
1
160
[
(
U
2
*
U
3
*
)
m
-
2
+
26
(
U
2
*
U
3
*
)
m
-
1
+
66
(
U
2
*
U
3
*
)
m
+
26
(
U
2
*
U
3
*
)
m
+
1
+
(
U
2
*
U
3
*
)
m
+
2
]
=
0
ϵ
80
k
{
[
U
2
,
m
-
2
n
+
1
+
26
U
2
,
m
-
1
n
+
1
+
66
U
2
,
m
n
+
1
+
26
U
2
,
m
+
1
n
+
1
+
U
2
,
m
+
2
n
+
1
]
-
[
U
2
,
m
-
2
n
+
26
U
2
,
m
-
1
n
+
66
U
2
,
m
n
+
26
U
2
,
m
+
1
n
+
U
2
,
m
+
2
n
]
}
-
3
8
h
2
[
U
1
,
m
+
2
*
+
2
U
1
,
m
+
1
*
-
6
U
1
,
m
*
+
2
U
1
,
m
-
1
*
+
U
1
,
m
-
2
*
]
+
1
160
[
(
U
1
*
U
3
*
)
m
-
2
+
26
(
U
1
*
U
3
*
)
m
-
1
+
66
(
U
1
*
U
3
*
)
m
+
26
(
U
1
*
U
3
*
)
m
+
1
+
(
U
1
*
U
3
*
)
m
+
2
]
=
0
1
80
k
{
[
U
3
,
m
-
2
n
+
1
+
26
U
3
,
m
-
1
n
+
1
+
66
U
3
,
m
n
+
1
+
26
U
3
,
m
+
1
n
+
1
+
U
3
,
m
+
2
n
+
1
]
-
[
U
3
,
m
-
2
n
+
26
U
3
,
m
-
1
n
+
66
U
3
,
m
n
+
26
U
3
,
m
+
1
n
+
U
3
,
m
+
2
n
]
}
+
3
8
h
3
[
U
3
,
m
+
2
*
-
2
U
3
,
m
+
1
*
+
2
U
3
,
m
-
1
*
-
U
3
,
m
-
2
*
]
+
1
32
h
[
G
(
U
1
,
m
+
2
*
,
U
2
,
m
+
2
*
,
U
3
,
m
+
2
*
)
+
10
G
(
U
1
,
m
+
1
*
,
U
2
,
m
+
1
*
,
U
3
,
m
+
1
*
)
-
10
G
(
U
1
,
m
-
1
*
,
U
2
,
m
-
1
*
,
U
3
,
m
-
1
*
)
-
G
(
U
1
,
m
-
2
*
,
U
2
,
m
-
2
*
,
U
3
,
m
-
2
*
)
]
=
0
.The numerical solution of the nonlinear system (19) can be obtained by many iterative methods like Newton’s method, fixed point method, and in this work we will adopt the fixed point method and this will be given next.


## 2.2. Fixed Point Method

The fixed point method for solving the nonlinear system (19) can be given as follows:
(20)
ϵ
80
k
{
[
U
1
,
m
-
2
n
+
1
,
s
+
1
+
26
U
1
,
m
-
1
n
+
1
,
s
+
1
+
66
U
1
,
m
n
+
1
,
s
+
1
+
26
U
1
,
m
+
1
n
+
1
,
s
+
1
+
U
1
,
m
+
2
n
+
1
,
s
+
1
]
}
+
3
8
h
2
[
U
2
,
m
+
2
n
+
1
,
s
+
1
+
2
U
2
,
m
+
1
n
+
1
,
s
+
1
-
6
U
2
,
m
n
+
1
,
s
+
1
+
2
U
2
,
m
-
1
n
+
1
,
s
+
1
+
U
2
,
m
-
2
n
+
1
,
s
+
1
]
-
1
160
[
(
U
2
n
+
1
,
s
+
1
U
~
3
)
m
-
2
+
26
(
U
2
n
+
1
,
s
+
1
U
~
3
)
m
-
1
+
66
(
U
2
n
+
1
,
s
+
1
U
~
3
)
m
+
26
(
U
2
n
+
1
,
s
+
1
U
~
3
)
m
+
1
+
(
U
2
n
+
1
,
s
+
1
U
~
3
)
m
+
2
]
=
ϵ
80
k
{
[
U
1
,
m
-
2
n
+
26
U
1
,
m
-
1
n
+
66
U
1
,
m
n
+
26
U
1
,
m
+
1
n
+
U
1
,
m
+
2
n
]
}
+
3
8
h
2
[
U
2
,
m
+
2
n
+
2
U
2
,
m
+
1
n
-
6
U
2
,
m
n
+
2
U
2
,
m
-
1
n
+
U
2
,
m
-
2
n
]
ϵ
80
k
{
[
U
2
,
m
-
2
n
+
1
,
s
+
1
+
26
U
2
,
m
-
1
n
+
1
,
s
+
1
+
66
U
2
,
m
n
+
1
,
s
+
1
+
26
U
2
,
m
+
1
n
+
1
,
s
+
1
+
U
2
,
m
+
2
n
+
1
,
s
+
1
]
}
+
3
8
h
2
[
U
1
,
m
+
2
n
+
1
,
s
+
1
+
2
U
1
,
m
+
1
n
+
1
,
s
+
1
-
6
U
1
,
m
n
+
1
,
s
+
1
+
2
U
1
,
m
-
1
n
+
1
,
s
+
1
+
U
1
,
m
-
2
n
+
1
,
s
+
1
]
-
1
160
[
(
U
1
n
+
1
,
s
+
1
U
~
3
)
m
-
2
+
26
(
U
1
n
+
1
,
s
+
1
U
~
3
)
m
-
1
+
66
(
U
1
n
+
1
,
s
+
1
U
~
3
)
m
+
26
(
U
1
n
+
1
,
s
+
1
U
~
3
)
m
+
1
+
(
U
1
n
+
1
,
s
+
1
U
~
3
)
m
+
2
]
=
ϵ
80
k
{
[
U
2
,
m
-
2
n
+
26
U
2
,
m
-
1
n
+
66
U
2
,
m
n
+
26
U
2
,
m
+
1
n
+
U
2
,
m
+
2
n
]
}
-
3
8
h
2
[
U
1
,
m
+
2
n
+
2
U
1
,
m
+
1
n
-
6
U
1
,
m
n
+
2
U
1
,
m
-
1
n
+
U
1
,
m
-
2
n
]
,
1
80
k
{
[
U
3
,
m
-
2
n
+
1
,
s
+
1
+
26
U
3
,
m
-
1
n
+
1
,
s
+
1
+
66
U
3
,
m
n
+
1
,
s
+
1
+
26
U
3
,
m
+
1
n
+
1
,
s
+
1
+
U
3
,
m
+
2
n
+
1
,
s
+
1
]
}
+
3
8
h
3
[
U
3
,
m
+
2
n
+
1
,
s
+
1
-
2
U
3
,
m
+
1
n
+
1
,
s
+
1
+
2
U
3
,
m
-
1
n
+
1
,
s
+
1
-
U
3
,
m
-
2
n
+
1
,
s
+
1
]
=
1
80
k
{
[
U
3
,
m
-
2
n
+
26
U
3
,
m
-
1
n
+
66
U
3
,
m
n
+
26
U
3
,
m
+
1
n
+
U
3
,
m
+
2
n
]
}
+
3
8
h
3
[
U
3
,
m
+
2
n
-
2
U
3
,
m
+
1
n
+
2
U
3
,
m
-
1
n
-
U
3
,
m
-
2
n
]
+
1
32
h
[
G
(
U
~
1
,
m
+
2
,
U
~
2
,
m
+
2
,
U
~
3
,
m
+
2
)
+
10
G
(
U
~
1
,
m
+
1
,
U
~
2
,
m
+
1
,
U
~
3
,
m
+
1
)
-
10
G
(
U
~
1
,
m
-
1
,
U
~
2
,
m
-
1
,
U
~
3
,
m
-
1
)
-
G
(
U
~
1
,
m
-
2
,
U
~
2
,
m
-
2
,
U
~
3
,
m
-
2
)
]
=
0
,
where
(21)
U
~
j
,
m
=
(
U
j
,
m
n
+
1
,
s
+
U
j
,
m
n
)
2
,
j
=
1,2
,
3
.
U
j
,
m
n
+
1,0
=
U
j
,
m
n
,
j
=
1,2
,
3
.We apply the iterative scheme (20) until the condition
(22)
∥
U
n
+
1
,
s
+
1
-
U
n
+
1
,
s
∥
∞
≤
1
0
-
8
,
s
=
0,1
,
…
is satisfied. Concerning the accuracy of the scheme, the scheme is of second-order in time and fourth-order in space.


## 2.3. Stability of the Scheme

To study the stability of the resulting scheme (19), we use von Neumann stability method. This method can only be applied for linear schemes, so we consider the linear version of the proposed scheme (19) by freezing all terms which make the scheme nonlinear. The linear version of the proposed method can be displayed as follows:
(23)
[
U
1
,
m
-
2
n
+
1
+
26
U
1
,
m
-
1
n
+
1
+
66
U
1
,
m
n
+
1
+
26
U
1
,
m
+
1
n
+
1
+
U
1
,
m
+
2
n
+
1
]
-
[
U
1
,
m
-
2
n
+
26
U
1
,
m
-
1
n
+
66
U
1
,
m
n
+
26
U
1
,
m
+
1
n
+
U
1
,
m
+
2
n
]
+
p
1
[
U
2
,
m
+
2
*
+
2
U
2
,
m
+
1
*
-
6
U
2
,
m
*
+
2
U
2
,
m
-
1
*
+
U
2
,
m
-
2
*
]
-
p
2
[
U
2
,
m
-
2
*
+
26
U
2
,
m
-
1
*
+
66
U
2
,
m
*
+
26
U
2
,
m
+
1
*
+
U
2
,
m
+
2
*
]
=
0
[
U
2
,
m
-
2
n
+
1
+
26
U
2
,
m
-
1
n
+
1
+
66
U
2
,
m
n
+
1
+
26
U
2
,
m
+
1
n
+
1
+
U
2
,
m
+
2
n
+
1
]
-
[
U
2
,
m
-
2
n
+
26
U
2
,
m
-
1
n
+
66
U
2
,
m
n
+
26
U
2
,
m
+
1
n
+
U
2
,
m
+
2
n
]
-
p
1
[
U
1
,
m
+
2
*
+
2
U
1
,
m
+
1
*
-
6
U
1
,
m
*
+
2
U
1
,
m
-
1
*
+
U
1
,
m
-
2
*
]
+
p
2
[
U
1
,
m
-
2
*
+
26
U
1
,
m
-
1
*
+
66
U
1
,
m
*
+
26
U
1
,
m
+
1
*
+
U
1
,
m
+
2
*
]
=
0
1
80
k
{
[
U
3
,
m
-
2
n
+
1
+
26
U
3
,
m
-
1
n
+
1
+
66
U
3
,
m
n
+
1
+
26
U
3
,
m
+
1
n
+
1
+
U
3
,
m
+
2
n
+
1
]
-
[
U
3
,
m
-
2
n
+
26
U
3
,
m
-
1
n
+
66
U
3
,
m
n
+
26
U
3
,
m
+
1
n
+
U
3
,
m
+
2
n
]
}
+
3
8
h
3
[
U
3
,
m
+
2
*
-
2
U
3
,
m
+
1
*
+
2
U
3
,
m
-
1
*
-
U
3
,
m
-
2
*
]
=
0
,
where p
1
=
15
k
/
ϵ
h
2,  p
2
=
k
U
-
/
4
ϵ, and U
- is assumed to be constant on the whole range.We consider the first two equations and the third equation can be done in the similar way. We assume that the solution of (23) to be of the form
(24)
U
1
,
m
n
=
W
1
n
e
i
β
m
h
,
U
2
,
m
n
=
W
2
n
e
i
β
m
h
.By substituting this into (23), this will lead us after some manipulations to the following system:
(25)
γ
1
W
1
n
+
1
+
(
p
1
γ
2
-
p
2
γ
1
)
W
2
n
+
1
=
γ
1
W
1
n
-
(
p
1
γ
2
-
p
2
γ
1
)
W
2
n
-
(
p
1
γ
2
-
p
2
γ
1
)
W
1
n
+
1
+
γ
1
W
2
n
+
1
=
(
p
1
γ
2
-
p
2
γ
1
)
W
1
n
+
γ
1
W
1
n
,
where
(26)
γ
1
=
66
+
2
cos
⁡
(
2
β
h
)
+
52
cos
⁡
(
β
h
)
,
γ
2
=
2
cos
⁡
(
2
β
h
)
+
4
cos
⁡
(
β
h
)
-
6
.
The system (25) can be written in a matrix-vector form as
(27)
Ψ
n
+
1
=
B
Ψ
n
,
where Ψ
n
=
[
W
1
n
,
W
2
n
]
t and B is the (
2
×
2
) matrix
(28)
B
=
[
c
d
-
d
c
]
-
1
[
c
-
d
d
c
]
,
where
(29)
c
=
γ
1
,
d
=
p
1
γ
2
-
p
2
γ
1
.The von Neumann stability condition for the system (27) is the maximum modulus of the eigenvalues of the matrix B are to be less than or equal to one. The eigenvalues of the matrix B are
(30)
λ
1
=
c
+
i
d
c
+
i
d
,
λ
2
=
c
-
i
d
c
+
i
d
,
with modulus equal to one. This means the scheme which we have derived is unconditionally stable according to von Neumann stability analysis.The same methodology of stability can be applied for the third equation.


## 2.4. Runge-Kutta Method of Fourth Order (RK4)

Another approach for solving the semidiscrete system (16) is to use the Runge-Kutta method of fourth-order. The semidiscrete system (16) can be written in a matrix vector form as
(31)
M
d
U
1
d
t
=
F
1
(
U
1
,
U
2
,
U
3
)
,
M
d
U
2
d
t
=
F
2
(
U
1
,
U
2
,
U
3
)
,
M
d
U
3
d
t
=
F
3
(
U
1
,
U
2
,
U
3
)
,
where M is a pentadiagonal matrix of the form
(32)
M
=
1
80
(
66
26
1
0
⋯
⋯
⋯
0
26
66
26
1
0
⋯
⋯
0
1
26
66
26
1
0
⋯
0
0
⋱
⋱
⋱
⋱
⋱
⋱
⋮
⋮
⋱
⋱
⋱
⋱
⋱
⋱
0
⋮
⋱
⋱
⋱
⋱
⋱
⋱
1
⋮
⋱
⋱
⋱
1
26
66
26
0
⋯
⋯
⋯
0
1
26
66
)
,
which can be in a compact form as
(33)
W
˙
=
M
-
1
G
,
where
(34)
W
=
[
U
1
,
U
2
,
U
3
]
t
,
G
=
[
F
1
(
W
)
,
F
2
(
W
)
,
F
3
(
W
)
]
t
.
The RK4 method for the block system (33) can be given as follows:
(35)
W
n
+
1
=
W
n
+
1
6
(
K
1
+
2
K
2
+
2
K
3
+
K
4
)
,
where
(36)
K
1
=
k
M
-
1
G
(
W
n
)
K
2
=
k
M
-
1
G
(
W
n
+
1
2
K
1
)
K
3
=
k
M
-
1
G
(
W
n
+
1
2
K
2
)
K
4
=
k
M
-
1
G
(
W
n
+
K
3
)
.
So in order to calculate K
1
,
K
2
,
K
3, and K
4, we need to solve a linear pentadiagonal system, and this can be done by Crout’s method as follows.At the beginning of the calculations, we factor the matrixM into L
U, where L and U are lower and upper triangular matrices one time only, and then we left only with a lower and an upper triangular systems which can be solved easily by forward and backward substitution process, respectively. This can be easily accomplished with minimum cost. The RK4 method is conditionally stable. In order to get stable results we choose the time step k as (
k
≈
h
3
). Regarding the accuracy in this case the scheme is of fourth-order in time and space directions.


## 3. Numerical Results

In this section, we present some numerical results for the proposed schemes. The accuracy of our methods is tested by calculating theL
∞ error norms
(37)
L
∞
u
=
max
⁡
m
|
U
m
n
-
u
(
x
m
,
t
n
)
|
,
L
∞
v
=
max
⁡
m
|
V
m
n
-
v
(
x
m
,
t
n
)
|
.The conserved quantities of the coupled Schrödinger-KdV equation are calculated using trapezoidal rule. In all calculations, we choose the initial condition(38)
u
(
x
,
0
)
=
-
6
5
3
α
tanh
ξ
cosh
⁡
ξ
exp
⁡
{
i
α
[
-
ϵ
x
3
]
}
v
=
-
9
5
α
1
cosh
⁡
2
ξ
,
where ξ
=
(
α
/
10
)
x; α and ϵ are free positive parameters. The following tests will be discussed(i) Implicit Midpoint Rule with
ϵ
=
1. To compute the numerical solution using the system obtained by using the implicit midpoint rule (19), the following parameters are used:
(39)
x
l
=
-
50
,
x
r
=
50
,
h
=
0.1
,
k
=
0.1
,
α
=
0.45
,
ϵ
=
1.0
.In Tables1 and 2, we display the errors and the conserved quantities. In Figures 1 and 2, we display the numerical solution of |
u
| and v for t
=
0,0.5,1
,
…
,
10. The maximum of |
u
| and v is constant and equal, respectively, to 0.4676
,
-
0.8100 for t
=
0,1
,
…
,
10. It is very easy to see the accuracy and how the conserved quantities are almost conserved.Table 1
Errors using implicit midpoint rule.
Time
L
∞
u
L
∞
v
1
5.3165142
E
-
05
4.7715757
E
-
06
2
5.5607597
E
-
05
9.5039732
E
-
06
5
7.150639
E
-
05
2.1101674
E
-
05
8
9.4115007
E
-
05
2.9205086
E
-
05Table 2
Conservative quantities.
Time
I
1
I
2
I
3
0
2.749231
-
7.636753
-
1.151015
2
2.749229
-
7.636751
-
1.151014
5
2.749226
-
7.636747
-
1.151012
8
2.749223
-
7.636743
-
1.151010Figure 1
The evolution of the numerical solution of|
u
| with ϵ
=
1.Figure 2
The evolution of the numerical solution ofv with ϵ
=
1.0.(ii) Implicit Midpoint Rule with
ϵ
=
8. We compute the numerical solution using the implicit midpoint rule. We choose in this case ϵ
=
8, we display the numerical solution of |
u
| and v for t
=
0,1
,
…
10. In Figures 3 and 4, we display the real (U
1
n
) and imaginary (U
2
n
) parts of the numerical solution, the same maximum values obtained as the case of ϵ
=
1.0. This is compatible with the exact solution, and we have noticed some differences between our results in Figures 3 and 4 with the one in [2]. We think this occurrence is due to their selected interval [
-
30,30
] where boundary values are not zero. In Figures 5 and 6, we display the simulation of the numerical solution of the modulus of u and U
3
,
m
n, respectively.Figure 3
The evolution of the numerical solution of real part ofu and ϵ
=
8.0.Figure 4
The evolution of the numerical solution of imaginary l part ofu and ϵ
=
8.0.Figure 5
The evolution of the numerical solution of|
u
| with ϵ
=
8.0.Figure 6
The evolution of the numerical solution ofv with ϵ
=
8.0.(iii) Runge-Kutta Method of Fourth-Order. Now we will test the method obtained by using the RK4. We choose the following set of parameters:
(40)
x
l
=
-
50
,
x
r
=
50
,
h
=
0.1
,
k
=
0.0001
,
α
=
0.45
,
ϵ
=
1.0
.The errors and conserved quantities are displayed in Tables3 and 4; we can easily see the high accuracy of the proposed method and how the conserved quantities are exactly conserved.Table 3
Errors using Runge-Kutta method.
Time
L
∞
u
L
∞
v
0.2
4.1779
E
-
05
6.5270
E
-
08
0.6
4.5066
E
-
05
5.1653
E
-
08
0.8
4.6236
E
-
05
4.9147
E
-
08
1.0
4.7321
E
-
05
5.6538
E
-
08Table 4
Conservative quantities.
Time
I
1
I
2
I
3
0
2.749231
-
7.636753
-
1.151015
2
2.749231
-
7.636753
-
1.151015
5
2.749231
-
7.636753
-
1.151015
8
2.749231
-
7.636753
-
1.151015(iv) Comparison of Some Existing Methods. A comparative study has been conducted with some existing methods and the results are displayed in Table 5. From the numerical results displayed for different methods, we can see the superiority of our proposed method comparing to the previous methods which are mentioned in [1, 2].Table 5
The comparison of maximum errors using different schemes att
=
0.1 with k
=
0.0001.
Method
L
∞
u
(
h
=
1
)
L
∞
v
(
h
=
1
)
Present (with midpoint rule)
1.69295
E
-
5
2.468298
E
-
05
Present (with RK4)
8.3382
E
-
06
4.2857
E
-
06
SSQBS FEM
8.831056
E
-
04
3.213056
E
-
04
SSCBS FEM
8.093879
E
-
04
2.723093
E
-
04
Semidiscrete FEM
8.881089
E
-
04
3.221997
E
-
04
Crank Nicolson
8.928318
E
-
04
5.431211
E
-
04


## 4. Conclusion

In this work, we derived highly accurate numerical schemes for the SKdV equation, using finite element method (Petrov-Galerkin) with product approximation technique, where we have used the linear hat functions as trial function and the cubic B-spline as test functions. The differential systems obtained are solved by using implicit midpoint rule. Also we have solved the first-order ordinary differential system using the RK4, and in this case the numerical scheme is conditionally stable and of fourth-order in space and time directions. We held a comparison with some existing methods, and we found that our scheme produced highly accurate results and conserved the conserved quantities almost exactly.



---

*Source: 705204-2014-06-15.xml*