# Nonlinear Stability and Convergence of Two-Step Runge-Kutta Methods for Volterra Delay Integro-Differential Equations

Authors: Haiyan Yuan; Cheng Song
Journal: Abstract and Applied Analysis
 (2013)

---

## Abstract

This paper introduces the stability and convergence of two-step Runge-Kutta methods with compound quadrature formula for solving nonlinear Volterra delay integro-differential equations. First, the definitions of(
k
,
l
)-algebraically stable and asymptotically stable are introduced; then the asymptotical stability of a (
k
,
l
)-algebraically stable two-step Runge-Kutta method with 0
<
k
<
1 is proved. For the convergence, the concepts of D-convergence, diagonally stable, and generalized stage order are firstly introduced; then it is proved by some theorems that if a two-step Runge-Kutta method is algebraically stable and diagonally stable and its generalized stage order is p, then the method with compound quadrature formula is D-convergent of order at least min
{
p
,
ν
}, where ν depends on the compound quadrature formula.


---

## Body

## 1. Introduction

Volterra delay integro-differential equations (VDIDEs) arise widely in the mathematical modeling of physical and biological phenomena. Significant advances in the theoretical analysis and in the numerical analysis for these problems have been made in the last few decades (see, e.g., [1, 2]). For the case of linear stability and convergence for these equations, numerical time-integration techniques of one-step collocation and Runge-Kutta type were investigated in [3–8]. Linear multistep-based methods were studied in [9–12]. De la Sen and Luo studied the uniform exponential stability of a wide class of linear time-delay systems in [13]; De la Sen considered the stability of impulsive time-varying systems in [14].For the case of nonlinear stability and convergence, stability results were obtained in [15, 16], where the authors investigated the nonlinear stability of continuous Runge-Kutta methods, discrete Runge-Kutta methods, and backward differentiation (BDF) methods, respectively. However, most of these important results are based on the classical Lipschitz conditions, while the classical Lipschitz conditions are so strong that there are few equations satisfying them. Most of the Volterra delay integro-differential equations satisfy the one-sided Lipschitz condition, while the studies focusing on the stability and convergence of the numerical method for nonlinear VDIDEs based on a one-sided Lipschitz condition have not yet been seen in the literature. By means of a one-sided Lipschitz condition, we will discuss the stability and convergence of two-step Runge-Kutta (TSRK) methods for nonlinear VDIDEs in the present paper.The paper is organized as follows. In Section2, a fairly general class of VDIDEs is defined. We present a stability criterion for such problems, which generalizes the criteria in the above references. A class of two-step Runge-Kutta methods is also derived for solving VDIDEs. They are obtained by compound quadrature rules. In Sections 3 and 4, nonlinear stability and convergence of TSRK method for NDDEs are derived and proved. In Section 5 we present some numerical examples in order to illustrate the nonlinear stability and convergence of a two-step Runge-Kutta method. These numerical results show that the new methods are quite effective.


## 2. A Class of VDIDES and the Two-Step Runge-Kutta Methods

It is the purpose of this paper to investigate the nonlinear stability and convergence properties of the following initial-value problem VDIDEs:(1)
y
′
(
t
)
=
f
(
t
,
y
(
t
)
,
y
(
t
-
τ
)
,
∫
t
-
τ
t
g
(
t
,
v
,
y
(
v
)
)
d
v
)
,
t
∈
[
0
,
Τ
]
,
y
(
t
)
=
φ
(
t
)
,
t
∈
[
-
τ
,
0
]
,
where f, g, and φ are smooth enough such that (1) has a unique solution y
(
t
) and τ is a positive delay term. We assume that there exist some inner product 〈
·
,
·
〉 and the induced norm |
|
·
|
| in C
N, such that φ
:
[
-
τ
,
0
]
→
C
N and f
:
[
0
,
Τ
]
×
C
N
×
C
N
×
C
N
→
C
N satisfy the following conditions:
(2)
Re
〈
f
(
t
,
u
1
,
v
,
w
)
-
f
(
t
,
u
2
,
v
,
w
)
,
u
1
-
u
2
〉
≤
α
∥
u
1
-
u
2
∥
2
,
(3)
∥
f
(
t
,
u
,
v
1
,
w
1
)
-
f
(
t
,
u
,
v
2
,
w
2
)
∥
≤
β
∥
v
1
-
v
2
∥
+
σ
∥
w
1
-
w
2
∥
,
(4)
∥
g
(
t
,
v
,
u
1
)
-
g
(
t
,
v
,
u
2
)
∥
≤
γ
∥
u
1
-
u
2
∥
,
(
t
,
v
)
∈
D
,
for t
∈
[
0
,
Τ
], D
=
{
(
t
,
v
)
:
t
∈
[
0
,
Τ
)
,
v
∈
[
t
-
τ
,
t
]
}, and for all u
1
,
u
2
,
v
1
,
v
2
,
w
1
,
w
2
∈
C
N
(
-
α
)
,
β
,
σ, and γ are nonnegative constants.Throughout this paper, we assume that (1) has unique solution y
(
t
) and denote the problem class R
(
α
,
β
,
σ
,
γ
) that consist problem of type (1) with (2)–(4). In order to make the error analysis feasible, we always assume that (1) has a unique solution y
(
t
) which is sufficiently differentiable and satisfies
(5)
∥
d
i
y
(
t
)
d
t
i
∥
≤
M
i
,
i
=
1,2
,
…
,
t
∈
[
-
τ
,
T
]
.
For function g
(
t
,
θ
,
y
(
θ
)
), we make the following assumptions, unless otherwise stated; all its partial derivatives used later exist and satisfy the following:
(6)
∥
∂
i
g
(
t
,
θ
,
y
(
θ
)
)
∂
θ
i
∥
≤
N
i
,
t
∈
[
0
,
T
]
,
θ
∈
[
-
τ
,
T
]
.
Many numerical methods have been proposed for the numerical solution of (1).In this paper, we are concerned with two-step Runge-Kutta (TSRK) methods of the form(7a)
Y
i
(
n
)
=
y
n
+
h
∑
j
=
1
s
a
i
j
f
(
t
n
+
c
j
h
,
Y
j
(
n
)
)
,
i
=
1,2
,
…
,
s
,
(7b)
Y
i
(
n
-
1
)
=
y
n
-
1
+
h
∑
j
=
1
s
a
i
j
f
(
t
n
-
1
+
c
j
h
,
Y
j
(
n
-
1
)
)
,
i
=
1,2
,
…
,
s
,
(7c)
y
n
+
1
=
(
1
-
θ
)
y
n
+
θ
y
n
-
1
+
h
∑
i
=
1
s
b
i
f
(
t
n
+
c
i
h
,
Y
i
(
n
)
)
+
h
∑
i
=
1
s
b
~
i
f
(
t
n
-
1
+
c
i
h
,
Y
i
(
n
-
1
)
)
,where ∑
i
=
1
s
b
i
+
∑
i
=
1
s
b
~
i
=
1
+
θ, c
i
=
∑
j
=
1
s
a
i
j, h
=
τ
/
m is a stepsize, m is an arbitrarily given positive integer, and 0
≤
θ
≤
1. The above methods are studied in [17]. Now we consider the adaptation of the two-step Runge-Kutta method to (1):(8a)
Y
i
(
n
)
=
y
n
+
h
∑
j
=
1
s
a
i
j
f
(
t
n
+
c
j
h
,
Y
j
(
n
)
,
Y
-
j
(
n
)
,
Y
~
j
(
n
)
)
,
i
=
1,2
,
…
,
s
,
(8b)
Y
i
(
n
-
1
)
=
y
n
-
1
+
h
∑
j
=
1
s
a
i
j
f
(
t
n
-
1
+
c
j
h
,
Y
j
(
n
-
1
)
,
Y
-
j
(
n
-
1
)
,
Y
~
j
(
n
-
1
)
)
,
i
=
1,2
,
…
,
s
,
(8c)
y
n
+
1
=
(
1
-
θ
)
y
n
+
θ
y
n
-
1
+
h
∑
i
=
1
s
b
i
f
(
t
n
+
c
i
h
,
Y
i
(
n
)
,
Y
-
i
(
n
)
,
Y
~
i
(
n
)
)
+
h
∑
i
=
1
s
b
~
i
f
(
t
n
-
1
+
c
i
h
,
Y
i
(
n
-
1
)
,
Y
-
i
(
n
-
1
)
,
Y
~
i
(
n
-
1
)
)
.In particular, y
0
=
φ
(
0
), y
n is the numerical approximation at t
n
=
n
h to the analytic solution y
(
t
n
), the argument Y
-
j
(
n
)
=
Y
j
(
n
-
m
) denotes an approximation to y
(
t
n
+
c
j
h
-
τ
), and the argument Y
~
j
(
n
) denotes an approximation to ∫
t
n
+
c
j
h
-
τ
t
n
+
c
j
h
g
(
t
n
+
c
j
h
,
ξ
,
y
(
ξ
)
)
d
ξ which are obtained by a convergent compound quadrature (CQ) formula:
(9)
Y
~
j
(
n
)
=
h
∑
q
=
0
m
p
~
q
g
(
t
n
+
c
j
h
,
t
n
-
q
+
c
j
h
,
Y
j
(
n
-
q
)
)
,
j
=
1,2
,
…
,
s
,
using values Y
j
(
n
-
q
)
=
φ
(
t
n
-
q
+
c
j
h
) with n
≤
q, t
n
-
q
+
c
j
h
≤
0.The class of Runge-Kutta methods with CQ formula has been applied to delay-integro-differential equations by many authors (c.f. [18, 19]). For the CQ formula (9), we usually adopt the repeated trapezoidal rule, the repeated Simpson’s rule, or the repeated Newton-cotes rule, and so forth, denote η
=
max
{
p
~
0
,
p
~
1
,
…
,
p
~
m
}. It should be pointed out that the adopted quadrature formula (9) is only a class of quadrature formula for Y
~
i
(
n
), there also exist some other types of quadrature formula, such as Pouzet quadrature (PQ) formula and the quadrature formula based on Laguerre-Radau interpolations [20, 21]. It is the aim of our future research to investigate the adaptation of PQ formula and the quadrature formula based on Laguerre-Radau interpolations to VDIDEs.


## 3. The Nonlinear Stability Analysis

In this section, we will investigate the stability of the two-step Runge-Kutta methods for VDIDEs. In order to consider the stability property, we also need to consider the perturbed problem of (1):
(10)
z
′
(
t
)
=
f
(
t
,
z
(
t
)
,
z
(
t
-
τ
)
,
∫
t
-
τ
t
g
(
t
,
ξ
,
z
(
ξ
)
)
d
ξ
)
,
t
∈
[
0
,
Τ
]
,
z
(
t
)
=
ψ
(
t
)
,
t
∈
[
-
τ
,
0
]
,
where ψ
:
[
-
τ
,
0
]
→
C
N is a given continuous function. The unique exact solution of the problem (10) is denoted as z
(
t
).Applying the two-step Runge-Kutta method (7a)–(7c) to (10) leads to(11a)
Z
i
(
n
)
=
z
n
+
h
∑
j
=
1
s
a
i
j
f
(
t
n
+
c
j
h
,
Z
j
(
n
)
,
Z
-
j
(
n
)
,
Z
~
j
(
n
)
)
,
i
=
1,2
,
…
,
s
,
(11b)
Z
i
(
n
-
1
)
=
z
n
-
1
+
h
∑
j
=
1
s
a
i
j
f
(
t
n
-
1
+
c
j
h
,
Z
j
(
n
-
1
)
,
Z
-
j
(
n
-
1
)
,
Z
~
j
(
n
-
1
)
)
,
i
=
1,2
,
…
,
s
,
(11c)
z
n
+
1
=
(
1
-
θ
)
z
n
+
θ
z
n
-
1
+
h
∑
i
=
1
s
b
i
f
(
t
n
+
c
i
h
,
Z
i
(
n
)
,
Z
-
i
(
n
)
,
Z
~
i
(
n
)
)
+
h
∑
i
=
1
s
b
~
i
f
(
t
n
-
1
+
c
i
h
,
Z
i
(
n
-
1
)
,
Z
-
i
(
n
-
1
)
,
Z
~
i
(
n
-
1
)
)
,where z
n and Z
i
(
n
) denote approximations to z
(
t
n
) and z
(
t
n
+
c
i
h
), respectively; the argument Z
~
j
(
n
) denotes an approximation to ∫
t
n
+
c
j
h
-
τ
t
n
+
c
j
h
g
(
t
n
+
c
j
h
,
ξ
,
z
(
ξ
)
)
d
ξ which are obtained by a convergent compound quadrature (CQ) formula:
(12)
Z
~
i
(
n
)
=
h
∑
q
=
0
m
p
~
q
g
(
t
n
+
c
i
h
,
t
n
-
q
+
c
i
h
,
Z
i
(
n
-
q
)
)
,
i
=
1,2
,
…
,
s
,
and Z
-
j
(
n
)
=
Z
j
(
n
-
m
) denotes an approximation to z
(
t
n
+
c
j
h
-
τ
), and z
n
=
ψ
(
t
n
) for n
=
0,
(13)
Z
i
(
n
-
q
)
=
ψ
(
t
n
-
q
+
c
i
h
)
,
for
t
n
-
q
+
c
i
h
≤
0
.


### 3.1. Some Concepts

For the stability analysis, we need the compound quadrature formula (9) that satisfies the following condition:
(14)
h
2
(
m
+
1
)
2
∑
q
=
0
m
|
p
~
q
|
2
<
p
~
2
,
where p
~ is a positive constant.Let(15)
w
n
=
y
n
-
z
n
,
W
j
(
n
)
=
Y
j
(
n
)
-
Z
j
(
n
)
,
W
~
j
(
n
)
=
Y
~
j
(
n
)
-
Z
~
j
(
n
)
,
Q
j
(
n
)
=
h
[
f
(
t
n
+
c
j
h
,
Y
j
(
n
)
,
Y
-
j
(
n
)
,
Y
~
j
(
n
)
)
-
f
(
t
n
+
c
j
h
,
Z
j
(
n
)
,
Z
-
j
(
n
)
,
Z
~
j
(
n
)
)
]
,
j
=
1,2
,
…
,
s
.
It follows from (8a)–(8c) and (11a)–(11c) that(16a)
W
i
(
n
)
=
w
n
+
∑
j
=
1
s
a
i
j
Q
j
(
n
)
,
i
=
1,2
,
…
,
s
,
(16b)
W
i
(
n
-
1
)
=
w
n
-
1
+
∑
j
=
1
s
a
i
j
Q
j
(
n
-
1
)
,
i
=
1,2
,
…
,
s
,
(16c)
w
n
+
1
=
(
1
-
θ
)
w
n
+
θ
w
n
-
1
+
∑
i
=
1
s
b
i
Q
i
(
n
)
+
∑
i
=
1
s
b
~
i
Q
i
(
n
-
1
)
.Now we will write the s-stage TSRK methods (7a)–(7c) as a general linear method.LetV
i
(
n
)
=
(
W
i
(
n
)
Τ
,
W
i
(
n
-
1
)
Τ
)
Τ be the internal stages and μ
n
=
(
w
n
+
1
Τ
,
w
n
Τ
)
Τ the external vectors and P
i
(
n
)
=
(
Q
i
(
n
)
Τ
,
Q
i
(
n
-
1
)
Τ
)
Τ. Then we have a 2
(
s
+
1
)-stage partitioned general linear method:(17a)
V
i
(
n
)
=
∑
j
=
1
s
C
i
j
11
P
i
(
n
)
+
∑
j
=
1
s
C
i
j
12
μ
n
,
i
=
1,2
,
…
,
s
,
(17b)
μ
n
+
1
=
∑
j
=
1
s
C
i
j
21
P
i
(
n
)
+
∑
j
=
1
s
C
i
j
22
μ
n
,
i
=
1,2
,
…
,
s
,where C
i
j
11
=
a
i
j
,
i
,
j
=
1,2
,
⋯
,
s, C
i
j
11
=
0
,
i
=
1,2
,
⋯
,
s
,
j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, C
i
j
11
=
0, i
=
s
+
1
,
s
+
2
,
⋯
,
2
s, j
=
1,2
,
⋯
,
s, C
i
j
11
=
a
i
j, i
,
j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, C
i
j
12
=
1, i
,
j
=
1,2
,
⋯
,
s, C
i
j
12
=
0, i
=
1,2
,
⋯
,
s, j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, C
i
j
12
=
0, i
=
s
+
1
,
s
+
2
,
⋯
,
2
s, j
=
1,2
,
⋯
,
s, C
i
j
12
=
1, i
,
j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, C
i
j
21
=
b
j, i
,
j
=
1,2
,
⋯
,
s, C
i
j
21
=
b
~
j, i
=
1,2
,
⋯
,
s, j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, C
i
j
21
=
0, i
=
s
+
1
,
s
+
2
,
⋯
,
2
s, j
=
1,2
,
⋯
,
s, C
i
j
21
=
0, i
,
j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, C
i
j
22
=
1
-
θ, i
,
j
=
1,2
,
⋯
,
s, C
i
j
22
=
θ, i
=
1,2
,
⋯
,
s, j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, C
i
j
22
=
1, i
=
s
+
1
,
s
+
2
,
⋯
,
2
s, j
=
1,2
,
⋯
,
s, C
i
j
22
=
0, i
,
j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, and A
=
(
a
i
j
), e
=
(
1,1
,
…
1
), b
=
(
b
1
,
b
2
,
…
,
b
s
), b
~
=
(
b
~
1
,
b
~
2
,
…
,
b
~
s
).We can identify the coefficient matrices(18)
C
11
=
(
C
i
j
11
)
=
(
A
0
0
A
)
,
C
12
=
(
C
i
j
12
)
=
(
e
0
0
e
)
,
C
21
=
(
C
i
j
21
)
=
(
b
b
~
0
0
)
,
C
22
=
(
C
i
j
22
)
=
(
1
-
θ
θ
1
0
)
.Definition 1 (see [22]).
Letk
,
l be real constants, a TSRK method is said to be (
k
,
l
) algebraically stable if there exists a diagonal matrix D
=
diag
(
d
1
,
d
2
,
…
,
d
2
s
) and a diagonal nonnegative matrix G such that M
=
(
m
i
j
) is nonnegative, where(19)
M
(
k
,
l
)
=
(
k
G
-
C
22
Τ
G
C
22
-
2
l
C
12
Τ
D
C
12
C
12
Τ
D
-
C
22
Τ
G
C
21
-
2
l
C
12
Τ
D
C
11
D
C
12
-
C
21
Τ
G
C
22
-
2
l
C
11
Τ
D
C
12
C
11
Τ
D
+
D
C
11
-
C
21
Τ
G
C
21
-
2
l
C
11
Τ
D
C
11
)
.In particular, a (
1,0
)-algebraically stable method is called algebraically stable.Definition 2.
The TSRK method (7a)–(7c) is called asymptotically stable for (1) with (2)–(4) if lim
n
→
∞
∥
w
n
∥
=
0.


### 3.2. Numerical Stability of the Methods

Numerical stability is an important feature of an effective numerical method. An unstable numerical method may be consistent of high order, yet arbitrarily small perturbations will eventually cause large deviations from the true solution. In this section, we will focus on the asymptotic stability of the TSRK method.Theorem 3.
Assume that the TSRK method (7a)–(7c) is (
k
,
l
)-algebraically stable with 0
<
k
≤
1, suppose that the quadrature formula (9) satisfies the condition (14) and the conditions (2)–(4) hold, then, method (17a) and (17b) satisfies the following:
(20)
∥
y
n
-
z
n
∥
≤
C
max
t
∈
[
-
τ
,
0
]
∥
φ
(
t
)
-
ψ
(
t
)
∥
,
n
=
1,2
,
…
,
when the following condition holds:
(21)
h
(
α
+
β
+
σ
+
σ
γ
2
p
~
2
)
≤
l
,
where C depends only on the method, the parameters α, β, σ, γ, η, and τ.Proof.
It follows from a fairly straightforward (but tedious) computation and(
k
,
l
) algebraically stability that (compare also [22, 23])
(22)
∥
μ
n
+
1
∥
2
-
k
∥
μ
n
∥
2
-
2
∑
j
=
1
2
s
d
j
Re
〈
V
j
(
n
)
,
P
j
(
n
)
-
l
V
j
(
n
)
〉
=
-
∑
i
=
1
2
s
+
2
∑
j
=
1
2
s
+
2
M
i
j
〈
r
i
,
r
j
〉
,
where
(23)
r
1
=
w
n
+
1
,
r
2
=
w
n
,
r
j
=
Q
j
-
2
(
n
)
,
j
=
3,4
,
…
,
s
+
2
,
r
j
=
Q
j
-
s
-
2
(
n
-
1
)
,
j
=
s
+
3
,
s
+
4
,
…
,
2
s
+
2
.
By means of (
k
,
l
) algebraical stability of the method, we have the following:
(24)
∥
μ
n
+
1
∥
2
≤
k
∥
μ
n
∥
2
+
2
∑
j
=
1
2
s
d
j
Re
〈
V
j
(
n
)
,
P
j
(
n
)
-
l
V
j
(
n
)
〉
.
It follows from (2)–(4) that
(25)
2
Re
〈
W
j
(
n
)
,
Q
j
(
n
)
〉
=
2
h
Re
〈
W
j
(
n
)
,
f
(
t
n
+
c
j
h
,
Y
j
(
n
)
,
Y
-
j
(
n
)
,
Y
~
j
(
n
)
)
-
f
(
t
n
+
c
j
h
,
Z
j
(
n
)
,
Z
-
j
(
n
)
,
Z
~
j
(
n
)
)
〉
=
2
h
Re
〈
W
j
(
n
)
,
f
(
t
n
+
c
j
h
,
Y
j
(
n
)
,
Y
-
j
(
n
)
,
Y
~
j
(
n
)
)
-
f
(
t
n
+
c
j
h
,
Z
j
(
n
)
,
Y
-
j
(
n
)
,
Y
~
j
(
n
)
)
〉
+
2
h
Re
〈
W
j
(
n
)
,
f
(
t
n
+
c
j
h
,
Z
j
(
n
)
,
Y
-
j
(
n
)
,
Y
~
j
(
n
)
)
-
f
(
t
n
+
c
j
h
,
Z
j
(
n
)
,
Z
-
j
(
n
)
,
Z
~
j
(
n
)
)
〉
≤
2
h
α
∥
W
j
(
n
)
∥
2
+
2
h
β
∥
W
j
(
n
)
∥
·
∥
W
-
j
(
n
)
∥
+
2
h
σ
∥
W
j
(
n
)
∥
·
∥
W
~
j
(
n
)
∥
≤
2
h
α
∥
W
j
(
n
)
∥
2
+
h
β
(
∥
W
j
(
n
)
∥
2
+
∥
W
-
j
(
n
)
∥
2
)
+
h
σ
(
∥
W
j
(
n
)
∥
2
+
∥
W
~
j
(
n
)
∥
2
)
=
(
2
h
α
+
h
β
+
h
σ
)
∥
W
j
(
n
)
∥
2
+
h
β
∥
W
-
j
(
n
)
∥
2
+
h
σ
∥
W
~
j
(
n
)
∥
2
=
(
2
h
α
+
h
β
+
h
σ
)
∥
W
j
(
n
)
∥
2
+
h
β
∥
W
-
j
(
n
)
∥
2
+
h
σ
∥
h
∑
q
=
0
m
p
~
q
g
(
t
n
+
c
j
h
,
t
n
-
q
+
c
j
h
,
Y
j
(
n
-
q
)
)
-
h
∑
q
=
0
m
p
~
q
g
(
t
n
+
c
j
h
,
t
n
-
q
+
c
j
h
,
Z
j
(
n
-
q
)
)
∥
2
≤
(
2
h
α
+
h
β
+
h
σ
)
∥
W
j
(
n
)
∥
2
+
h
β
∥
W
j
(
n
-
m
)
∥
2
+
h
σ
(
m
+
1
)
h
2
p
~
q
2
×
∑
q
=
0
m
∥
g
(
t
n
+
c
j
h
,
t
n
-
q
+
c
j
h
,
Y
j
(
n
-
q
)
)
-
g
(
t
n
+
c
j
h
,
t
n
-
q
+
c
j
h
,
Z
j
(
n
-
q
)
)
∥
2
≤
(
2
h
α
+
h
β
+
h
σ
)
∥
W
j
(
n
)
∥
2
+
h
β
∥
W
j
(
n
-
m
)
∥
2
+
h
σ
(
m
+
1
)
h
2
γ
2
×
∑
q
=
0
m
p
~
q
2
·
∥
W
j
(
n
-
q
)
∥
2
.
Substituting (25) into (24), using (14) we get the following:
(26)
∥
μ
n
+
1
∥
2
≤
∥
μ
n
∥
2
+
∑
j
=
1
2
s
d
j
{
∑
q
=
0
m
∥
W
j
(
n
-
q
)
∥
2
(
2
h
α
+
h
β
+
h
σ
-
2
l
)
∥
W
j
(
n
)
∥
2
+
h
β
∥
W
j
(
n
-
m
)
∥
+
h
σ
γ
2
p
~
2
∑
q
=
0
m
∥
W
j
(
n
-
q
)
∥
2
}
.
By induction, we get the following:
(27)
∥
μ
n
+
1
∥
2
≤
∥
μ
-
1
∥
2
+
∑
i
=
-
1
n
∑
j
=
1
2
s
d
j
{
(
2
h
α
+
2
h
β
+
h
σ
-
2
l
+
h
σ
γ
2
p
~
2
)
∥
W
j
(
i
)
∥
2
}
+
∑
i
=
-
m
-
2
∑
j
=
1
2
s
d
j
{
(
h
β
+
h
σ
γ
2
p
~
2
)
∥
W
j
(
i
)
∥
2
}
≤
∥
μ
-
1
∥
2
+
2
∑
i
=
-
1
n
∑
j
=
1
2
s
d
j
{
h
(
α
+
β
+
σ
+
σ
γ
2
p
~
2
-
l
)
∥
W
j
(
i
)
∥
2
}
+
∑
i
=
-
m
-
1
-
2
∑
j
=
1
2
s
d
j
{
(
h
β
+
h
σ
γ
2
p
~
2
)
∥
W
j
(
i
)
∥
2
}
≤
∥
μ
-
1
∥
2
+
m
h
∑
j
=
1
2
s
d
j
(
β
+
σ
γ
2
p
~
2
)
max
-
m
-
1
≤
i
≤
-
2
∥
W
j
(
i
)
∥
2
≤
∥
w
0
∥
2
+
∥
w
-
1
∥
2
+
m
h
∑
j
=
1
2
s
d
j
(
β
+
σ
γ
2
p
~
2
)
max
-
m
-
1
≤
i
≤
-
2
∥
W
j
(
i
)
∥
2
≤
(
2
+
τ
(
β
+
σ
γ
2
p
~
2
)
∑
j
=
1
2
s
d
j
)
max
t
∈
[
-
τ
,
0
]
∥
φ
(
t
)
-
ψ
(
t
)
∥
2
.
Hence, ∥
μ
n
+
1
∥
≤
C
max
t
∈
[
-
τ
,
0
]
∥
φ
(
t
)
-
ψ
(
t
)
∥, where C
=
2
+
τ
(
β
+
σ
γ
2
p
~
2
)
∑
j
=
1
2
s
d
j.
The proof of Theorem3 is completed.Theorem 4.
Assume that a TSRK method (7a)–(7c) is (
k
,
l
)-algebraically stable with 0
<
k
<
1, then the TSRK method (7a)–(7c) with (9), (12) and (14) is asymptotically stable for the problem (1) with (2)–(4), when the following condition holds:
(28)
h
(
α
+
β
+
1
2
σ
+
1
2
σ
γ
2
p
~
2
)
<
l
.Proof.
Let(29)
u
=
(
2
α
+
β
+
σ
)
h
-
2
l
,
k
-
=
max
{
k
,
(
β
+
σ
γ
2
p
~
2
)
h
(
-
u
)
1
/
m
}
.
Then when h
(
α
+
β
+
(
1
/
2
)
σ
+
(
1
/
2
)
σ
γ
2
p
~
2
)
<
l, we have u
<
-
(
β
+
σ
γ
2
p
~
2
)
h and 0
<
k
-
<
1.
The application of Theorem3 yields the following:
(30)
∥
μ
n
+
1
∥
2
≤
k
-
∥
μ
n
∥
2
+
∑
j
=
1
2
s
d
j
{
∑
q
=
0
m
∥
W
j
(
n
-
q
)
∥
2
(
2
h
α
+
h
β
+
h
σ
-
2
l
)
∥
W
j
(
n
)
∥
2
+
h
β
∥
W
j
(
n
-
m
)
∥
+
h
σ
γ
2
p
~
2
∑
q
=
0
m
∥
W
j
(
n
-
q
)
∥
2
}
≤
k
-
∥
μ
n
∥
2
+
∑
j
=
1
2
s
d
j
{
∑
q
=
0
m
∥
W
j
(
n
-
q
)
∥
2
u
∥
W
j
(
n
)
∥
2
+
h
β
∥
W
j
(
n
-
m
)
∥
+
h
σ
γ
2
p
~
2
∑
q
=
0
m
∥
W
j
(
n
-
q
)
∥
2
}
.
By induction, we get the following:
(31)
∥
μ
n
+
1
∥
2
≤
k
-
n
+
2
∥
μ
-
1
∥
2
+
∑
i
=
-
1
n
k
-
n
-
m
-
i
×
∑
j
=
1
2
s
d
j
(
u
k
-
m
+
h
β
+
h
σ
γ
2
p
~
2
)
·
∥
W
j
(
i
)
∥
2
+
∑
i
=
-
m
-
1
-
2
∑
j
=
1
2
s
d
j
k
-
n
-
m
-
i
(
h
β
+
h
σ
γ
2
p
~
2
)
·
∥
W
j
(
i
)
∥
2
≤
k
-
n
+
2
∥
μ
-
1
∥
2
+
k
-
n
-
m
+
2
·
m
h
·
∑
j
=
1
2
s
d
j
(
β
+
σ
γ
2
p
~
2
)
·
max
-
m
-
1
≤
i
≤
-
2
∥
W
j
(
i
)
∥
2
(32)
≤
k
-
n
+
2
(
∥
w
-
1
∥
2
+
∥
w
0
∥
2
)
+
k
-
n
-
m
+
2
·
τ
·
∑
j
=
1
2
s
d
j
(
β
+
σ
γ
2
p
~
2
)
·
max
-
m
-
1
≤
i
≤
-
2
∥
W
j
(
i
)
∥
2
≤
{
∑
2
2
2
2
k
-
n
+
2
+
k
-
n
-
m
+
2
·
τ
·
∑
j
=
1
2
s
d
j
(
β
+
σ
γ
2
p
~
2
)
∑
2
2
2
}
·
max
t
∈
[
-
τ
,
0
]
∥
φ
(
t
)
-
ψ
(
t
)
∥
2
.
The inequality together with the knowledge 0
<
k
-
<
1 leads to the following:
(33)
lim
n
→
∞
∥
μ
n
+
1
∥
=
0
.
Because μ
n
=
(
w
n
Τ
,
w
n
-
1
Τ
)
Τ, we can get from (33) that
(34)
lim
n
→
∞
∥
w
n
∥
=
0
.
The proof of the theorem is completed.


## 3.1. Some Concepts

For the stability analysis, we need the compound quadrature formula (9) that satisfies the following condition:
(14)
h
2
(
m
+
1
)
2
∑
q
=
0
m
|
p
~
q
|
2
<
p
~
2
,
where p
~ is a positive constant.Let(15)
w
n
=
y
n
-
z
n
,
W
j
(
n
)
=
Y
j
(
n
)
-
Z
j
(
n
)
,
W
~
j
(
n
)
=
Y
~
j
(
n
)
-
Z
~
j
(
n
)
,
Q
j
(
n
)
=
h
[
f
(
t
n
+
c
j
h
,
Y
j
(
n
)
,
Y
-
j
(
n
)
,
Y
~
j
(
n
)
)
-
f
(
t
n
+
c
j
h
,
Z
j
(
n
)
,
Z
-
j
(
n
)
,
Z
~
j
(
n
)
)
]
,
j
=
1,2
,
…
,
s
.
It follows from (8a)–(8c) and (11a)–(11c) that(16a)
W
i
(
n
)
=
w
n
+
∑
j
=
1
s
a
i
j
Q
j
(
n
)
,
i
=
1,2
,
…
,
s
,
(16b)
W
i
(
n
-
1
)
=
w
n
-
1
+
∑
j
=
1
s
a
i
j
Q
j
(
n
-
1
)
,
i
=
1,2
,
…
,
s
,
(16c)
w
n
+
1
=
(
1
-
θ
)
w
n
+
θ
w
n
-
1
+
∑
i
=
1
s
b
i
Q
i
(
n
)
+
∑
i
=
1
s
b
~
i
Q
i
(
n
-
1
)
.Now we will write the s-stage TSRK methods (7a)–(7c) as a general linear method.LetV
i
(
n
)
=
(
W
i
(
n
)
Τ
,
W
i
(
n
-
1
)
Τ
)
Τ be the internal stages and μ
n
=
(
w
n
+
1
Τ
,
w
n
Τ
)
Τ the external vectors and P
i
(
n
)
=
(
Q
i
(
n
)
Τ
,
Q
i
(
n
-
1
)
Τ
)
Τ. Then we have a 2
(
s
+
1
)-stage partitioned general linear method:(17a)
V
i
(
n
)
=
∑
j
=
1
s
C
i
j
11
P
i
(
n
)
+
∑
j
=
1
s
C
i
j
12
μ
n
,
i
=
1,2
,
…
,
s
,
(17b)
μ
n
+
1
=
∑
j
=
1
s
C
i
j
21
P
i
(
n
)
+
∑
j
=
1
s
C
i
j
22
μ
n
,
i
=
1,2
,
…
,
s
,where C
i
j
11
=
a
i
j
,
i
,
j
=
1,2
,
⋯
,
s, C
i
j
11
=
0
,
i
=
1,2
,
⋯
,
s
,
j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, C
i
j
11
=
0, i
=
s
+
1
,
s
+
2
,
⋯
,
2
s, j
=
1,2
,
⋯
,
s, C
i
j
11
=
a
i
j, i
,
j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, C
i
j
12
=
1, i
,
j
=
1,2
,
⋯
,
s, C
i
j
12
=
0, i
=
1,2
,
⋯
,
s, j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, C
i
j
12
=
0, i
=
s
+
1
,
s
+
2
,
⋯
,
2
s, j
=
1,2
,
⋯
,
s, C
i
j
12
=
1, i
,
j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, C
i
j
21
=
b
j, i
,
j
=
1,2
,
⋯
,
s, C
i
j
21
=
b
~
j, i
=
1,2
,
⋯
,
s, j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, C
i
j
21
=
0, i
=
s
+
1
,
s
+
2
,
⋯
,
2
s, j
=
1,2
,
⋯
,
s, C
i
j
21
=
0, i
,
j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, C
i
j
22
=
1
-
θ, i
,
j
=
1,2
,
⋯
,
s, C
i
j
22
=
θ, i
=
1,2
,
⋯
,
s, j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, C
i
j
22
=
1, i
=
s
+
1
,
s
+
2
,
⋯
,
2
s, j
=
1,2
,
⋯
,
s, C
i
j
22
=
0, i
,
j
=
s
+
1
,
s
+
2
,
⋯
,
2
s, and A
=
(
a
i
j
), e
=
(
1,1
,
…
1
), b
=
(
b
1
,
b
2
,
…
,
b
s
), b
~
=
(
b
~
1
,
b
~
2
,
…
,
b
~
s
).We can identify the coefficient matrices(18)
C
11
=
(
C
i
j
11
)
=
(
A
0
0
A
)
,
C
12
=
(
C
i
j
12
)
=
(
e
0
0
e
)
,
C
21
=
(
C
i
j
21
)
=
(
b
b
~
0
0
)
,
C
22
=
(
C
i
j
22
)
=
(
1
-
θ
θ
1
0
)
.Definition 1 (see [22]).
Letk
,
l be real constants, a TSRK method is said to be (
k
,
l
) algebraically stable if there exists a diagonal matrix D
=
diag
(
d
1
,
d
2
,
…
,
d
2
s
) and a diagonal nonnegative matrix G such that M
=
(
m
i
j
) is nonnegative, where(19)
M
(
k
,
l
)
=
(
k
G
-
C
22
Τ
G
C
22
-
2
l
C
12
Τ
D
C
12
C
12
Τ
D
-
C
22
Τ
G
C
21
-
2
l
C
12
Τ
D
C
11
D
C
12
-
C
21
Τ
G
C
22
-
2
l
C
11
Τ
D
C
12
C
11
Τ
D
+
D
C
11
-
C
21
Τ
G
C
21
-
2
l
C
11
Τ
D
C
11
)
.In particular, a (
1,0
)-algebraically stable method is called algebraically stable.Definition 2.
The TSRK method (7a)–(7c) is called asymptotically stable for (1) with (2)–(4) if lim
n
→
∞
∥
w
n
∥
=
0.


## 3.2. Numerical Stability of the Methods

Numerical stability is an important feature of an effective numerical method. An unstable numerical method may be consistent of high order, yet arbitrarily small perturbations will eventually cause large deviations from the true solution. In this section, we will focus on the asymptotic stability of the TSRK method.Theorem 3.
Assume that the TSRK method (7a)–(7c) is (
k
,
l
)-algebraically stable with 0
<
k
≤
1, suppose that the quadrature formula (9) satisfies the condition (14) and the conditions (2)–(4) hold, then, method (17a) and (17b) satisfies the following:
(20)
∥
y
n
-
z
n
∥
≤
C
max
t
∈
[
-
τ
,
0
]
∥
φ
(
t
)
-
ψ
(
t
)
∥
,
n
=
1,2
,
…
,
when the following condition holds:
(21)
h
(
α
+
β
+
σ
+
σ
γ
2
p
~
2
)
≤
l
,
where C depends only on the method, the parameters α, β, σ, γ, η, and τ.Proof.
It follows from a fairly straightforward (but tedious) computation and(
k
,
l
) algebraically stability that (compare also [22, 23])
(22)
∥
μ
n
+
1
∥
2
-
k
∥
μ
n
∥
2
-
2
∑
j
=
1
2
s
d
j
Re
〈
V
j
(
n
)
,
P
j
(
n
)
-
l
V
j
(
n
)
〉
=
-
∑
i
=
1
2
s
+
2
∑
j
=
1
2
s
+
2
M
i
j
〈
r
i
,
r
j
〉
,
where
(23)
r
1
=
w
n
+
1
,
r
2
=
w
n
,
r
j
=
Q
j
-
2
(
n
)
,
j
=
3,4
,
…
,
s
+
2
,
r
j
=
Q
j
-
s
-
2
(
n
-
1
)
,
j
=
s
+
3
,
s
+
4
,
…
,
2
s
+
2
.
By means of (
k
,
l
) algebraical stability of the method, we have the following:
(24)
∥
μ
n
+
1
∥
2
≤
k
∥
μ
n
∥
2
+
2
∑
j
=
1
2
s
d
j
Re
〈
V
j
(
n
)
,
P
j
(
n
)
-
l
V
j
(
n
)
〉
.
It follows from (2)–(4) that
(25)
2
Re
〈
W
j
(
n
)
,
Q
j
(
n
)
〉
=
2
h
Re
〈
W
j
(
n
)
,
f
(
t
n
+
c
j
h
,
Y
j
(
n
)
,
Y
-
j
(
n
)
,
Y
~
j
(
n
)
)
-
f
(
t
n
+
c
j
h
,
Z
j
(
n
)
,
Z
-
j
(
n
)
,
Z
~
j
(
n
)
)
〉
=
2
h
Re
〈
W
j
(
n
)
,
f
(
t
n
+
c
j
h
,
Y
j
(
n
)
,
Y
-
j
(
n
)
,
Y
~
j
(
n
)
)
-
f
(
t
n
+
c
j
h
,
Z
j
(
n
)
,
Y
-
j
(
n
)
,
Y
~
j
(
n
)
)
〉
+
2
h
Re
〈
W
j
(
n
)
,
f
(
t
n
+
c
j
h
,
Z
j
(
n
)
,
Y
-
j
(
n
)
,
Y
~
j
(
n
)
)
-
f
(
t
n
+
c
j
h
,
Z
j
(
n
)
,
Z
-
j
(
n
)
,
Z
~
j
(
n
)
)
〉
≤
2
h
α
∥
W
j
(
n
)
∥
2
+
2
h
β
∥
W
j
(
n
)
∥
·
∥
W
-
j
(
n
)
∥
+
2
h
σ
∥
W
j
(
n
)
∥
·
∥
W
~
j
(
n
)
∥
≤
2
h
α
∥
W
j
(
n
)
∥
2
+
h
β
(
∥
W
j
(
n
)
∥
2
+
∥
W
-
j
(
n
)
∥
2
)
+
h
σ
(
∥
W
j
(
n
)
∥
2
+
∥
W
~
j
(
n
)
∥
2
)
=
(
2
h
α
+
h
β
+
h
σ
)
∥
W
j
(
n
)
∥
2
+
h
β
∥
W
-
j
(
n
)
∥
2
+
h
σ
∥
W
~
j
(
n
)
∥
2
=
(
2
h
α
+
h
β
+
h
σ
)
∥
W
j
(
n
)
∥
2
+
h
β
∥
W
-
j
(
n
)
∥
2
+
h
σ
∥
h
∑
q
=
0
m
p
~
q
g
(
t
n
+
c
j
h
,
t
n
-
q
+
c
j
h
,
Y
j
(
n
-
q
)
)
-
h
∑
q
=
0
m
p
~
q
g
(
t
n
+
c
j
h
,
t
n
-
q
+
c
j
h
,
Z
j
(
n
-
q
)
)
∥
2
≤
(
2
h
α
+
h
β
+
h
σ
)
∥
W
j
(
n
)
∥
2
+
h
β
∥
W
j
(
n
-
m
)
∥
2
+
h
σ
(
m
+
1
)
h
2
p
~
q
2
×
∑
q
=
0
m
∥
g
(
t
n
+
c
j
h
,
t
n
-
q
+
c
j
h
,
Y
j
(
n
-
q
)
)
-
g
(
t
n
+
c
j
h
,
t
n
-
q
+
c
j
h
,
Z
j
(
n
-
q
)
)
∥
2
≤
(
2
h
α
+
h
β
+
h
σ
)
∥
W
j
(
n
)
∥
2
+
h
β
∥
W
j
(
n
-
m
)
∥
2
+
h
σ
(
m
+
1
)
h
2
γ
2
×
∑
q
=
0
m
p
~
q
2
·
∥
W
j
(
n
-
q
)
∥
2
.
Substituting (25) into (24), using (14) we get the following:
(26)
∥
μ
n
+
1
∥
2
≤
∥
μ
n
∥
2
+
∑
j
=
1
2
s
d
j
{
∑
q
=
0
m
∥
W
j
(
n
-
q
)
∥
2
(
2
h
α
+
h
β
+
h
σ
-
2
l
)
∥
W
j
(
n
)
∥
2
+
h
β
∥
W
j
(
n
-
m
)
∥
+
h
σ
γ
2
p
~
2
∑
q
=
0
m
∥
W
j
(
n
-
q
)
∥
2
}
.
By induction, we get the following:
(27)
∥
μ
n
+
1
∥
2
≤
∥
μ
-
1
∥
2
+
∑
i
=
-
1
n
∑
j
=
1
2
s
d
j
{
(
2
h
α
+
2
h
β
+
h
σ
-
2
l
+
h
σ
γ
2
p
~
2
)
∥
W
j
(
i
)
∥
2
}
+
∑
i
=
-
m
-
2
∑
j
=
1
2
s
d
j
{
(
h
β
+
h
σ
γ
2
p
~
2
)
∥
W
j
(
i
)
∥
2
}
≤
∥
μ
-
1
∥
2
+
2
∑
i
=
-
1
n
∑
j
=
1
2
s
d
j
{
h
(
α
+
β
+
σ
+
σ
γ
2
p
~
2
-
l
)
∥
W
j
(
i
)
∥
2
}
+
∑
i
=
-
m
-
1
-
2
∑
j
=
1
2
s
d
j
{
(
h
β
+
h
σ
γ
2
p
~
2
)
∥
W
j
(
i
)
∥
2
}
≤
∥
μ
-
1
∥
2
+
m
h
∑
j
=
1
2
s
d
j
(
β
+
σ
γ
2
p
~
2
)
max
-
m
-
1
≤
i
≤
-
2
∥
W
j
(
i
)
∥
2
≤
∥
w
0
∥
2
+
∥
w
-
1
∥
2
+
m
h
∑
j
=
1
2
s
d
j
(
β
+
σ
γ
2
p
~
2
)
max
-
m
-
1
≤
i
≤
-
2
∥
W
j
(
i
)
∥
2
≤
(
2
+
τ
(
β
+
σ
γ
2
p
~
2
)
∑
j
=
1
2
s
d
j
)
max
t
∈
[
-
τ
,
0
]
∥
φ
(
t
)
-
ψ
(
t
)
∥
2
.
Hence, ∥
μ
n
+
1
∥
≤
C
max
t
∈
[
-
τ
,
0
]
∥
φ
(
t
)
-
ψ
(
t
)
∥, where C
=
2
+
τ
(
β
+
σ
γ
2
p
~
2
)
∑
j
=
1
2
s
d
j.
The proof of Theorem3 is completed.Theorem 4.
Assume that a TSRK method (7a)–(7c) is (
k
,
l
)-algebraically stable with 0
<
k
<
1, then the TSRK method (7a)–(7c) with (9), (12) and (14) is asymptotically stable for the problem (1) with (2)–(4), when the following condition holds:
(28)
h
(
α
+
β
+
1
2
σ
+
1
2
σ
γ
2
p
~
2
)
<
l
.Proof.
Let(29)
u
=
(
2
α
+
β
+
σ
)
h
-
2
l
,
k
-
=
max
{
k
,
(
β
+
σ
γ
2
p
~
2
)
h
(
-
u
)
1
/
m
}
.
Then when h
(
α
+
β
+
(
1
/
2
)
σ
+
(
1
/
2
)
σ
γ
2
p
~
2
)
<
l, we have u
<
-
(
β
+
σ
γ
2
p
~
2
)
h and 0
<
k
-
<
1.
The application of Theorem3 yields the following:
(30)
∥
μ
n
+
1
∥
2
≤
k
-
∥
μ
n
∥
2
+
∑
j
=
1
2
s
d
j
{
∑
q
=
0
m
∥
W
j
(
n
-
q
)
∥
2
(
2
h
α
+
h
β
+
h
σ
-
2
l
)
∥
W
j
(
n
)
∥
2
+
h
β
∥
W
j
(
n
-
m
)
∥
+
h
σ
γ
2
p
~
2
∑
q
=
0
m
∥
W
j
(
n
-
q
)
∥
2
}
≤
k
-
∥
μ
n
∥
2
+
∑
j
=
1
2
s
d
j
{
∑
q
=
0
m
∥
W
j
(
n
-
q
)
∥
2
u
∥
W
j
(
n
)
∥
2
+
h
β
∥
W
j
(
n
-
m
)
∥
+
h
σ
γ
2
p
~
2
∑
q
=
0
m
∥
W
j
(
n
-
q
)
∥
2
}
.
By induction, we get the following:
(31)
∥
μ
n
+
1
∥
2
≤
k
-
n
+
2
∥
μ
-
1
∥
2
+
∑
i
=
-
1
n
k
-
n
-
m
-
i
×
∑
j
=
1
2
s
d
j
(
u
k
-
m
+
h
β
+
h
σ
γ
2
p
~
2
)
·
∥
W
j
(
i
)
∥
2
+
∑
i
=
-
m
-
1
-
2
∑
j
=
1
2
s
d
j
k
-
n
-
m
-
i
(
h
β
+
h
σ
γ
2
p
~
2
)
·
∥
W
j
(
i
)
∥
2
≤
k
-
n
+
2
∥
μ
-
1
∥
2
+
k
-
n
-
m
+
2
·
m
h
·
∑
j
=
1
2
s
d
j
(
β
+
σ
γ
2
p
~
2
)
·
max
-
m
-
1
≤
i
≤
-
2
∥
W
j
(
i
)
∥
2
(32)
≤
k
-
n
+
2
(
∥
w
-
1
∥
2
+
∥
w
0
∥
2
)
+
k
-
n
-
m
+
2
·
τ
·
∑
j
=
1
2
s
d
j
(
β
+
σ
γ
2
p
~
2
)
·
max
-
m
-
1
≤
i
≤
-
2
∥
W
j
(
i
)
∥
2
≤
{
∑
2
2
2
2
k
-
n
+
2
+
k
-
n
-
m
+
2
·
τ
·
∑
j
=
1
2
s
d
j
(
β
+
σ
γ
2
p
~
2
)
∑
2
2
2
}
·
max
t
∈
[
-
τ
,
0
]
∥
φ
(
t
)
-
ψ
(
t
)
∥
2
.
The inequality together with the knowledge 0
<
k
-
<
1 leads to the following:
(33)
lim
n
→
∞
∥
μ
n
+
1
∥
=
0
.
Because μ
n
=
(
w
n
Τ
,
w
n
-
1
Τ
)
Τ, we can get from (33) that
(34)
lim
n
→
∞
∥
w
n
∥
=
0
.
The proof of the theorem is completed.


## 4. The Convergence of TSRK Method for NDDEs

### 4.1. Some Concepts

In order to study the convergence of the method, we define the following:(35)
Y
(
n
)
=
(
Y
1
(
n
)
Τ
,
…
,
Y
s
(
n
)
Τ
,
Y
1
(
n
-
1
)
Τ
,
…
,
Y
s
(
n
-
1
)
Τ
)
Τ
∈
C
2
s
N
,
Y
-
(
n
)
=
(
Y
-
1
(
n
)
Τ
,
…
,
Y
-
s
(
n
)
Τ
,
Y
-
1
(
n
-
1
)
Τ
,
…
,
Y
-
s
(
n
-
1
)
Τ
)
Τ
∈
C
2
s
N
,
Y
~
(
n
)
=
(
Y
~
1
(
n
)
Τ
,
…
,
Y
~
s
(
n
)
Τ
,
Y
~
1
(
n
-
1
)
Τ
,
…
,
Y
~
s
(
n
-
1
)
Τ
)
Τ
∈
C
2
s
N
,
ζ
(
n
+
1
)
=
(
y
(
n
+
1
)
Τ
,
y
(
n
)
Τ
)
Τ
∈
C
2
N
,
F
(
t
n
,
Y
(
n
)
,
Y
-
(
n
)
,
Y
~
(
n
)
)
=
(
f
(
t
n
+
c
1
h
,
Y
1
(
n
)
,
Y
-
1
(
n
)
,
Y
~
1
(
n
)
)
Τ
,
…
,
f
(
t
n
+
c
s
h
,
Y
s
(
n
)
,
Y
-
s
(
n
)
,
Y
~
s
(
n
)
)
Τ
,
f
(
t
n
-
1
+
c
1
h
,
Y
1
(
n
-
1
)
,
Y
-
1
(
n
-
1
)
,
Y
~
1
(
n
-
1
)
)
Τ
,
…
,
f
(
t
n
-
1
+
c
s
h
,
Y
s
(
n
-
1
)
,
Y
-
s
(
n
-
1
)
,
Y
~
s
(
n
-
1
)
)
Τ
)
Τ
∈
C
2
s
N
.
Thus, process (8a)–(8c) can be written in the more compact form:(36a)
Y
(
n
)
=
h
C
11
F
(
t
n
,
Y
(
n
)
,
Y
-
(
n
)
,
Y
~
(
n
)
)
+
C
12
ζ
(
n
)
,
(36b)
ζ
(
n
+
1
)
=
h
C
21
F
(
t
n
,
Y
(
n
)
,
Y
-
(
n
)
,
Y
~
(
n
)
)
+
C
22
ζ
(
n
)
.Definition 5.
Method (7a)–(7c) with an approximation procedure (9) is said to be D convergent of order p if the global error satisfies a bound of the form
(37)
∥
H
(
t
n
)
-
ζ
(
n
)
∥
≤
ρ
1
(
t
n
)
(
∥
H
(
t
0
)
-
ζ
(
0
)
∥
+
max
k
≤
0
∥
Y
(
t
k
)
-
Y
(
k
)
∥
+
h
p
)
,
n
≥
1
,
h
∈
(
0
,
h
0
]
,
where H
(
t
) and Y
(
t
) are defined by
(38)
H
(
t
)
=
(
y
(
t
+
h
)
,
y
(
t
)
)
∈
C
2
N
,
Y
(
t
)
=
(
y
(
t
+
c
1
h
)
,
…
,
y
(
t
+
c
s
h
)
,
y
(
t
-
h
+
c
1
h
)
,
…
,
y
(
t
-
h
+
c
s
h
)
)
∈
C
2
s
N
,
where ρ
1
(
t
) and h
0 depend on M
i, α
,
β
,
γ
,
σ, and τ.Definition 6.
TSRK Method (7a)–(7c) is said to be diagonally stable if there exist a 2
s
×
2
s diagonal matrix Q
>
0 such that the matrix Q
C
11
+
C
11
Τ
Q is positive definite.Remark 7.
The concepts of algebraic stability and diagonal stability of TSRK method are the generalizations of corresponding concepts of Runge-Kutta methods. Although it is difficult to examine these conditions, many results have been found, especially, there exist algebraically stable and diagonally stable multistep formulas of arbitrarily high order (cf. [24]).Definition 8.
TSRK Method (7a)–(7c) is said to have generalized stage order p if p is the largest integer which possesses the following properties.
For any given problem (1) and h
∈
[
0
,
h
-
0
], there exists an abstract function H
h
(
t
),
(39)
H
h
(
t
)
=
(
H
1
h
(
t
)
,
H
2
h
(
t
)
)
∈
C
2
N
such that
(40)
∥
H
(
t
)
=
H
h
(
t
)
∥
≤
p
1
h
p
,
∥
Δ
h
(
t
)
∥
≤
p
1
h
p
+
1
,
∥
δ
h
(
t
)
∥
≤
p
1
h
p
+
1
,
where the maximum stepsize h
-
0
>
0 and the constant p
1 depend only on the method and the bounds M
i, Δ
h
(
t
) and δ
h
(
t
), they are defined by the following equations:(41a)
Y
(
t
)
=
h
C
11
Y
′
(
t
)
+
C
12
H
h
(
t
)
+
Δ
h
(
t
)
,
(41b)
H
h
(
t
+
h
)
=
h
C
21
Y
′
(
t
)
+
C
22
H
h
(
t
)
+
δ
h
(
t
)
.
The function Y
′
(
t
) is defined by the following:
(41c)
Y
(
t
)
=
(
y
′
(
t
+
c
1
h
)
,
…
,
y
′
(
t
+
c
s
h
)
,
y
′
(
t
-
h
+
c
1
h
)
,
…
,
y
′
(
t
-
h
+
c
s
h
)
)
∈
C
2
s
N
.Particularly, when H
(
t
)
=
H
h
(
t
), generalized stage order is called stage order.


### 4.2.D-Convergence and Proofs

In this section, we focus on the error analysis of TSRK method for (1). For the sake of simplicity, we always assume that all constants h
i, γ
i, d
i, and L
i are dependent on the method, the bounds M
i, N
i, the parameters α
,
β
,
γ
,
σ, and τ.First, we give a preliminary result which will later be used several times. To simplify, we denoteY
,
Y
-
,
Y
~
,
Z
,
Z
-
,
Z
~
∈
C
2
s
N for Y
(
n
)
,
Y
-
(
n
)
,
Y
~
(
n
)
,
Z
(
n
)
,
Z
-
(
n
)  Z
~
(
n
), and ζ
,
ϖ
∈
C
2
N for ζ
(
n
+
1
), ϖ
(
n
+
1
), where ϖ
(
n
+
1
)
=
(
z
(
n
+
1
)
,
z
(
n
)
). Define Δ
~ and δ
~ by the following:(42a)
Δ
~
=
Y
-
Z
-
h
C
11
(
F
(
t
,
Y
,
Y
-
,
Y
~
)
-
F
(
t
,
Z
,
Z
-
,
Z
~
)
)
,
(42b)
δ
~
=
ζ
-
ϖ
-
h
C
21
(
F
(
t
,
Y
,
Y
-
,
Y
~
)
-
F
(
t
,
Z
,
Z
-
,
Z
~
)
)
.Theorem 9.
Suppose method (7a)–(7c) is diagonally stable, then there exist constants h
1, p
2, and p
3 such that(43a)
∥
Y
-
Z
∥
≤
p
2
(
∥
Δ
~
∥
+
h
∥
Y
-
-
Z
-
∥
+
h
∥
Y
~
-
Z
~
∥
)
,
h
∈
(
0
,
h
1
]
,
(43b)
∥
ζ
-
ϖ
∥
≤
p
3
(
∥
Δ
~
∥
+
∥
δ
~
∥
+
h
∥
Y
-
-
Z
-
∥
+
h
∥
Y
~
-
Z
~
∥
)
,
h
∈
(
0
,
h
1
]
.Proof.
Since the method (7a)–(7c) is diagonally stable, there exists a positive definite diagonal matrix Q such that the matrix E
=
Q
C
11
+
C
11
Τ
Q is positive definite. Therefore, the matrix C
11 is obviously nonsingular and there exists a l
>
0 depends only on the method such that the matrix
(44)
E
l
=
C
11
-
Τ
E
C
11
-
2
l
Q
is also positive definite.
Define(45a)
W
=
Y
-
Z
,
W
-
=
Y
-
-
Z
-
,
W
~
=
Y
~
-
Z
~
,
(45b)
K
1
=
h
(
F
(
t
,
Y
,
Y
-
,
Y
~
)
-
F
(
t
,
Z
,
Y
-
,
Y
~
)
)
,
(45c)
K
2
=
h
(
F
(
t
,
Z
,
Y
-
,
Y
~
)
-
F
(
t
,
Z
,
Z
-
,
Y
~
)
)
,
(45d)
K
3
=
h
(
F
(
t
,
Z
,
Z
-
,
Y
~
)
-
F
(
t
,
Z
,
Z
-
,
Z
~
)
)
;then(46a)
Δ
~
=
W
-
C
11
(
K
1
+
K
2
+
K
3
)
,
(46b)
δ
~
=
ζ
-
ϖ
-
C
21
(
K
1
+
K
2
+
K
3
)
.Using (2)–(4), (44), (46a), and (46b), we have for h
α
≤
l,
(47)
0
≤
〈
W
,
2
l
Q
W
〉
+
2

Re

〈
K
1
,
-
Q
W
〉
=
-
〈
W
,
E
l
W
〉
+
2

Re
〈
W
,
Q
(
C
11
-
1
W
-
K
1
)
〉
≤
-
λ
1
∥
W
∥
2
+
2

Re
〈
W
,
Q
(
C
11
-
1
Δ
~
+
K
2
+
K
3
)
〉
≤
-
λ
l
∥
W
∥
2
+
2
∥
Q
C
11
-
1
∥
∥
W
∥
·
∥
Δ
~
∥
+
2
h
β
∥
Q
∥
·
∥
W
∥
·
∥
W
-
∥
+
2
h
σ
∥
Q
∥
·
∥
W
∥
·
∥
W
~
∥
,
where λ
l is the minimum eigenvalue of E
l. Therefore,
(48)
∥
Y
-
Z
∥
≤
p
2
(
∥
Δ
~
∥
+
h
∥
Y
-
-
Z
-
∥
+
h
∥
Y
~
-
Z
~
∥
)
,
h
∈
(
0
,
h
1
]
,
where(49a)
p
2
=
2
λ
1
max
{
∥
Q
C
11
-
1
∥
,
β
∥
Q
∥
,
σ
∥
Q
∥
}
,
(49b)
h
1
=
{
1
α
≤
0
,
min
(
1
,
l
α
)
α
>
0
.From (43a), (49a), and (49b), it follows that
(50)
∥
ζ
-
ϖ
∥
=
∥
δ
~
+
C
21
C
11
-
1
(
W
-
Δ
~
)
∥
≤
∥
δ
~
∥
+
p
1
h
∥
C
21
C
11
-
1
∥
(
∥
W
-
∥
+
∥
W
~
∥
)
+
(
1
+
p
1
)
∥
C
21
C
11
-
1
∥
·
∥
Δ
~
∥
≤
p
3
(
∥
Δ
~
∥
+
∥
δ
~
∥
+
h
∥
Y
-
-
Z
-
∥
+
h
∥
Y
~
-
Z
~
∥
)
,
h
∈
(
0
,
h
1
]
,
where p
3
=
max
{
1
,
(
1
+
p
2
)
∥
C
21
C
11
-
1
∥
}, which completes the proof of Theorem 9.Consider the compact form of (11a)–(11c):(51a)
Z
(
n
)
=
h
C
11
F
(
t
n
,
Z
(
n
)
,
Z
-
(
n
)
,
Z
~
(
n
)
)
+
C
12
ϖ
(
n
)
,
(51b)
ϖ
(
n
+
1
)
=
h
C
21
F
(
t
n
,
Z
(
n
)
,
Z
-
(
n
)
,
Z
~
(
n
)
)
+
C
22
ϖ
(
n
)
,where
(52)
Z
(
n
)
=
(
Z
1
(
n
)
,
…
,
Z
s
(
n
)
,
Z
1
(
n
-
1
)
,
…
,
Z
s
(
n
-
1
)
)
∈
C
2
s
N
,
Z
-
(
n
)
=
(
Z
-
1
(
n
)
,
…
,
Z
-
s
(
n
)
,
Z
-
1
(
n
-
1
)
,
…
,
Z
-
s
(
n
-
1
)
)
∈
C
2
s
N
,
Z
~
(
n
)
=
(
Z
~
1
(
n
)
,
…
,
Z
~
s
(
n
)
,
Z
~
1
(
n
-
1
)
,
…
,
Z
~
s
(
n
-
1
)
)
∈
C
2
s
N
,
ϖ
(
n
+
1
)
=
(
z
(
n
+
1
)
,
z
(
n
)
)
∈
C
2
N
.Theorem 10.
Suppose the method (7a)–(7c) is algebraically stable for the matrices G and D, then for (36a), (36b), (51a), and (51b), one has the following:
(53)
∥
ζ
(
n
+
1
)
-
ϖ
(
n
+
1
)
∥
G
2
≤
∥
ζ
(
n
)
-
ϖ
(
n
)
∥
G
2
+
p
4
h
(
∥
Y
~
(
n
)
-
Z
~
(
n
)
∥
2
∥
Y
(
n
)
-
Z
(
n
)
∥
2
+
∥
Y
(
n
-
m
)
-
Z
(
n
-
m
)
∥
2
+
∥
Y
~
(
n
)
-
Z
~
(
n
)
∥
2
)
,
where p
4
=
∥
D
∥
·
max
{
(
2
α
+
β
+
σ
)
,
β
,
σ
}, ∥
·
∥ is a norm on C
2
s
N defined by the following:
(54)
∥
U
∥
G
=
〈
U
,
G
U
〉
1
/
2
=
(
∑
i
,
j
=
1
2
s
g
i
j
〈
u
i
,
u
j
〉
)
1
/
2
,
U
=
(
u
1
,
u
2
,
…
,
u
2
s
)
∈
C
2
s
N
,
u
i
∈
C
N
.Proof.
Defineu
(
n
)
=
ζ
(
n
)
-
ϖ
(
n
), we get from (45a)–(45d) that
(55)
K
(
n
)
=
h
(
F
(
t
n
,
Y
(
n
)
,
Y
-
(
n
)
,
Y
~
(
n
)
)
-
F
(
t
n
,
Z
(
n
)
,
Z
-
(
n
)
,
Z
~
(
n
)
)
)
,
W
(
n
)
=
C
11
K
(
n
)
+
C
12
u
(
n
)
,
u
(
n
+
1
)
=
C
21
K
(
n
)
+
C
22
u
(
n
)
.
With algebraic stability, the matrix
(56)
M
=
[
G
-
C
22
Τ
G
C
22
C
12
Τ
D
-
C
22
Τ
G
C
21
D
C
12
-
C
21
Τ
G
C
22
D
C
11
+
C
11
Τ
D
-
C
21
Τ
G
C
21
]
is nonnegative definite. As in [25], we have
(57)
〈
u
(
n
+
1
)
,
G
u
(
n
+
1
)
〉
-
〈
u
(
n
)
,
G
u
(
n
)
〉
-
2

Re
〈
W
(
n
)
,
D
K
(
n
)
〉
=
-
〈
〈
u
(
n
)
,
K
(
n
)
〉

,
M
〈
u
(
n
)
,
K
(
n
)
〉
〉
≤
0
.
Using (2)–(4), we further obtain the following:
(58)
〈
u
(
n
+
1
)
,
G
u
(
n
+
1
)
〉
≤
〈
u
(
n
)
,
G
u
(
n
)
〉
+
2

Re
〈
W
(
n
)
,
D
K
(
n
)
〉
≤
〈
u
(
n
)
,
G
u
(
n
)
〉
+
2

Re
〈
W
(
n
)
,
D
K
1
(
n
)
〉
+
2

Re
〈
W
(
n
)
,
D
K
2
(
n
)
+
D
K
3
(
n
)
〉
≤
〈
u
(
n
)
,
G
u
(
n
)
〉
+
2
h
α
〈
W
(
n
)
,
D
W
(
n
)
〉
+
2
h
β
∥
D
1
/
2
W
(
n
)
∥
·
∥
D
1
/
2
W
-
(
n
)
∥
+
2
h
σ
∥
D
1
/
2
W
(
n
)
∥
·
∥
D
1
/
2
W
~
(
n
)
∥
≤
〈
u
(
n
)
,
G
u
(
n
)
〉
+
2
h
α
〈
W
(
n
)
,
D
W
(
n
)
〉
+
h
β
(
∥
D
∥
·
∥
W
(
n
)
∥
2
+
∥
D
∥
·
∥
W
-
(
n
)
∥
2
)
+
h
σ
(
∥
D
∥
·
∥
W
(
n
)
∥
2
+
∥
D
∥
·
∥
W
~
(
n
)
∥
2
)
≤
〈
u
(
n
)
,
G
u
(
n
)
〉
+
(
2
h
α
+
h
β
+
h
σ
)
∥
D
∥
·
∥
W
(
n
)
∥
2
+
∥
D
∥
(
h
β
∥
W
(
n
-
m
)
∥
2
+
h
σ
∥
W
~
(
n
)
∥
2
)
,
which gives (53). The proof is completed.In the following, we assume that the method (7a)–(7c) has generalized stage order p; that is, there exists a function H
h
(
t
) such that (40) holds. For any >
0, we define Y
^
(
n
) and y
^
(
n
+
1
) by the following:(59a)
Y
^
(
n
)
=
h
C
11
F
(
t
n
,
Y
^
(
n
)
,
Y
(
t
n
-
τ
)
,
Y
~
^
(
n
)
)
+
C
12
H
h
(
t
n
)
,
(59b)
y
^
(
n
+
1
)
=
h
C
21
F
(
t
n
,
Y
^
(
n
)
,
Y
(
t
n
-
τ
)
,
Y
~
^
(
n
)
)
+
C
22
H
h
(
t
n
)
,where
(60)
Y
~
^
(
n
)
=
(
∫
t
n
-
τ
+
c
1
h
t
n
+
c
1
h
g
(
t
,
ξ
,
H
1
h
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
τ
+
c
s
h
t
n
+
c
s
h
g
(
t
,
ξ
,
H
1
h
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
τ
+
c
1
h
t
n
+
c
1
h
g
(
t
,
ξ
,
H
2
h
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
τ
+
c
s
h
t
n
+
c
s
h
g
(
t
,
ξ
,
H
2
h
(
ξ
)
d
ξ
)
=
(
∫
t
n
-
τ
+
c
1
h
t
n
+
c
1
h
g
(
t
,
ξ
,
H
1
h
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
τ
+
c
s
h
t
n
+
c
s
h
g
(
t
,
ξ
,
H
1
h
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
τ
-
h
+
c
1
h
t
n
-
h
+
c
1
h
g
(
t
,
ξ
,
H
1
h
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
h
-
τ
+
c
s
h
t
n
-
h
+
c
s
h
g
(
t
,
ξ
,
H
1
h
(
ξ
)
d
ξ
)
.Theorem 11.
Suppose the method (7a)–(7c) is diagonally stable and its generalized stage order is p, then there exist constants p
5 and h
2 such that
(61)
∑
k
=
1
n
∥
Y
-
(
k
)
-
Y
(
t
k
-
τ
)
∥
2
=
∑
k
=
1
n
∥
Y
(
k
-
m
)
-
Y
(
t
k
-
τ
)
∥
2
≤
p
5
(
∑
k
=
1
n
∥
ζ
(
k
)
-
H
h
(
t
k
)
∥
2
+
C
1
h
(
p
+
1
)
+
C
2
h
2
(
ν
+
1
)
+
C
3
∑
k
=
-
m
+
1
0
∥
Y
(
k
)
-
Y
(
t
k
)
∥
2
)
,
h
∈
(
0
,
h
2
]
,
where p
5
,
C
1
,
C
2, and C
3 depend on p
2
,
s
,
d
1
,
h
1, and ∥
C
12
∥.Proof.
A combination of (36a) and (41a) leads to the following:
(62)
Y
(
k
)
-
Y
(
t
k
)
=
h
C
11
[
F
(
t
k
,
Y
(
k
)
,
Y
-
(
k
)
,
Y
~
(
k
)
)
-
F
(
t
k
,
Y
(
t
k
)
,
Y
-
(
t
k
)
,
Y
~
(
t
k
)
)
(
t
k
,
Y
(
k
)
,
Y
-
(
k
)
,
Y
~
(
k
)
)
]
+
C
12
(
ζ
(
k
)
-
H
h
(
t
k
)
)
-
Δ
h
(
t
)
.
It follows from Theorem 9 that
(63)
∥
Y
(
k
)
-
Y
(
t
k
)
∥
≤
p
2
(
h
∥
Y
(
k
-
m
)
-
Y
(
t
k
-
τ
)
∥
+
h
∥
Y
~
(
k
)
-
Y
~
(
t
k
)
∥
+
∥
C
12
(
ζ
(
k
)
-
H
h
(
t
k
)
)
-
Δ
h
(
t
)
∥
h
∥
Y
(
k
-
m
)
-
Y
(
t
k
-
τ
)
∥
∥
Y
~
(
k
)
-
Y
~
(
t
k
)
∥
)
.
A combination of (4) and (9) gives the following:
(64)
∥
Y
~
(
k
)
-
Y
~
(
t
k
)
∥
≤
∥
(
h
∑
q
=
0
m
p
~
q
g
(
t
k
+
c
1
h
,
t
k
-
q
+
c
1
h
,
Y
1
(
k
-
q
)
)
,
…
,
h
∑
q
=
0
m
p
~
q
g
(
t
k
+
c
s
h
,
t
k
-
q
+
c
s
h
,
Y
s
(
k
-
q
)
)
,
h
∑
q
=
0
m
p
~
q
g
(
t
k
-
1
+
c
1
h
,
t
k
-
q
-
1
+
c
1
h
,
Y
1
(
k
-
q
-
1
)
)
,
…
,
h
∑
q
=
0
m
p
~
q
g
(
t
k
-
1
+
c
s
h
,
t
k
-
q
-
1
+
c
s
h
,
Y
s
(
k
-
q
-
1
)
)
)
-
(
∫
t
n
-
τ
+
c
1
h
t
n
+
c
1
h
g
(
t
,
ξ
,
y
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
τ
+
c
s
h
t
n
+
c
s
h
g
(
t
,
ξ
,
y
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
τ
+
c
1
h
t
n
+
c
1
h
g
(
t
,
ξ
,
y
(
ξ
-
h
)
d
ξ
,
…
,
∫
t
n
-
τ
+
c
s
h
t
n
+
c
s
h
g
(
t
,
ξ
,
y
(
ξ
-
h
)
d
ξ
)
∑
q
=
0
m
∥
≤
2
s
d
~
1
h
ν
,
where d
~
1 and ν depend only on the compound quadrature (CQ) formula (9), which on substitution into (63) gives the following:
(65)
∑
k
=
1
n
∥
Y
(
k
)
-
Y
(
t
k
)
∥
2
≤
3
(
p
1
2
h
2
∑
k
=
1
n
∥
Y
(
k
-
m
)
-
Y
(
t
k
-
τ
)
∥
2
+
p
2
2
∥
C
12
∥
2
∑
k
=
1
n
∥
ζ
(
k
)
-
H
h
(
t
k
)
∥
2
+
p
1
2
·
∥
Δ
h
(
t
k
)
∥
2
+
p
2
2
4
s
2
d
~
1
2
h
2
ν
+
2
∑
k
=
1
n
)
,
∑
k
=
1
n
∥
Y
(
k
-
m
)
-
Y
(
t
k
-
τ
)
∥
2
=
∑
k
=
1
n
∥
Y
(
k
-
m
)
-
Y
(
t
k
-
m
h
)
∥
2
=
∑
i
=
1
-
m
n
-
m
∥
Y
(
i
)
-
Y
(
t
i
)
∥
2
≤
∑
k
=
1
n
∥
Y
(
k
)
-
Y
(
t
k
)
∥
2
+
∑
k
=
1
-
m
0
∥
Y
(
k
)
-
Y
(
t
k
)
∥
2
.
Therefore, there exist p
5 and h
2 such that (61) holds. The proof of Theorem 11 is completed.Theorem 12.
Suppose method (7a)–(7c) is algebraically stable and diagonally stable and its generalized stage order is p. Then the method with quadrature formula (9) is D-convergent of order at least min
{
p
,
ν
}.Proof.
A combination of (36b), (59b), and (53) leads to the following:
(66)
∥
ζ
(
n
+
1
)
-
y
^
(
n
+
1
)
∥
G
2
≤
∥
ζ
(
n
)
-
y
^
(
n
)
∥
G
2
+
p
4
h
×
[
∥
Y
(
n
)
-
Y
^
(
n
)
∥
2
+
∥
Y
(
n
-
m
)
-
Y
(
t
n
-
m
)
∥
2
+
∥
Y
~
(
n
)
-
Y
~
^
(
n
)
∥
2
]
.
Using Theorem 9, we have the following:
(67)
∥
Y
(
n
)
-
Y
^
(
n
)
∥
≤
p
2
[
∥
Y
(
n
-
m
)
-
Y
(
t
n
-
m
)
∥
h
∥
Y
~
(
n
)
-
Y
~
^
(
n
)
∥
∥
C
12
∥
·
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
+
h
∥
Y
(
n
-
m
)
-
Y
(
t
n
-
m
)
∥
+
h
∥
Y
~
(
n
)
-
Y
~
^
(
n
)
∥
]
,
which on substitution into (66) gives the following:
(68)
∥
ζ
(
n
+
1
)
-
y
^
(
n
+
1
)
∥
G
(
2
)
≤
(
1
+
3
h
p
4
p
2
2
∥
C
12
∥
2
λ
2
)
·
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
G
2
+
p
4
(
1
+
3
h
2
p
2
2
)
h
∥
Y
(
n
-
m
)
-
Y
(
t
n
-
m
)
∥
2
+
p
4
(
1
+
3
h
2
p
2
2
)
h
2
s
d
~
2
h
2
ν
,
where λ
2 is the minimum characteristic value of G. On the other hand,
(69)
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
G
2
≤
∥
ζ
(
n
)
-
y
^
(
n
)
∥
G
2
+
∥
y
^
(
n
)
-
H
h
(
t
n
)
∥
G
2
+
2
∥
ζ
(
n
)
-
y
^
(
n
)
∥
G
·
∥
y
^
(
n
)
-
H
h
(
t
n
)
∥
G
≤
(
1
+
h
)
∥
ζ
(
n
)
-
y
^
(
n
)
∥
G
2
+
(
1
+
1
h
)
∥
y
^
(
n
)
-
H
h
(
t
n
)
∥
G
2
.
In view of (41a)–(41c) and (59a) and (59b), the application of Theorem 9 leads to the following:
(70)
∥
y
^
(
n
)
-
H
h
(
t
n
)
∥
≤
p
3
(
∥
Δ
h
(
t
n
)
∥
+
∥
δ
h
(
t
n
)
∥
)
,
h
∈
(
0
,
h
-
0
]
which gives
(71)
∥
y
^
(
n
)
-
H
h
(
t
n
)
∥
G
2
≤
λ
1
p
3
2
(
∥
Δ
h
(
t
n
)
∥
+
∥
δ
h
(
t
n
)
∥
)
2
,
h
∈
(
0
,
h
-
0
]
,
where λ
1 denotes the maximum eigenvalue of the matrix G. A combination of (40) and (66)–(71) leads to the following:
(72)
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
G
2
≤
(
1
+
h
)
[
γ
1
∥
ζ
(
n
-
1
)
-
H
h
(
t
n
-
1
)
∥
G
2
+
γ
2
∥
Y
(
n
-
m
-
1
)
-
Y
(
t
n
=
m
-
1
)
∥
2
+
γ
3
h
2
ν
+
1
∥
ζ
(
n
-
1
)
-
H
h
(
t
n
-
1
)
∥
G
2
]
+
4
(
1
+
1
h
)
λ
1
p
3
2
p
1
2
h
2
p
+
2
,
h
∈
(
0
,
h
3
]
,
where
(73)
h
3
=
min
{
h
-
0
,
h
1
,
h
2
}
≤
1
,
γ
1
=
1
+
3
h
p
4
p
2
2
∥
C
12
∥
2
λ
2
,
γ
2
=
p
4
h
(
1
+
3
p
2
2
)
h
2
,
γ
3
=
2
s
d
~
2
p
4
(
1
+
3
h
2
p
2
2
)
.
Considering Theorem 11, we further obtain the following:
(74)
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
G
2
≤
(
1
+
h
)
n
γ
1
n
∥
ζ
(
0
)
-
H
h
(
t
0
)
∥
G
2
+
(
1
+
h
)
γ
2
×
∑
i
=
1
n
(
1
+
h
)
i
-
1
γ
1
i
-
1
∥
Y
(
n
-
m
-
i
)
-
Y
(
t
n
=
m
-
i
)
∥
G
2
+
[
(
1
+
h
)
γ
3
h
2
ν
+
1
+
4
(
1
+
1
h
)
λ
1
p
3
2
p
1
2
h
2
p
+
2
]
×
∑
i
=
1
n
(
1
+
h
)
i
-
1
γ
1
i
-
1
.
Using discrete Bellman inequality, we have the following:
(75)
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
G
2
≤
{
∑
i
=
1
n
(
1
+
h
)
n
γ
1
n
∥
ζ
(
0
)
-
H
h
(
t
0
)
∥
G
2
+
(
1
+
h
)
γ
2
×
∑
i
=
1
n
(
1
+
h
)
i
-
1
γ
1
i
-
1
∑
j
=
1
-
m
0
∥
Y
(
j
)
-
Y
(
t
j
)
∥
G
2
+
[
(
1
+
h
)
γ
3
h
2
ν
+
1
+
4
(
1
+
1
h
)
λ
1
p
3
2
p
1
2
h
2
p
+
2
]
×
∑
i
=
1
n
(
1
+
h
)
i
-
1
γ
1
i
-
1
}
×
exp
[
(
1
+
h
)
γ
2
∑
i
=
1
n
(
1
+
h
)
i
-
1
γ
1
i
-
1
]
,
where h
∈
(
0
,
h
0
], and h
0
=
min
{
h
-
0
,
h
1
,
h
2
,
h
3
}
≤
1.
Considering∥
H
(
t
)
-
H
h
(
t
n
)
∥
≤
p
0
h
p, we obtain the following:
(76)
∥
ζ
(
n
)
-
H
(
t
n
)
∥
≤
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
+
∥
H
h
(
t
n
)
-
H
(
t
n
)
∥
≤
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
+
p
0
h
p
.
Considering (75) and (76), we can easily conclude that method (7a)–(7c) with quadrature formula (9) is D-convergent of order at least >min
{
p
,
ν
}. The proof is completed.


## 4.1. Some Concepts

In order to study the convergence of the method, we define the following:(35)
Y
(
n
)
=
(
Y
1
(
n
)
Τ
,
…
,
Y
s
(
n
)
Τ
,
Y
1
(
n
-
1
)
Τ
,
…
,
Y
s
(
n
-
1
)
Τ
)
Τ
∈
C
2
s
N
,
Y
-
(
n
)
=
(
Y
-
1
(
n
)
Τ
,
…
,
Y
-
s
(
n
)
Τ
,
Y
-
1
(
n
-
1
)
Τ
,
…
,
Y
-
s
(
n
-
1
)
Τ
)
Τ
∈
C
2
s
N
,
Y
~
(
n
)
=
(
Y
~
1
(
n
)
Τ
,
…
,
Y
~
s
(
n
)
Τ
,
Y
~
1
(
n
-
1
)
Τ
,
…
,
Y
~
s
(
n
-
1
)
Τ
)
Τ
∈
C
2
s
N
,
ζ
(
n
+
1
)
=
(
y
(
n
+
1
)
Τ
,
y
(
n
)
Τ
)
Τ
∈
C
2
N
,
F
(
t
n
,
Y
(
n
)
,
Y
-
(
n
)
,
Y
~
(
n
)
)
=
(
f
(
t
n
+
c
1
h
,
Y
1
(
n
)
,
Y
-
1
(
n
)
,
Y
~
1
(
n
)
)
Τ
,
…
,
f
(
t
n
+
c
s
h
,
Y
s
(
n
)
,
Y
-
s
(
n
)
,
Y
~
s
(
n
)
)
Τ
,
f
(
t
n
-
1
+
c
1
h
,
Y
1
(
n
-
1
)
,
Y
-
1
(
n
-
1
)
,
Y
~
1
(
n
-
1
)
)
Τ
,
…
,
f
(
t
n
-
1
+
c
s
h
,
Y
s
(
n
-
1
)
,
Y
-
s
(
n
-
1
)
,
Y
~
s
(
n
-
1
)
)
Τ
)
Τ
∈
C
2
s
N
.
Thus, process (8a)–(8c) can be written in the more compact form:(36a)
Y
(
n
)
=
h
C
11
F
(
t
n
,
Y
(
n
)
,
Y
-
(
n
)
,
Y
~
(
n
)
)
+
C
12
ζ
(
n
)
,
(36b)
ζ
(
n
+
1
)
=
h
C
21
F
(
t
n
,
Y
(
n
)
,
Y
-
(
n
)
,
Y
~
(
n
)
)
+
C
22
ζ
(
n
)
.Definition 5.
Method (7a)–(7c) with an approximation procedure (9) is said to be D convergent of order p if the global error satisfies a bound of the form
(37)
∥
H
(
t
n
)
-
ζ
(
n
)
∥
≤
ρ
1
(
t
n
)
(
∥
H
(
t
0
)
-
ζ
(
0
)
∥
+
max
k
≤
0
∥
Y
(
t
k
)
-
Y
(
k
)
∥
+
h
p
)
,
n
≥
1
,
h
∈
(
0
,
h
0
]
,
where H
(
t
) and Y
(
t
) are defined by
(38)
H
(
t
)
=
(
y
(
t
+
h
)
,
y
(
t
)
)
∈
C
2
N
,
Y
(
t
)
=
(
y
(
t
+
c
1
h
)
,
…
,
y
(
t
+
c
s
h
)
,
y
(
t
-
h
+
c
1
h
)
,
…
,
y
(
t
-
h
+
c
s
h
)
)
∈
C
2
s
N
,
where ρ
1
(
t
) and h
0 depend on M
i, α
,
β
,
γ
,
σ, and τ.Definition 6.
TSRK Method (7a)–(7c) is said to be diagonally stable if there exist a 2
s
×
2
s diagonal matrix Q
>
0 such that the matrix Q
C
11
+
C
11
Τ
Q is positive definite.Remark 7.
The concepts of algebraic stability and diagonal stability of TSRK method are the generalizations of corresponding concepts of Runge-Kutta methods. Although it is difficult to examine these conditions, many results have been found, especially, there exist algebraically stable and diagonally stable multistep formulas of arbitrarily high order (cf. [24]).Definition 8.
TSRK Method (7a)–(7c) is said to have generalized stage order p if p is the largest integer which possesses the following properties.
For any given problem (1) and h
∈
[
0
,
h
-
0
], there exists an abstract function H
h
(
t
),
(39)
H
h
(
t
)
=
(
H
1
h
(
t
)
,
H
2
h
(
t
)
)
∈
C
2
N
such that
(40)
∥
H
(
t
)
=
H
h
(
t
)
∥
≤
p
1
h
p
,
∥
Δ
h
(
t
)
∥
≤
p
1
h
p
+
1
,
∥
δ
h
(
t
)
∥
≤
p
1
h
p
+
1
,
where the maximum stepsize h
-
0
>
0 and the constant p
1 depend only on the method and the bounds M
i, Δ
h
(
t
) and δ
h
(
t
), they are defined by the following equations:(41a)
Y
(
t
)
=
h
C
11
Y
′
(
t
)
+
C
12
H
h
(
t
)
+
Δ
h
(
t
)
,
(41b)
H
h
(
t
+
h
)
=
h
C
21
Y
′
(
t
)
+
C
22
H
h
(
t
)
+
δ
h
(
t
)
.
The function Y
′
(
t
) is defined by the following:
(41c)
Y
(
t
)
=
(
y
′
(
t
+
c
1
h
)
,
…
,
y
′
(
t
+
c
s
h
)
,
y
′
(
t
-
h
+
c
1
h
)
,
…
,
y
′
(
t
-
h
+
c
s
h
)
)
∈
C
2
s
N
.Particularly, when H
(
t
)
=
H
h
(
t
), generalized stage order is called stage order.


## 4.2.D-Convergence and Proofs

In this section, we focus on the error analysis of TSRK method for (1). For the sake of simplicity, we always assume that all constants h
i, γ
i, d
i, and L
i are dependent on the method, the bounds M
i, N
i, the parameters α
,
β
,
γ
,
σ, and τ.First, we give a preliminary result which will later be used several times. To simplify, we denoteY
,
Y
-
,
Y
~
,
Z
,
Z
-
,
Z
~
∈
C
2
s
N for Y
(
n
)
,
Y
-
(
n
)
,
Y
~
(
n
)
,
Z
(
n
)
,
Z
-
(
n
)  Z
~
(
n
), and ζ
,
ϖ
∈
C
2
N for ζ
(
n
+
1
), ϖ
(
n
+
1
), where ϖ
(
n
+
1
)
=
(
z
(
n
+
1
)
,
z
(
n
)
). Define Δ
~ and δ
~ by the following:(42a)
Δ
~
=
Y
-
Z
-
h
C
11
(
F
(
t
,
Y
,
Y
-
,
Y
~
)
-
F
(
t
,
Z
,
Z
-
,
Z
~
)
)
,
(42b)
δ
~
=
ζ
-
ϖ
-
h
C
21
(
F
(
t
,
Y
,
Y
-
,
Y
~
)
-
F
(
t
,
Z
,
Z
-
,
Z
~
)
)
.Theorem 9.
Suppose method (7a)–(7c) is diagonally stable, then there exist constants h
1, p
2, and p
3 such that(43a)
∥
Y
-
Z
∥
≤
p
2
(
∥
Δ
~
∥
+
h
∥
Y
-
-
Z
-
∥
+
h
∥
Y
~
-
Z
~
∥
)
,
h
∈
(
0
,
h
1
]
,
(43b)
∥
ζ
-
ϖ
∥
≤
p
3
(
∥
Δ
~
∥
+
∥
δ
~
∥
+
h
∥
Y
-
-
Z
-
∥
+
h
∥
Y
~
-
Z
~
∥
)
,
h
∈
(
0
,
h
1
]
.Proof.
Since the method (7a)–(7c) is diagonally stable, there exists a positive definite diagonal matrix Q such that the matrix E
=
Q
C
11
+
C
11
Τ
Q is positive definite. Therefore, the matrix C
11 is obviously nonsingular and there exists a l
>
0 depends only on the method such that the matrix
(44)
E
l
=
C
11
-
Τ
E
C
11
-
2
l
Q
is also positive definite.
Define(45a)
W
=
Y
-
Z
,
W
-
=
Y
-
-
Z
-
,
W
~
=
Y
~
-
Z
~
,
(45b)
K
1
=
h
(
F
(
t
,
Y
,
Y
-
,
Y
~
)
-
F
(
t
,
Z
,
Y
-
,
Y
~
)
)
,
(45c)
K
2
=
h
(
F
(
t
,
Z
,
Y
-
,
Y
~
)
-
F
(
t
,
Z
,
Z
-
,
Y
~
)
)
,
(45d)
K
3
=
h
(
F
(
t
,
Z
,
Z
-
,
Y
~
)
-
F
(
t
,
Z
,
Z
-
,
Z
~
)
)
;then(46a)
Δ
~
=
W
-
C
11
(
K
1
+
K
2
+
K
3
)
,
(46b)
δ
~
=
ζ
-
ϖ
-
C
21
(
K
1
+
K
2
+
K
3
)
.Using (2)–(4), (44), (46a), and (46b), we have for h
α
≤
l,
(47)
0
≤
〈
W
,
2
l
Q
W
〉
+
2

Re

〈
K
1
,
-
Q
W
〉
=
-
〈
W
,
E
l
W
〉
+
2

Re
〈
W
,
Q
(
C
11
-
1
W
-
K
1
)
〉
≤
-
λ
1
∥
W
∥
2
+
2

Re
〈
W
,
Q
(
C
11
-
1
Δ
~
+
K
2
+
K
3
)
〉
≤
-
λ
l
∥
W
∥
2
+
2
∥
Q
C
11
-
1
∥
∥
W
∥
·
∥
Δ
~
∥
+
2
h
β
∥
Q
∥
·
∥
W
∥
·
∥
W
-
∥
+
2
h
σ
∥
Q
∥
·
∥
W
∥
·
∥
W
~
∥
,
where λ
l is the minimum eigenvalue of E
l. Therefore,
(48)
∥
Y
-
Z
∥
≤
p
2
(
∥
Δ
~
∥
+
h
∥
Y
-
-
Z
-
∥
+
h
∥
Y
~
-
Z
~
∥
)
,
h
∈
(
0
,
h
1
]
,
where(49a)
p
2
=
2
λ
1
max
{
∥
Q
C
11
-
1
∥
,
β
∥
Q
∥
,
σ
∥
Q
∥
}
,
(49b)
h
1
=
{
1
α
≤
0
,
min
(
1
,
l
α
)
α
>
0
.From (43a), (49a), and (49b), it follows that
(50)
∥
ζ
-
ϖ
∥
=
∥
δ
~
+
C
21
C
11
-
1
(
W
-
Δ
~
)
∥
≤
∥
δ
~
∥
+
p
1
h
∥
C
21
C
11
-
1
∥
(
∥
W
-
∥
+
∥
W
~
∥
)
+
(
1
+
p
1
)
∥
C
21
C
11
-
1
∥
·
∥
Δ
~
∥
≤
p
3
(
∥
Δ
~
∥
+
∥
δ
~
∥
+
h
∥
Y
-
-
Z
-
∥
+
h
∥
Y
~
-
Z
~
∥
)
,
h
∈
(
0
,
h
1
]
,
where p
3
=
max
{
1
,
(
1
+
p
2
)
∥
C
21
C
11
-
1
∥
}, which completes the proof of Theorem 9.Consider the compact form of (11a)–(11c):(51a)
Z
(
n
)
=
h
C
11
F
(
t
n
,
Z
(
n
)
,
Z
-
(
n
)
,
Z
~
(
n
)
)
+
C
12
ϖ
(
n
)
,
(51b)
ϖ
(
n
+
1
)
=
h
C
21
F
(
t
n
,
Z
(
n
)
,
Z
-
(
n
)
,
Z
~
(
n
)
)
+
C
22
ϖ
(
n
)
,where
(52)
Z
(
n
)
=
(
Z
1
(
n
)
,
…
,
Z
s
(
n
)
,
Z
1
(
n
-
1
)
,
…
,
Z
s
(
n
-
1
)
)
∈
C
2
s
N
,
Z
-
(
n
)
=
(
Z
-
1
(
n
)
,
…
,
Z
-
s
(
n
)
,
Z
-
1
(
n
-
1
)
,
…
,
Z
-
s
(
n
-
1
)
)
∈
C
2
s
N
,
Z
~
(
n
)
=
(
Z
~
1
(
n
)
,
…
,
Z
~
s
(
n
)
,
Z
~
1
(
n
-
1
)
,
…
,
Z
~
s
(
n
-
1
)
)
∈
C
2
s
N
,
ϖ
(
n
+
1
)
=
(
z
(
n
+
1
)
,
z
(
n
)
)
∈
C
2
N
.Theorem 10.
Suppose the method (7a)–(7c) is algebraically stable for the matrices G and D, then for (36a), (36b), (51a), and (51b), one has the following:
(53)
∥
ζ
(
n
+
1
)
-
ϖ
(
n
+
1
)
∥
G
2
≤
∥
ζ
(
n
)
-
ϖ
(
n
)
∥
G
2
+
p
4
h
(
∥
Y
~
(
n
)
-
Z
~
(
n
)
∥
2
∥
Y
(
n
)
-
Z
(
n
)
∥
2
+
∥
Y
(
n
-
m
)
-
Z
(
n
-
m
)
∥
2
+
∥
Y
~
(
n
)
-
Z
~
(
n
)
∥
2
)
,
where p
4
=
∥
D
∥
·
max
{
(
2
α
+
β
+
σ
)
,
β
,
σ
}, ∥
·
∥ is a norm on C
2
s
N defined by the following:
(54)
∥
U
∥
G
=
〈
U
,
G
U
〉
1
/
2
=
(
∑
i
,
j
=
1
2
s
g
i
j
〈
u
i
,
u
j
〉
)
1
/
2
,
U
=
(
u
1
,
u
2
,
…
,
u
2
s
)
∈
C
2
s
N
,
u
i
∈
C
N
.Proof.
Defineu
(
n
)
=
ζ
(
n
)
-
ϖ
(
n
), we get from (45a)–(45d) that
(55)
K
(
n
)
=
h
(
F
(
t
n
,
Y
(
n
)
,
Y
-
(
n
)
,
Y
~
(
n
)
)
-
F
(
t
n
,
Z
(
n
)
,
Z
-
(
n
)
,
Z
~
(
n
)
)
)
,
W
(
n
)
=
C
11
K
(
n
)
+
C
12
u
(
n
)
,
u
(
n
+
1
)
=
C
21
K
(
n
)
+
C
22
u
(
n
)
.
With algebraic stability, the matrix
(56)
M
=
[
G
-
C
22
Τ
G
C
22
C
12
Τ
D
-
C
22
Τ
G
C
21
D
C
12
-
C
21
Τ
G
C
22
D
C
11
+
C
11
Τ
D
-
C
21
Τ
G
C
21
]
is nonnegative definite. As in [25], we have
(57)
〈
u
(
n
+
1
)
,
G
u
(
n
+
1
)
〉
-
〈
u
(
n
)
,
G
u
(
n
)
〉
-
2

Re
〈
W
(
n
)
,
D
K
(
n
)
〉
=
-
〈
〈
u
(
n
)
,
K
(
n
)
〉

,
M
〈
u
(
n
)
,
K
(
n
)
〉
〉
≤
0
.
Using (2)–(4), we further obtain the following:
(58)
〈
u
(
n
+
1
)
,
G
u
(
n
+
1
)
〉
≤
〈
u
(
n
)
,
G
u
(
n
)
〉
+
2

Re
〈
W
(
n
)
,
D
K
(
n
)
〉
≤
〈
u
(
n
)
,
G
u
(
n
)
〉
+
2

Re
〈
W
(
n
)
,
D
K
1
(
n
)
〉
+
2

Re
〈
W
(
n
)
,
D
K
2
(
n
)
+
D
K
3
(
n
)
〉
≤
〈
u
(
n
)
,
G
u
(
n
)
〉
+
2
h
α
〈
W
(
n
)
,
D
W
(
n
)
〉
+
2
h
β
∥
D
1
/
2
W
(
n
)
∥
·
∥
D
1
/
2
W
-
(
n
)
∥
+
2
h
σ
∥
D
1
/
2
W
(
n
)
∥
·
∥
D
1
/
2
W
~
(
n
)
∥
≤
〈
u
(
n
)
,
G
u
(
n
)
〉
+
2
h
α
〈
W
(
n
)
,
D
W
(
n
)
〉
+
h
β
(
∥
D
∥
·
∥
W
(
n
)
∥
2
+
∥
D
∥
·
∥
W
-
(
n
)
∥
2
)
+
h
σ
(
∥
D
∥
·
∥
W
(
n
)
∥
2
+
∥
D
∥
·
∥
W
~
(
n
)
∥
2
)
≤
〈
u
(
n
)
,
G
u
(
n
)
〉
+
(
2
h
α
+
h
β
+
h
σ
)
∥
D
∥
·
∥
W
(
n
)
∥
2
+
∥
D
∥
(
h
β
∥
W
(
n
-
m
)
∥
2
+
h
σ
∥
W
~
(
n
)
∥
2
)
,
which gives (53). The proof is completed.In the following, we assume that the method (7a)–(7c) has generalized stage order p; that is, there exists a function H
h
(
t
) such that (40) holds. For any >
0, we define Y
^
(
n
) and y
^
(
n
+
1
) by the following:(59a)
Y
^
(
n
)
=
h
C
11
F
(
t
n
,
Y
^
(
n
)
,
Y
(
t
n
-
τ
)
,
Y
~
^
(
n
)
)
+
C
12
H
h
(
t
n
)
,
(59b)
y
^
(
n
+
1
)
=
h
C
21
F
(
t
n
,
Y
^
(
n
)
,
Y
(
t
n
-
τ
)
,
Y
~
^
(
n
)
)
+
C
22
H
h
(
t
n
)
,where
(60)
Y
~
^
(
n
)
=
(
∫
t
n
-
τ
+
c
1
h
t
n
+
c
1
h
g
(
t
,
ξ
,
H
1
h
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
τ
+
c
s
h
t
n
+
c
s
h
g
(
t
,
ξ
,
H
1
h
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
τ
+
c
1
h
t
n
+
c
1
h
g
(
t
,
ξ
,
H
2
h
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
τ
+
c
s
h
t
n
+
c
s
h
g
(
t
,
ξ
,
H
2
h
(
ξ
)
d
ξ
)
=
(
∫
t
n
-
τ
+
c
1
h
t
n
+
c
1
h
g
(
t
,
ξ
,
H
1
h
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
τ
+
c
s
h
t
n
+
c
s
h
g
(
t
,
ξ
,
H
1
h
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
τ
-
h
+
c
1
h
t
n
-
h
+
c
1
h
g
(
t
,
ξ
,
H
1
h
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
h
-
τ
+
c
s
h
t
n
-
h
+
c
s
h
g
(
t
,
ξ
,
H
1
h
(
ξ
)
d
ξ
)
.Theorem 11.
Suppose the method (7a)–(7c) is diagonally stable and its generalized stage order is p, then there exist constants p
5 and h
2 such that
(61)
∑
k
=
1
n
∥
Y
-
(
k
)
-
Y
(
t
k
-
τ
)
∥
2
=
∑
k
=
1
n
∥
Y
(
k
-
m
)
-
Y
(
t
k
-
τ
)
∥
2
≤
p
5
(
∑
k
=
1
n
∥
ζ
(
k
)
-
H
h
(
t
k
)
∥
2
+
C
1
h
(
p
+
1
)
+
C
2
h
2
(
ν
+
1
)
+
C
3
∑
k
=
-
m
+
1
0
∥
Y
(
k
)
-
Y
(
t
k
)
∥
2
)
,
h
∈
(
0
,
h
2
]
,
where p
5
,
C
1
,
C
2, and C
3 depend on p
2
,
s
,
d
1
,
h
1, and ∥
C
12
∥.Proof.
A combination of (36a) and (41a) leads to the following:
(62)
Y
(
k
)
-
Y
(
t
k
)
=
h
C
11
[
F
(
t
k
,
Y
(
k
)
,
Y
-
(
k
)
,
Y
~
(
k
)
)
-
F
(
t
k
,
Y
(
t
k
)
,
Y
-
(
t
k
)
,
Y
~
(
t
k
)
)
(
t
k
,
Y
(
k
)
,
Y
-
(
k
)
,
Y
~
(
k
)
)
]
+
C
12
(
ζ
(
k
)
-
H
h
(
t
k
)
)
-
Δ
h
(
t
)
.
It follows from Theorem 9 that
(63)
∥
Y
(
k
)
-
Y
(
t
k
)
∥
≤
p
2
(
h
∥
Y
(
k
-
m
)
-
Y
(
t
k
-
τ
)
∥
+
h
∥
Y
~
(
k
)
-
Y
~
(
t
k
)
∥
+
∥
C
12
(
ζ
(
k
)
-
H
h
(
t
k
)
)
-
Δ
h
(
t
)
∥
h
∥
Y
(
k
-
m
)
-
Y
(
t
k
-
τ
)
∥
∥
Y
~
(
k
)
-
Y
~
(
t
k
)
∥
)
.
A combination of (4) and (9) gives the following:
(64)
∥
Y
~
(
k
)
-
Y
~
(
t
k
)
∥
≤
∥
(
h
∑
q
=
0
m
p
~
q
g
(
t
k
+
c
1
h
,
t
k
-
q
+
c
1
h
,
Y
1
(
k
-
q
)
)
,
…
,
h
∑
q
=
0
m
p
~
q
g
(
t
k
+
c
s
h
,
t
k
-
q
+
c
s
h
,
Y
s
(
k
-
q
)
)
,
h
∑
q
=
0
m
p
~
q
g
(
t
k
-
1
+
c
1
h
,
t
k
-
q
-
1
+
c
1
h
,
Y
1
(
k
-
q
-
1
)
)
,
…
,
h
∑
q
=
0
m
p
~
q
g
(
t
k
-
1
+
c
s
h
,
t
k
-
q
-
1
+
c
s
h
,
Y
s
(
k
-
q
-
1
)
)
)
-
(
∫
t
n
-
τ
+
c
1
h
t
n
+
c
1
h
g
(
t
,
ξ
,
y
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
τ
+
c
s
h
t
n
+
c
s
h
g
(
t
,
ξ
,
y
(
ξ
)
d
ξ
,
…
,
∫
t
n
-
τ
+
c
1
h
t
n
+
c
1
h
g
(
t
,
ξ
,
y
(
ξ
-
h
)
d
ξ
,
…
,
∫
t
n
-
τ
+
c
s
h
t
n
+
c
s
h
g
(
t
,
ξ
,
y
(
ξ
-
h
)
d
ξ
)
∑
q
=
0
m
∥
≤
2
s
d
~
1
h
ν
,
where d
~
1 and ν depend only on the compound quadrature (CQ) formula (9), which on substitution into (63) gives the following:
(65)
∑
k
=
1
n
∥
Y
(
k
)
-
Y
(
t
k
)
∥
2
≤
3
(
p
1
2
h
2
∑
k
=
1
n
∥
Y
(
k
-
m
)
-
Y
(
t
k
-
τ
)
∥
2
+
p
2
2
∥
C
12
∥
2
∑
k
=
1
n
∥
ζ
(
k
)
-
H
h
(
t
k
)
∥
2
+
p
1
2
·
∥
Δ
h
(
t
k
)
∥
2
+
p
2
2
4
s
2
d
~
1
2
h
2
ν
+
2
∑
k
=
1
n
)
,
∑
k
=
1
n
∥
Y
(
k
-
m
)
-
Y
(
t
k
-
τ
)
∥
2
=
∑
k
=
1
n
∥
Y
(
k
-
m
)
-
Y
(
t
k
-
m
h
)
∥
2
=
∑
i
=
1
-
m
n
-
m
∥
Y
(
i
)
-
Y
(
t
i
)
∥
2
≤
∑
k
=
1
n
∥
Y
(
k
)
-
Y
(
t
k
)
∥
2
+
∑
k
=
1
-
m
0
∥
Y
(
k
)
-
Y
(
t
k
)
∥
2
.
Therefore, there exist p
5 and h
2 such that (61) holds. The proof of Theorem 11 is completed.Theorem 12.
Suppose method (7a)–(7c) is algebraically stable and diagonally stable and its generalized stage order is p. Then the method with quadrature formula (9) is D-convergent of order at least min
{
p
,
ν
}.Proof.
A combination of (36b), (59b), and (53) leads to the following:
(66)
∥
ζ
(
n
+
1
)
-
y
^
(
n
+
1
)
∥
G
2
≤
∥
ζ
(
n
)
-
y
^
(
n
)
∥
G
2
+
p
4
h
×
[
∥
Y
(
n
)
-
Y
^
(
n
)
∥
2
+
∥
Y
(
n
-
m
)
-
Y
(
t
n
-
m
)
∥
2
+
∥
Y
~
(
n
)
-
Y
~
^
(
n
)
∥
2
]
.
Using Theorem 9, we have the following:
(67)
∥
Y
(
n
)
-
Y
^
(
n
)
∥
≤
p
2
[
∥
Y
(
n
-
m
)
-
Y
(
t
n
-
m
)
∥
h
∥
Y
~
(
n
)
-
Y
~
^
(
n
)
∥
∥
C
12
∥
·
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
+
h
∥
Y
(
n
-
m
)
-
Y
(
t
n
-
m
)
∥
+
h
∥
Y
~
(
n
)
-
Y
~
^
(
n
)
∥
]
,
which on substitution into (66) gives the following:
(68)
∥
ζ
(
n
+
1
)
-
y
^
(
n
+
1
)
∥
G
(
2
)
≤
(
1
+
3
h
p
4
p
2
2
∥
C
12
∥
2
λ
2
)
·
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
G
2
+
p
4
(
1
+
3
h
2
p
2
2
)
h
∥
Y
(
n
-
m
)
-
Y
(
t
n
-
m
)
∥
2
+
p
4
(
1
+
3
h
2
p
2
2
)
h
2
s
d
~
2
h
2
ν
,
where λ
2 is the minimum characteristic value of G. On the other hand,
(69)
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
G
2
≤
∥
ζ
(
n
)
-
y
^
(
n
)
∥
G
2
+
∥
y
^
(
n
)
-
H
h
(
t
n
)
∥
G
2
+
2
∥
ζ
(
n
)
-
y
^
(
n
)
∥
G
·
∥
y
^
(
n
)
-
H
h
(
t
n
)
∥
G
≤
(
1
+
h
)
∥
ζ
(
n
)
-
y
^
(
n
)
∥
G
2
+
(
1
+
1
h
)
∥
y
^
(
n
)
-
H
h
(
t
n
)
∥
G
2
.
In view of (41a)–(41c) and (59a) and (59b), the application of Theorem 9 leads to the following:
(70)
∥
y
^
(
n
)
-
H
h
(
t
n
)
∥
≤
p
3
(
∥
Δ
h
(
t
n
)
∥
+
∥
δ
h
(
t
n
)
∥
)
,
h
∈
(
0
,
h
-
0
]
which gives
(71)
∥
y
^
(
n
)
-
H
h
(
t
n
)
∥
G
2
≤
λ
1
p
3
2
(
∥
Δ
h
(
t
n
)
∥
+
∥
δ
h
(
t
n
)
∥
)
2
,
h
∈
(
0
,
h
-
0
]
,
where λ
1 denotes the maximum eigenvalue of the matrix G. A combination of (40) and (66)–(71) leads to the following:
(72)
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
G
2
≤
(
1
+
h
)
[
γ
1
∥
ζ
(
n
-
1
)
-
H
h
(
t
n
-
1
)
∥
G
2
+
γ
2
∥
Y
(
n
-
m
-
1
)
-
Y
(
t
n
=
m
-
1
)
∥
2
+
γ
3
h
2
ν
+
1
∥
ζ
(
n
-
1
)
-
H
h
(
t
n
-
1
)
∥
G
2
]
+
4
(
1
+
1
h
)
λ
1
p
3
2
p
1
2
h
2
p
+
2
,
h
∈
(
0
,
h
3
]
,
where
(73)
h
3
=
min
{
h
-
0
,
h
1
,
h
2
}
≤
1
,
γ
1
=
1
+
3
h
p
4
p
2
2
∥
C
12
∥
2
λ
2
,
γ
2
=
p
4
h
(
1
+
3
p
2
2
)
h
2
,
γ
3
=
2
s
d
~
2
p
4
(
1
+
3
h
2
p
2
2
)
.
Considering Theorem 11, we further obtain the following:
(74)
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
G
2
≤
(
1
+
h
)
n
γ
1
n
∥
ζ
(
0
)
-
H
h
(
t
0
)
∥
G
2
+
(
1
+
h
)
γ
2
×
∑
i
=
1
n
(
1
+
h
)
i
-
1
γ
1
i
-
1
∥
Y
(
n
-
m
-
i
)
-
Y
(
t
n
=
m
-
i
)
∥
G
2
+
[
(
1
+
h
)
γ
3
h
2
ν
+
1
+
4
(
1
+
1
h
)
λ
1
p
3
2
p
1
2
h
2
p
+
2
]
×
∑
i
=
1
n
(
1
+
h
)
i
-
1
γ
1
i
-
1
.
Using discrete Bellman inequality, we have the following:
(75)
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
G
2
≤
{
∑
i
=
1
n
(
1
+
h
)
n
γ
1
n
∥
ζ
(
0
)
-
H
h
(
t
0
)
∥
G
2
+
(
1
+
h
)
γ
2
×
∑
i
=
1
n
(
1
+
h
)
i
-
1
γ
1
i
-
1
∑
j
=
1
-
m
0
∥
Y
(
j
)
-
Y
(
t
j
)
∥
G
2
+
[
(
1
+
h
)
γ
3
h
2
ν
+
1
+
4
(
1
+
1
h
)
λ
1
p
3
2
p
1
2
h
2
p
+
2
]
×
∑
i
=
1
n
(
1
+
h
)
i
-
1
γ
1
i
-
1
}
×
exp
[
(
1
+
h
)
γ
2
∑
i
=
1
n
(
1
+
h
)
i
-
1
γ
1
i
-
1
]
,
where h
∈
(
0
,
h
0
], and h
0
=
min
{
h
-
0
,
h
1
,
h
2
,
h
3
}
≤
1.
Considering∥
H
(
t
)
-
H
h
(
t
n
)
∥
≤
p
0
h
p, we obtain the following:
(76)
∥
ζ
(
n
)
-
H
(
t
n
)
∥
≤
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
+
∥
H
h
(
t
n
)
-
H
(
t
n
)
∥
≤
∥
ζ
(
n
)
-
H
h
(
t
n
)
∥
+
p
0
h
p
.
Considering (75) and (76), we can easily conclude that method (7a)–(7c) with quadrature formula (9) is D-convergent of order at least >min
{
p
,
ν
}. The proof is completed.


## 5. Numerical Experiments

Consider the following nonlinear Volterra delay integro-differential equations:(77)
y
′
(
t
)
=
-
(
6
+
sin
t
)
y
(
t
)
+
y
(
t
-
π
4
)
-
∫
t
-
π
/
4
t
sin
(
v
)
y
(
v
)
d
v
+
5
exp
(
cos
t
)
,
t
≥
0
,
y
(
t
)
=
exp
(
cos
t
)
,
-
π
4
≤
t
≤
0
and its perturbed problem
(78)
z
′
(
t
)
=
-
(
6
+
sin
t
)
z
(
t
)
+
z
(
t
-
π
4
)
-
∫
t
-
π
/
4
t
sin
(
v
)
z
(
v
)
d
v
+
5
exp
(
cos
t
)
,
t
≥
0
,
z
(
t
)
=
2
,
-
π
4
≤
t
≤
0
.
One may check that this system has the solution y
(
t
)
=
exp
(
cos
t
) and a
=
-
4, β
=
1
/
2, σ
=
1
/
2, γ
=
1. With Theorem 9, we conclude that system (77) satisfies stability and convergence properties.Apply the two-step Runge-Kutta method induced by the GL method in [12] (79a)to (77) and (78), Y
~
j
(
n
) is computed by the compound Simpson rule with an even integer m
≥
4 formula:
(79b)
Y
~
j
(
n
)
=
h
3
[
sin
(
t
n
+
c
j
h
)
Y
j
(
n
)
+
4
∑
q
=
1
m
/
2
sin
(
t
n
-
2
q
+
1
+
c
j
h
)
Y
j
(
n
-
2
q
+
1
)
+
2
∑
q
=
1
m
-
2
/
2
sin
(
t
n
-
2
q
+
c
j
h
)
Y
j
(
n
-
2
q
)
+
sin
(
t
n
-
m
+
c
j
h
)
Y
j
(
n
-
m
)
+
2
∑
q
=
1
m
-
2
/
2
sin
(
t
n
-
2
q
+
c
j
h
)
Y
j
(
n
-
2
q
)
]Using the condition (28), we know this method is stable; the convergence order is shown in Table 1.For Figures1 and 2, it is obvious that the corresponding method for VDIDEs is stable and convergent, and the convergent order is min
{
4,4
}.Figure 1
The numerical solution of (79a) and (79b) for (77) with h
=
0.1.Figure 2
The numerical solution of (79a) and (79b) for (78) with h
=
0.1.


## 6. Conclusions

This paper is devoted to the stability and convergence analysis of the two-step Runge-Kutta (TSRK) methods with compound quadrature formula for the numerical solution for a nonlinear Volterra delay integro-differential equations. Nonlinear stability andD-convergence are introduced and proved.We believe that the results presented in this paper can be extended to other general DIDEs and NDIDEs. However, it is difficult to extend the results presented in this paper to more general delay differential equations; the discussion will be discussed later.



---

*Source: 679075-2013-05-19.xml*