# The Stability of Two-Step Runge-Kutta Methods for Neutral Delay Integro Differential-Algebraic Equations with Many Delays

Authors: Haiyan Yuan; Jihong Shen
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

This paper studies the asymptotic stability of the two-step Runge-Kutta methods for neutral delay integro differential-algebraic equations with many delays. It proves that A-stable two-step Runge-Kutta methods are asymptotically stable for neutral delay integro differential-algebraic equations with many delays.


---

## Body

## 1. Introduction

The stability of numerical methods for delay differential equations has been intensively studied in [1–3] for many years. These equations appeared in a wide variety of scientific and engineering fields, such as circuit analysis, computer-aided design power systems, and optimal control. The structure for these, the order of convergence, and the asymptotic stability of numerical methods have been studied in [4–6]. Zhu and Petzold investigated the asymptotic stability of neutral delay differential equations with θ-methods, Runge-Kutta methods, BDF methods, and linear multistep methods [7]. Zhao et al. studied the stability of neutral delay differential equations with Rosenbrock methods [8]. Yu et al. studied the general neutral delay differential equations with multistep methods [9]. More recently, there is a growing interest in the analysis of delay integro differential equations. Baker and Ford [10] studied the asymptotic stability of a class of linear multistep (LM) methods for scalar linear delay integro differential equations; Koto [11] dealt with the linear stability of Runge-Kutta (R-K) methods for systems of delay integro differential equations; Huang and Vandewalle [12] gave sufficient and necessary stability conditions for exact and discrete solutions of linear scalar delay integro differential equations, and Luzyanina et al. [13] developed computational procedures for determining the stability of delay integro differential equations. Zhang and Vandewalle [14] gave the stability criteria for exact and discrete solution of neutral multidelay integro differential equations. Although the stability of numerical methods for delay integro differential equations has been very intensively studied, the stability of delay integro differential equations with many delays has not been studied so far.In this paper, we focus on the asymptotic stability of numerical methods for neutral delay integro differential-algebraic equations with many delays. This paper is structured as follows. In Section2 we give asymptotic stability of the analytical solution and introduce two-step Runge-Kutta methods and the stability region. In Section 3, we deal with the asymptotic stability of two-step Runge-Kutta method for neutral delay integro differential-algebraic equations with many delays; the theoretical results are proved. In Section 4, an example is given to illustrate the theoretical results.


## 2. Asymptotic Stability of the Analytical Solution

### 2.1. Asymptotic Stability of the Analytical Solution of Neutral Delay Integro Differential-Algebraic Equation with Many Delays

In this section, we consider the following linear system:(1)
A
u
′
(
t
)
+
B
u
(
t
)
+
∑
q
=
1
m
C
q
u
′
(
t
-
τ
q
)
+
∑
q
=
1
m
D
q
u
(
t
-
τ
q
)
+
∑
q
=
1
m
G
q
∫
t
-
τ
q
t
u
(
δ
)
d
δ
=
0
,
t
≥
0
,
u
(
t
)
=
φ
(
t
)
,
t
∈
[
-
τ
,
0
)
,
where A, B, C
q, D
q, G
q
∈
R
d
×
d, A is a singular matrix, τ
q is a given positive delay constant (q
=
1,2
,
…
,
m), and 0
<
τ
1
≤
τ
2
≤
⋯
≤
τ
m
=
τ. φ
(
t
) denotes a given vector-valued function and u
(
t
) is a vector-valued unknown function to be solved for t
≥
0.In order to obtain the characteristic equation of system (1), we focus on the exponential solutions u
(
t
)
=
e
s
t
x of (1); here x
=
(
x
1
,
x
2
,
…
,
x
d
)
Τ
∈
C
d denotes the unknown vector. Then we have (2a)
u
(
t
)
=
(
e
s
t
x
1
,
e
s
t
x
2
,
…
,
e
s
t
x
d
)
Τ
,
(2b)
u
′
(
t
)
=
s
e
s
t
x
,
(2c)
u
(
t
-
τ
q
)
=
(
e
s
(
t
-
τ
q
)
x
1
,
e
s
(
t
-
τ
q
)
x
2
,
…
,
e
s
(
t
-
τ
q
)
x
d
)
Τ
,
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
q
=
1,2
,
…
,
m
,
(2d)
u
′
(
t
-
τ
q
)
=
s
e
s
(
t
-
τ
q
)
x
,
q
=
1,2
,
…
,
m
.Substituting the above results into (1), we have the following equation:
(2e)
[
s
A
+
B
+
s
∑
q
=
1
m
C
q
e
-
s
τ
q
+
∑
q
=
1
m
D
q
e
-
s
τ
q
i
i
i
+
s
-
1
∑
q
=
1
m
G
q
(
1
-
e
-
s
τ
q
)
]
x
=
0
.The existence of a nonzerox in (2e) implies the characteristic equation of system (1) holds; that is, the following equation holds:
(3)
det
⁡
[
s
A
+
B
+
s
∑
q
=
1
m
C
q
e
-
s
τ
q
i
i
i
i
i
i
i
i
+
∑
q
=
1
m
D
q
e
-
s
τ
q
i
i
i
i
i
i
i
i
+
s
-
1
∑
q
=
1
m
G
q
(
1
-
e
-
s
τ
q
)
]
=
0
.Definition 1 (see [13]).
Equation (1) is said to be asymptotically stable, if for any continuous differential initial function and for any delay τ
q
>
0, q
=
1,2
,
…
,
m the analytical solution to (1) satisfies lim
⁡
t
→
∞
⁡
u
(
t
)
=
0.We know that the stability of analytical solution can be studied via the characteristic equation, so we give a criterion for the asymptotic stability of (1), which is based on the following lemmas.Lemma 2 (see [14]).
Assume(4)
Sup
⁡
{
Re
(
λ
)
:
p
(
λ
)
=
0
}
<
0
,
where p
(
λ
)
=
det
⁡
{
λ
A
+
B
+
λ
∑
q
=
1
m
C
q
e
-
λ
τ
q
+
∑
q
=
1
m
D
q
e
-
λ
τ
q
+
∑
q
=
1
m
G
q
τ
q
η
(
e
-
λ
τ
q
)
} is the characteristic polynomial of (1). Then, system (1) is asymptotically stable.Whereη
(
z
) is a complex function defined by
(5)
η
(
z
)
=
{
1
-
z
ln
⁡
z
,
z
∈
C
∖
{
0,1
}
0
,
z
=
0
-
1
,
z
=
1
.
And ln
⁡
z
=
ln
⁡
|
z
|
+
i
arg
z (z
=
0,1
;
-
π
<
arg
z
≤
π) is the principal branch of the multivalued complex natural logarithm.Lemma 3 (see [14]).
Functionη
(
z
) is analytic in C
∖
R
0
- and satisfies |
η
(
z
)
|
≤
1 for |
z
|
≤
1, where R
0
-
=
{
x
∈
R
:
x
≤
0
}.Lemma 4.
If the matrix(
A
+
∑
q
=
1
m
C
q
e
-
λ
τ
q
) is invertible for Re
(
λ
)
≥
r, where r
∈
R, then the function
(6)
p
~
(
λ
)
=
det
⁡
{
λ
2
I
d
+
(
A
+
∑
q
=
1
m
C
q
e
-
λ
τ
q
)
-
1
i
i
i
i
i
i
i
i
i
i
i
i
×
(
λ
B
+
λ
∑
q
=
1
m
D
q
e
-
λ
τ
q
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
+
∑
q
=
1
m
G
q
(
1
-
e
-
λ
τ
q
)
{
λ
2
I
d
+
(
A
+
∑
q
=
1
m
C
q
e
-
λ
τ
q
)
-
1
)
}
has at most a finite number of zeros for Re
(
λ
)
≥
r.Proof.
WhenRe
(
λ
)
≥
r, the function p
~
(
λ
) can be expanded into the following form:
(7)
p
~
(
λ
)
=
λ
2
d
+
ψ
2
d
-
1
(
e
-
λ
τ
1
,
e
-
λ
τ
2
,
…
,
e
-
λ
τ
m
1
)
λ
2
d
-
1
+
⋯
+
ψ
0
(
e
-
λ
τ
1
,
e
-
λ
τ
2
,
…
,
e
-
λ
τ
m
1
)
,
where ψ
i
(
e
-
λ
τ
1
,
e
-
λ
τ
2
,
…
,
e
-
λ
τ
m
1
), i
=
0,1
,
…
,
2
d
-
1, are rational functions for the expressions e
-
λ
τ
1
,
e
-
λ
τ
2
,
…
,
e
-
λ
τ
m
1, and they have no poles for Re
(
λ
)
≥
r.
Sinceτ
i
>
0, we have that
(8)
|
e
-
λ
τ
i
|
=
e
-
τ
i
Re
(
λ
)
≤
e
-
τ
i
r
,
for
Re
(
λ
)
≥
r
.
Hence, there exist constantsK
i
>
0 such that
(9)
|
ψ
i
(
e
-
λ
τ
1
,
e
-
λ
τ
2
,
…
,
e
-
λ
τ
m
1
)
|
≤
K
i
,
i
=
0,1
,
…
,
2
d
-
1
.
Let M be a positive number large enough such that
(10)
K
2
d
-
1
M
+
K
2
d
-
2
M
2
+
⋯
+
K
0
M
2
d
<
1
,
which implies that, forRe
(
λ
)
≥
r and |
λ
|
≥
M,
(11)
|
p
~
(
λ
)
|
≥
|
λ
|
2
d
[
1
-
K
2
d
-
1
M
-
K
2
d
-
2
M
2
-
⋯
-
K
0
M
2
d
]
>
0
.
That is, p
~
(
λ
)
≠
0 in the set {
λ
:
Re
λ
≥
r
,
|
λ
|
≥
M
}.
By the isolation property of the zeros for analytic functions,p
~
(
λ
) has at most a finite number of zeros in the set {
λ
:
Re
λ
≥
r
,
|
λ
|
<
M
}; this proves the lemma.In the following, we denote the spectrum of a square matrixA by σ
(
A
) and introduce the set
(12)
C
-
=
{
z
∈
C
:
Re
(
z
)
<
0
}
.Theorem 5.
System (1) is asymptotically stable if the following conditions are satisfied:(a)
det
⁡
(
A
+
∑
q
=
1
m
ξ
q
C
q
)
≠
0 for |
ξ
q
|
≤
1,
(b)
σ
(
G
(
ξ
)
)
⊆
C
- for ξ
=
(
ξ
1
,
ξ
2
,
…
,
ξ
m
)
Τwith |
ξ
q
|
≤
1, where
(13)
G
(
ξ
)
=
(
A
+
∑
q
=
1
m
ξ
q
C
q
)
-
1
×
(
-
B
-
∑
q
=
1
m
ξ
q
D
q
-
∑
q
=
1
m
η
(
ξ
q
)
G
q
τ
q
)
.Proof.
When|
ξ
q
|
≤
1, q
=
1,2
,
…
,
m, condition (a) leads to
(14)
P
^
(
λ
,
ξ
1
,
ξ
2
,
…
,
ξ
m
)
=
det
⁡
(
λ
A
+
λ
∑
q
=
1
m
ξ
q
C
q
+
B
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
+
∑
q
=
1
m
ξ
q
D
q
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
+
∑
q
=
1
m
η
(
ξ
q
)
G
q
τ
q
)
=
det
⁡
(
A
+
∑
q
=
1
m
ξ
q
C
q
)
det
⁡
(
λ
I
d
-
G
(
ξ
)
)
.
Condition (b) leads to(15)
P
(
λ
)
=
P
^
(
λ
,
e
-
λ
τ
1
,
e
-
λ
τ
2
,
…
,
e
-
λ
τ
m
)
≠
0
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
for
Re
(
λ
)
≥
0
.
Hence(16)
Sup
⁡
{
Re
(
λ
)
:
P
(
λ
)
=
0
}
≤
0
.
Now we will show that the strict inequality in (16) holds.
Define(17)
F
(
ξ
1
,
ξ
2
,
…
,
ξ
m
)
=
det
⁡
(
A
+
∑
q
=
1
m
ξ
q
C
q
)
,
and then F
(
ξ
1
,
ξ
2
,
…
,
ξ
m
) is a multivariate polynomial and is nonzero on the compact domain defined by |
ξ
q
|
≤
1, q
=
1,2
,
…
,
m, and equal to 1 at the origin. Hence, its modulus is bounded; that is,
(18)
|
F
(
ξ
1
,
ξ
2
,
…
,
ξ
m
)
|
≥
ε
>
0
,
when
|
ξ
q
|
≤
1
,
for
q
=
1,2
,
…
,
m
.
By the continuity ofF, there exists a δ
>
0 such that
(19)
|
F
(
ξ
1
,
ξ
2
,
…
,
ξ
m
)
|
>
0
,
when
|
ξ
q
|
≤
1
+
δ
,
for
q
=
1,2
,
…
,
m
.
It follows from this that(20)
det
⁡
(
A
+
∑
q
=
1
m
ξ
q
C
q
)
≠
0
,
when
|
e
-
λ
τ
q
|
≤
1
+
δ
,
for
q
=
1,2
,
…
,
m
.
Letr be the strictly positive number r
=
ln
⁡
(
1
+
δ
)
/
τ; then
(21)
det
⁡
(
A
+
∑
q
=
1
m
e
-
λ
τ
q
C
q
)
≠
0
for
Re
(
λ
)
≥
-
r
.
Thus, the equationP
~
(
λ
)
=
0 has only a finite number of roots when Re
(
λ
)
≥
-
r, and it holds true for the equation p
(
λ
)
=
0 by condition (a). Combined with (16) we get that the characteristic equation has at most a finite number of roots in the region {
λ
:
-
r
≤
Re
(
λ
)
<
0
}.
Let(22)
-
γ
=
max
⁡
-
r
≤
Re
(
λ
)
<
0
{
Re
(
λ
)
}
;
then γ
>
0.
WhenRe
(
λ
)
>
-
γ, the characteristic equation p
(
λ
)
=
0 has no root. Hence, a strict inequality holds in (16). By Lemma 2, the proof is completed.


### 2.2. The Two-Step Runge-Kutta Methods and the Stability Region

Consider the two-step Runge-Kutta method:(23)
Y
(
n
)
=
h
C
11
F
(
t
n
,
Y
(
n
)
)
+
C
12
y
(
n
-
1
)
y
(
n
)
=
h
C
21
F
(
t
n
,
Y
(
n
)
)
+
C
22
y
(
n
-
1
)
for solving the initial value problem (1).In order to simplify the analysis, we consider two-step Runge-Kutta method (TSRK) of the form(24a)
u
n
+
1
=
(
1
-
θ
)
u
n
+
θ
u
n
-
1
+
h
∑
j
=
1
s
b
^
j
f
(
t
j
,
U
n
j
)
+
h
∑
j
=
1
s
b
~
j
f
(
t
~
j
,
U
n
-
1
j
)
,
(24b)
U
n
i
=
u
n
+
h
∑
j
=
1
s
a
~
i
,
j
f
(
t
j
,
U
n
j
)
i
=
1,2
,
…
,
s
,where t
j
=
t
n
+
c
j
h, t
~
j
=
t
n
-
1
+
c
j
h, u
i is an approximation to u
(
t
i
), h is a fixed step-size, θ, b
^
j, b
~
j, a
~
i
,
j, and c
j are coefficients of the method, 0
≤
θ
≤
1.These methods are a subclass of general linear methods introduced by Butcher [15] and could be possibly also referred to as two-step hybrid methods. They generalize k-step collocation methods (with k
=
2) for ordinary differential equations (ODEs) studied by Lie and Nørsett [16] and Lie [17] and two-step Runge-Kutta methods for ODEs investigated by Byrne and Lambert [18]. The variable stepsize continuous two-step Runge-Kutta methods for ODEs were investigated by Jackiewicz and Tracogna [19]. Here we will represent (24a) and (24b) by the following table of the coefficients:
(25)
C
A
b
^
T
s
θ
b
~
T
=
c
1
a
~
11
a
~
12
⋯
a
~
1
s
c
2
a
~
21
a
~
22
⋯
a
~
2
s
⋮
⋮
⋮
⋱
⋮
c
s
a
~
s
1
a
~
s
2
⋯
a
~
s
s
b
^
1
b
^
2
⋯
b
^
s
s
θ
b
~
1
b
~
2
⋯
b
~
s
,
where c
i
=
∑
j
=
1
s
a
~
i
j and ∑
j
=
1
s
(
b
^
j
+
b
~
j
)
=
1
+
θ.Apply (24a) and (24b) to the basic test equation
(26)
u
′
(
t
)
=
a
u
(
t
)
t
≥
0
,
Re
a
<
0
,
which gives the following equations:
(27)
u
n
+
1
=
(
1
-
θ
)
u
n
+
θ
u
n
-
1
+
α
(
∑
j
=
1
s
b
^
j
U
n
j
+
b
~
j
U
n
-
1
j
)
U
n
i
=
u
n
+
α
∑
j
=
1
s
a
~
i
,
j
U
n
j
.Rewriting (27) we obtain
(28)
u
i
+
1
=
R
(
α
,
θ
)
u
i
+
S
(
α
,
θ
)
u
i
-
1
,
where
(29)
R
(
α
,
θ
)
=
1
-
θ
+
α
b
^
Τ
(
I
-
α
A
~
)
-
1
e
,
S
(
α
,
θ
)
=
θ
+
α
b
~
Τ
(
I
-
α
A
~
)
-
1
e
,
α
=
a
h
,
e
=
[
1,1
,
…
,
1
]
s
Τ
.
To investigate the stability properties of (24a) and (24b) with (26), we must investigate the asymptotic behaviors of the solution to (28). This is determined by the location of roots of the characteristic polynomial
(30)
φ
(
z
)
=
z
2
-
R
(
a
,
θ
)
z
-
S
(
a
,
θ
)
.The stability region of the two-step Runge-Kutta methods (24a) and (24b) is the set of all points α for which the roots of φ
(
z
) are inside or on the unit circle with those on the unit circle being simple. If φ
(
z
) is a Schur polynomial for any α with Re
α
<
0, the stability of the two-step Runge-Kutta method contains the negative half plane; the method is said to be A-stable for ODEs.


## 2.1. Asymptotic Stability of the Analytical Solution of Neutral Delay Integro Differential-Algebraic Equation with Many Delays

In this section, we consider the following linear system:(1)
A
u
′
(
t
)
+
B
u
(
t
)
+
∑
q
=
1
m
C
q
u
′
(
t
-
τ
q
)
+
∑
q
=
1
m
D
q
u
(
t
-
τ
q
)
+
∑
q
=
1
m
G
q
∫
t
-
τ
q
t
u
(
δ
)
d
δ
=
0
,
t
≥
0
,
u
(
t
)
=
φ
(
t
)
,
t
∈
[
-
τ
,
0
)
,
where A, B, C
q, D
q, G
q
∈
R
d
×
d, A is a singular matrix, τ
q is a given positive delay constant (q
=
1,2
,
…
,
m), and 0
<
τ
1
≤
τ
2
≤
⋯
≤
τ
m
=
τ. φ
(
t
) denotes a given vector-valued function and u
(
t
) is a vector-valued unknown function to be solved for t
≥
0.In order to obtain the characteristic equation of system (1), we focus on the exponential solutions u
(
t
)
=
e
s
t
x of (1); here x
=
(
x
1
,
x
2
,
…
,
x
d
)
Τ
∈
C
d denotes the unknown vector. Then we have (2a)
u
(
t
)
=
(
e
s
t
x
1
,
e
s
t
x
2
,
…
,
e
s
t
x
d
)
Τ
,
(2b)
u
′
(
t
)
=
s
e
s
t
x
,
(2c)
u
(
t
-
τ
q
)
=
(
e
s
(
t
-
τ
q
)
x
1
,
e
s
(
t
-
τ
q
)
x
2
,
…
,
e
s
(
t
-
τ
q
)
x
d
)
Τ
,
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
q
=
1,2
,
…
,
m
,
(2d)
u
′
(
t
-
τ
q
)
=
s
e
s
(
t
-
τ
q
)
x
,
q
=
1,2
,
…
,
m
.Substituting the above results into (1), we have the following equation:
(2e)
[
s
A
+
B
+
s
∑
q
=
1
m
C
q
e
-
s
τ
q
+
∑
q
=
1
m
D
q
e
-
s
τ
q
i
i
i
+
s
-
1
∑
q
=
1
m
G
q
(
1
-
e
-
s
τ
q
)
]
x
=
0
.The existence of a nonzerox in (2e) implies the characteristic equation of system (1) holds; that is, the following equation holds:
(3)
det
⁡
[
s
A
+
B
+
s
∑
q
=
1
m
C
q
e
-
s
τ
q
i
i
i
i
i
i
i
i
+
∑
q
=
1
m
D
q
e
-
s
τ
q
i
i
i
i
i
i
i
i
+
s
-
1
∑
q
=
1
m
G
q
(
1
-
e
-
s
τ
q
)
]
=
0
.Definition 1 (see [13]).
Equation (1) is said to be asymptotically stable, if for any continuous differential initial function and for any delay τ
q
>
0, q
=
1,2
,
…
,
m the analytical solution to (1) satisfies lim
⁡
t
→
∞
⁡
u
(
t
)
=
0.We know that the stability of analytical solution can be studied via the characteristic equation, so we give a criterion for the asymptotic stability of (1), which is based on the following lemmas.Lemma 2 (see [14]).
Assume(4)
Sup
⁡
{
Re
(
λ
)
:
p
(
λ
)
=
0
}
<
0
,
where p
(
λ
)
=
det
⁡
{
λ
A
+
B
+
λ
∑
q
=
1
m
C
q
e
-
λ
τ
q
+
∑
q
=
1
m
D
q
e
-
λ
τ
q
+
∑
q
=
1
m
G
q
τ
q
η
(
e
-
λ
τ
q
)
} is the characteristic polynomial of (1). Then, system (1) is asymptotically stable.Whereη
(
z
) is a complex function defined by
(5)
η
(
z
)
=
{
1
-
z
ln
⁡
z
,
z
∈
C
∖
{
0,1
}
0
,
z
=
0
-
1
,
z
=
1
.
And ln
⁡
z
=
ln
⁡
|
z
|
+
i
arg
z (z
=
0,1
;
-
π
<
arg
z
≤
π) is the principal branch of the multivalued complex natural logarithm.Lemma 3 (see [14]).
Functionη
(
z
) is analytic in C
∖
R
0
- and satisfies |
η
(
z
)
|
≤
1 for |
z
|
≤
1, where R
0
-
=
{
x
∈
R
:
x
≤
0
}.Lemma 4.
If the matrix(
A
+
∑
q
=
1
m
C
q
e
-
λ
τ
q
) is invertible for Re
(
λ
)
≥
r, where r
∈
R, then the function
(6)
p
~
(
λ
)
=
det
⁡
{
λ
2
I
d
+
(
A
+
∑
q
=
1
m
C
q
e
-
λ
τ
q
)
-
1
i
i
i
i
i
i
i
i
i
i
i
i
×
(
λ
B
+
λ
∑
q
=
1
m
D
q
e
-
λ
τ
q
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
+
∑
q
=
1
m
G
q
(
1
-
e
-
λ
τ
q
)
{
λ
2
I
d
+
(
A
+
∑
q
=
1
m
C
q
e
-
λ
τ
q
)
-
1
)
}
has at most a finite number of zeros for Re
(
λ
)
≥
r.Proof.
WhenRe
(
λ
)
≥
r, the function p
~
(
λ
) can be expanded into the following form:
(7)
p
~
(
λ
)
=
λ
2
d
+
ψ
2
d
-
1
(
e
-
λ
τ
1
,
e
-
λ
τ
2
,
…
,
e
-
λ
τ
m
1
)
λ
2
d
-
1
+
⋯
+
ψ
0
(
e
-
λ
τ
1
,
e
-
λ
τ
2
,
…
,
e
-
λ
τ
m
1
)
,
where ψ
i
(
e
-
λ
τ
1
,
e
-
λ
τ
2
,
…
,
e
-
λ
τ
m
1
), i
=
0,1
,
…
,
2
d
-
1, are rational functions for the expressions e
-
λ
τ
1
,
e
-
λ
τ
2
,
…
,
e
-
λ
τ
m
1, and they have no poles for Re
(
λ
)
≥
r.
Sinceτ
i
>
0, we have that
(8)
|
e
-
λ
τ
i
|
=
e
-
τ
i
Re
(
λ
)
≤
e
-
τ
i
r
,
for
Re
(
λ
)
≥
r
.
Hence, there exist constantsK
i
>
0 such that
(9)
|
ψ
i
(
e
-
λ
τ
1
,
e
-
λ
τ
2
,
…
,
e
-
λ
τ
m
1
)
|
≤
K
i
,
i
=
0,1
,
…
,
2
d
-
1
.
Let M be a positive number large enough such that
(10)
K
2
d
-
1
M
+
K
2
d
-
2
M
2
+
⋯
+
K
0
M
2
d
<
1
,
which implies that, forRe
(
λ
)
≥
r and |
λ
|
≥
M,
(11)
|
p
~
(
λ
)
|
≥
|
λ
|
2
d
[
1
-
K
2
d
-
1
M
-
K
2
d
-
2
M
2
-
⋯
-
K
0
M
2
d
]
>
0
.
That is, p
~
(
λ
)
≠
0 in the set {
λ
:
Re
λ
≥
r
,
|
λ
|
≥
M
}.
By the isolation property of the zeros for analytic functions,p
~
(
λ
) has at most a finite number of zeros in the set {
λ
:
Re
λ
≥
r
,
|
λ
|
<
M
}; this proves the lemma.In the following, we denote the spectrum of a square matrixA by σ
(
A
) and introduce the set
(12)
C
-
=
{
z
∈
C
:
Re
(
z
)
<
0
}
.Theorem 5.
System (1) is asymptotically stable if the following conditions are satisfied:(a)
det
⁡
(
A
+
∑
q
=
1
m
ξ
q
C
q
)
≠
0 for |
ξ
q
|
≤
1,
(b)
σ
(
G
(
ξ
)
)
⊆
C
- for ξ
=
(
ξ
1
,
ξ
2
,
…
,
ξ
m
)
Τwith |
ξ
q
|
≤
1, where
(13)
G
(
ξ
)
=
(
A
+
∑
q
=
1
m
ξ
q
C
q
)
-
1
×
(
-
B
-
∑
q
=
1
m
ξ
q
D
q
-
∑
q
=
1
m
η
(
ξ
q
)
G
q
τ
q
)
.Proof.
When|
ξ
q
|
≤
1, q
=
1,2
,
…
,
m, condition (a) leads to
(14)
P
^
(
λ
,
ξ
1
,
ξ
2
,
…
,
ξ
m
)
=
det
⁡
(
λ
A
+
λ
∑
q
=
1
m
ξ
q
C
q
+
B
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
+
∑
q
=
1
m
ξ
q
D
q
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
+
∑
q
=
1
m
η
(
ξ
q
)
G
q
τ
q
)
=
det
⁡
(
A
+
∑
q
=
1
m
ξ
q
C
q
)
det
⁡
(
λ
I
d
-
G
(
ξ
)
)
.
Condition (b) leads to(15)
P
(
λ
)
=
P
^
(
λ
,
e
-
λ
τ
1
,
e
-
λ
τ
2
,
…
,
e
-
λ
τ
m
)
≠
0
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
for
Re
(
λ
)
≥
0
.
Hence(16)
Sup
⁡
{
Re
(
λ
)
:
P
(
λ
)
=
0
}
≤
0
.
Now we will show that the strict inequality in (16) holds.
Define(17)
F
(
ξ
1
,
ξ
2
,
…
,
ξ
m
)
=
det
⁡
(
A
+
∑
q
=
1
m
ξ
q
C
q
)
,
and then F
(
ξ
1
,
ξ
2
,
…
,
ξ
m
) is a multivariate polynomial and is nonzero on the compact domain defined by |
ξ
q
|
≤
1, q
=
1,2
,
…
,
m, and equal to 1 at the origin. Hence, its modulus is bounded; that is,
(18)
|
F
(
ξ
1
,
ξ
2
,
…
,
ξ
m
)
|
≥
ε
>
0
,
when
|
ξ
q
|
≤
1
,
for
q
=
1,2
,
…
,
m
.
By the continuity ofF, there exists a δ
>
0 such that
(19)
|
F
(
ξ
1
,
ξ
2
,
…
,
ξ
m
)
|
>
0
,
when
|
ξ
q
|
≤
1
+
δ
,
for
q
=
1,2
,
…
,
m
.
It follows from this that(20)
det
⁡
(
A
+
∑
q
=
1
m
ξ
q
C
q
)
≠
0
,
when
|
e
-
λ
τ
q
|
≤
1
+
δ
,
for
q
=
1,2
,
…
,
m
.
Letr be the strictly positive number r
=
ln
⁡
(
1
+
δ
)
/
τ; then
(21)
det
⁡
(
A
+
∑
q
=
1
m
e
-
λ
τ
q
C
q
)
≠
0
for
Re
(
λ
)
≥
-
r
.
Thus, the equationP
~
(
λ
)
=
0 has only a finite number of roots when Re
(
λ
)
≥
-
r, and it holds true for the equation p
(
λ
)
=
0 by condition (a). Combined with (16) we get that the characteristic equation has at most a finite number of roots in the region {
λ
:
-
r
≤
Re
(
λ
)
<
0
}.
Let(22)
-
γ
=
max
⁡
-
r
≤
Re
(
λ
)
<
0
{
Re
(
λ
)
}
;
then γ
>
0.
WhenRe
(
λ
)
>
-
γ, the characteristic equation p
(
λ
)
=
0 has no root. Hence, a strict inequality holds in (16). By Lemma 2, the proof is completed.


## 2.2. The Two-Step Runge-Kutta Methods and the Stability Region

Consider the two-step Runge-Kutta method:(23)
Y
(
n
)
=
h
C
11
F
(
t
n
,
Y
(
n
)
)
+
C
12
y
(
n
-
1
)
y
(
n
)
=
h
C
21
F
(
t
n
,
Y
(
n
)
)
+
C
22
y
(
n
-
1
)
for solving the initial value problem (1).In order to simplify the analysis, we consider two-step Runge-Kutta method (TSRK) of the form(24a)
u
n
+
1
=
(
1
-
θ
)
u
n
+
θ
u
n
-
1
+
h
∑
j
=
1
s
b
^
j
f
(
t
j
,
U
n
j
)
+
h
∑
j
=
1
s
b
~
j
f
(
t
~
j
,
U
n
-
1
j
)
,
(24b)
U
n
i
=
u
n
+
h
∑
j
=
1
s
a
~
i
,
j
f
(
t
j
,
U
n
j
)
i
=
1,2
,
…
,
s
,where t
j
=
t
n
+
c
j
h, t
~
j
=
t
n
-
1
+
c
j
h, u
i is an approximation to u
(
t
i
), h is a fixed step-size, θ, b
^
j, b
~
j, a
~
i
,
j, and c
j are coefficients of the method, 0
≤
θ
≤
1.These methods are a subclass of general linear methods introduced by Butcher [15] and could be possibly also referred to as two-step hybrid methods. They generalize k-step collocation methods (with k
=
2) for ordinary differential equations (ODEs) studied by Lie and Nørsett [16] and Lie [17] and two-step Runge-Kutta methods for ODEs investigated by Byrne and Lambert [18]. The variable stepsize continuous two-step Runge-Kutta methods for ODEs were investigated by Jackiewicz and Tracogna [19]. Here we will represent (24a) and (24b) by the following table of the coefficients:
(25)
C
A
b
^
T
s
θ
b
~
T
=
c
1
a
~
11
a
~
12
⋯
a
~
1
s
c
2
a
~
21
a
~
22
⋯
a
~
2
s
⋮
⋮
⋮
⋱
⋮
c
s
a
~
s
1
a
~
s
2
⋯
a
~
s
s
b
^
1
b
^
2
⋯
b
^
s
s
θ
b
~
1
b
~
2
⋯
b
~
s
,
where c
i
=
∑
j
=
1
s
a
~
i
j and ∑
j
=
1
s
(
b
^
j
+
b
~
j
)
=
1
+
θ.Apply (24a) and (24b) to the basic test equation
(26)
u
′
(
t
)
=
a
u
(
t
)
t
≥
0
,
Re
a
<
0
,
which gives the following equations:
(27)
u
n
+
1
=
(
1
-
θ
)
u
n
+
θ
u
n
-
1
+
α
(
∑
j
=
1
s
b
^
j
U
n
j
+
b
~
j
U
n
-
1
j
)
U
n
i
=
u
n
+
α
∑
j
=
1
s
a
~
i
,
j
U
n
j
.Rewriting (27) we obtain
(28)
u
i
+
1
=
R
(
α
,
θ
)
u
i
+
S
(
α
,
θ
)
u
i
-
1
,
where
(29)
R
(
α
,
θ
)
=
1
-
θ
+
α
b
^
Τ
(
I
-
α
A
~
)
-
1
e
,
S
(
α
,
θ
)
=
θ
+
α
b
~
Τ
(
I
-
α
A
~
)
-
1
e
,
α
=
a
h
,
e
=
[
1,1
,
…
,
1
]
s
Τ
.
To investigate the stability properties of (24a) and (24b) with (26), we must investigate the asymptotic behaviors of the solution to (28). This is determined by the location of roots of the characteristic polynomial
(30)
φ
(
z
)
=
z
2
-
R
(
a
,
θ
)
z
-
S
(
a
,
θ
)
.The stability region of the two-step Runge-Kutta methods (24a) and (24b) is the set of all points α for which the roots of φ
(
z
) are inside or on the unit circle with those on the unit circle being simple. If φ
(
z
) is a Schur polynomial for any α with Re
α
<
0, the stability of the two-step Runge-Kutta method contains the negative half plane; the method is said to be A-stable for ODEs.


## 3. Asymptotic Stability of TSRK Methods for Neutral Delay Integro Differential-Algebraic Equation with Many Delays

In this section, we will confine our discussion to neutral delay integro differential-algebraic equation with commensurate delays, that is, systems of the form (1) with τ
q
=
q
τ, τ
=
M
h, M is a positive integer, q
=
1,2
,
…
,
m.Definition 6 (see [20]).
A numerical method for asymptotically stable system (1) is called asymptotically stable if the numerical solution satisfies
(31)
lim
⁡
n
→
∞
u
n
=
0
.Applying the two-step method (24a) and (24b) to (1), we have
(32)
u
n
+
1
=
(
1
-
θ
)
u
n
+
θ
u
n
-
1
+
∑
j
=
1
s
b
^
j
K
n
,
j
+
∑
j
=
1
s
b
~
j
K
n
-
1
,
j
,
(33)
A
K
n
,
i
+
h
B
(
u
n
+
∑
j
=
1
s
a
~
i
j
k
n
,
j
)
+
∑
q
=
1
m
C
q
K
n
-
q
M
,
i
+
h
∑
q
=
1
m
D
q
(
u
n
-
q
M
+
∑
j
=
1
s
a
~
i
j
K
n
-
q
M
,
j
)
+
h
∑
q
=
1
m
∑
γ
=
0
q
M
γ
r
G
q
·
(
u
n
-
r
+
∑
j
=
1
s
a
~
i
j
k
n
-
r
,
j
)
=
0
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
for
i
=
1,2
,
…
,
s
,
where K
n
,
i
=
[
K
n
,
i
1
,
K
n
,
i
2
,
…
,
K
n
,
i
d
]
Τ, i
=
1,2
,
…
,
s, are stage derivatives multiplied by h.Let(34)
b
^
T
=
[
b
^
1
,
b
^
2
,
…
,
b
^
s
]
,
b
~
T
=
[
b
~
1
,
b
~
2
,
…
,
b
~
s
]
,
A
~
=
(
a
~
i
j
)
.We assume that all the eigenvalues ofA
~ have positive real part. Rearrange the variables of the stage derivatives as
(35)
K
n
=
[
K
n
,
1
1
,
K
n
,
2
1
,
…
,
K
n
,
s
1
,
K
n
,
1
2
,
K
n
,
2
2
,
i
i
i
i
i
i
…
,
K
n
,
s
2
,
…
,
K
n
,
1
d
,
K
n
,
2
d
,
…
,
K
n
,
s
d
]
Τ
.
Define
(36)
Y
n
=
(
K
n
Τ
,
u
n
+
1
Τ
)
Τ
,
B
-
=
h
B
,
D
-
q
=
h
D
q
,
G
̿
q
=
h
2
G
q
.
Rewrite (32) and (33) in the form
(37)
[
A
⊗
I
s
+
B
-
⊗
A
~
0
-
I
d
⊗
b
^
Τ
I
d
]
Y
n
+
[
0
B
-
⊗
e
-
I
d
⊗
b
~
Τ
-
(
1
-
θ
)
I
d
]
Y
n
-
1
+
[
0
0
0
-
θ
I
d
]
Y
n
-
2
+
∑
q
=
1
m
[
C
q
⊗
I
s
+
D
-
q
⊗
A
~
0
0
0
]
Y
n
-
q
M
+
∑
q
=
1
m
[
0
D
-
q
⊗
e
0
0
]
Y
n
-
q
M
-
1
+
∑
q
=
1
m
∑
r
=
0
M
q
[
C
q
⊗
I
s
+
D
-
q
⊗
A
~
0
0
0
]
Y
n
-
r
+
∑
q
=
1
m
∑
r
=
0
M
q
[
0
γ
r
G
̿
q
⊗
e
0
0
]
Y
n
-
r
-
1
=
0
.The characteristic polynomial of (37) is given by
(38)
p
(
z
)
=
det
⁡
[
T
1
(
z
)
T
2
(
z
)
T
3
(
z
)
T
4
(
z
)
]
,
z
∈
C
,
where
(39)
T
1
(
z
)
=
z
2
[
+
∑
q
=
1
m
∑
r
=
0
M
q
γ
r
G
̿
q
⊗
A
~
z
-
r
]
(
A
⊗
I
s
+
B
-
⊗
A
~
)
i
i
i
i
i
i
i
+
∑
q
=
1
m
(
C
q
⊗
I
s
+
D
-
q
⊗
A
~
)
z
-
q
M
i
i
i
i
i
i
+
∑
q
=
1
m
∑
r
=
0
M
q
γ
r
G
̿
q
⊗
A
~
z
-
r
]
,
T
2
(
z
)
=
z
[
B
-
⊗
e
+
∑
q
=
1
m
D
-
q
⊗
e
z
-
q
M
i
i
i
i
i
i
+
∑
q
=
1
m
∑
r
=
0
M
q
γ
r
G
̿
q
⊗
e
r
-
r
]
,
T
3
(
z
)
=
-
z
2
I
d
⊗
b
^
Τ
-
z
I
d
⊗
b
~
Τ
,
T
4
(
z
)
=
z
2
I
d
-
z
(
1
-
θ
)
I
d
-
θ
I
d
.
Following from the theorem on difference equations, we get that if all the zeros z of (38) satisfy |
z
|
<
1, then
(40)
lim
⁡
n
→
∞
Y
n
=
0
.
Hence, we formulate the following lemmas.Lemma 7 (see [21]).
If all the zerosz of (38) satisfy |
z
|
<
1, then numerical method ((32) and (33)) satisfies
(41)
lim
⁡
n
→
∞
u
n
=
0
.Lemma 8.
Assume that condition (a) of Theorem5 holds and assume that [
I
s
-
λ
l
(
r
(
z
)
)
A
~
] are invertible for |
z
|
≥
1, where
(42)
r
(
z
)
=
(
A
+
∑
q
=
1
m
z
-
q
M
C
q
)
-
1
×
(
-
B
-
-
∑
q
=
1
m
z
-
q
M
D
-
q
-
∑
q
=
1
m
∑
r
=
0
q
M
γ
r
z
-
r
G
̿
q
)
;
then, det
⁡
[
T
1
(
z
)
]
≠
0, for |
z
|
≥
1.Proof.
Condition (a) in Theorem5 implies that the matrix (
A
+
∑
q
=
1
m
z
-
q
M
C
q
) is invertible for |
z
|
≥
1; then T
1
(
z
)
=
z
2
[
(
A
+
∑
q
=
1
m
z
-
q
M
C
q
)
⊗
I
s
]
[
I
d
⊗
I
s
-
r
(
z
)
⊗
A
~
].
We have that(43)
det
⁡
[
T
1
(
z
)
]
=
z
2
s
d
[
det
⁡
(
A
+
∑
q
=
1
m
z
-
q
M
C
q
)
]
s
×
det
⁡
[
I
d
⊗
I
s
-
r
(
z
)
⊗
A
~
]
=
z
2
s
d
[
det
⁡
(
A
+
∑
q
=
1
m
z
-
q
M
C
q
)
]
s
×
∏
l
=
1
d
∐
j
=
1
s
[
1
-
λ
l
(
r
(
z
)
)
λ
j
(
A
~
)
]
.
The matrix I
s
-
λ
l
(
r
(
z
)
)
A
~ is invertible meaning that λ
l
(
r
(
z
)
)
λ
j
(
A
~
)
≠
1 for all l
,
j.
Hence,det
⁡
[
T
1
(
z
)
]
≠
0, for |
z
|
≥
1.Theorem 9.
If the system ((32) and (33)) satisfies Lemma 8 and the following conditions,(
a
~
)
det
⁡
(
A
+
∑
q
=
1
m
ξ
q
C
q
)
≠
0, for |
ξ
q
|
≤
1;
(
b
~
)
|
1
-
θ
-
z
-
1
θ
-
(
b
^
Τ
+
z
-
1
b
~
Τ
)
[
I
d
⊗
I
s
-
r
(
z
)
⊗
A
~
]
-
1
[
e
s
⊗
r
(
z
)
]
|
<
1, for |
z
|
≥
1,then the solution of the TSRK methods for (1) is asymptotically stable.Proof.
By Lemma7, we need to prove that all the zeros of (38) satisfy |
z
|
<
1.
If these were not true, there would exist az
0
∈
C with |
z
0
|
≥
1, such that
(44)
det
⁡
[
T
1
(
z
0
)
T
2
(
z
0
)
T
3
(
z
0
)
T
4
(
z
0
)
]
=
0
.
By Lemma 8, we have that det
⁡
[
T
1
(
z
0
)
]
≠
0.
Hence, (44) is equivalent to
(45)
det
⁡
⌊
T
4
(
z
0
)
-
T
3
(
z
0
)
T
1
-
1
(
z
0
)
T
2
(
z
0
)
⌋
=
0
.
Using the Kronecker product [5, chapter 4], we have that
(46)
det
⁡
⌊
T
4
(
z
0
)
-
T
3
(
z
0
)
T
1
-
1
(
z
0
)
T
2
(
z
0
)
⌋
=
det
⁡
{
[
(
A
+
∑
q
=
1
m
z
0
-
q
M
C
q
)
⊗
I
s
]
-
1
z
0
2
I
d
-
z
0
(
1
-
θ
)
I
d
-
θ
I
d
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
-
(
z
0
2
I
d
⊗
b
^
Τ
+
z
0
I
d
⊗
b
~
Τ
)
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
×
z
0
-
2
[
(
A
+
∑
q
=
1
m
z
0
-
q
M
C
q
)
⊗
I
s
]
-
1
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
·
[
I
d
⊗
I
s
-
r
(
z
0
)
⊗
A
~
]
-
1
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
·
z
0
[
B
-
⊗
e
+
∑
q
=
1
m
D
-
q
⊗
e
z
0
-
q
M
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
+
∑
q
=
1
m
∑
r
=
0
M
q
γ
r
G
̿
q
⊗
e
z
0
-
r
]
[
(
A
+
∑
q
=
1
m
z
0
-
q
M
C
q
)
⊗
I
s
]
-
1
}
=
z
0
d
det
⁡
{
z
0
I
d
-
[
(
1
-
θ
)
I
d
-
z
0
-
1
θ
I
d
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
-
(
I
d
⊗
b
^
Τ
+
z
0
-
1
I
d
⊗
b
~
Τ
)
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
×
[
I
d
⊗
I
s
-
r
(
z
0
)
⊗
A
~
]
-
1
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
×
[
e
⊗
r
(
z
0
)
]
(
1
-
θ
)
I
d
-
z
0
-
1
θ
I
d
]
}
.
Combining (45) and (46) gives that
(47)
det
⁡
{
z
0
I
d
-
[
[
I
d
⊗
I
s
-
r
(
z
0
)
⊗
A
~
]
-
1
(
1
-
θ
)
I
d
-
z
0
-
1
θ
I
d
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
-
(
I
d
⊗
b
^
Τ
+
z
0
-
1
I
d
⊗
b
~
Τ
)
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
i
·
[
I
d
⊗
I
s
-
r
(
z
0
)
⊗
A
~
]
-
1
[
e
⊗
r
(
z
0
)
]
]
}
=
0
,
which implies |
1
-
θ
-
z
0
-
1
θ
-
(
b
^
Τ
+
z
0
-
1
b
~
Τ
)
[
I
d
⊗
I
s
-
r
(
z
0
)
⊗
A
~
]
-
1
[
e
⊗
r
(
z
0
)
]
|
≥
1.
This contradicts the assumption that(
b
~
) for |
z
0
|
≥
1. Hence, the theorem is proved.


## 4. Numerical Experiments

Example 1.
Consider the following linear system:(48)
A
u
′
(
t
)
+
B
u
(
t
)
+
∑
q
=
1
m
C
q
u
′
(
t
-
τ
q
)
+
∑
q
=
1
m
D
q
u
(
t
-
τ
q
)
+
∑
q
=
1
m
G
q
∫
t
-
τ
q
t
u
(
δ
)
d
δ
=
0
,
t
≥
0
,
u
(
t
)
=
φ
(
t
)
,
t
∈
[
-
τ
,
0
)
,
with m
=
2, τ
1
=
1, τ
2
=
2, and φ
(
t
)
=
(
cos
⁡
(
t
)
,
sin
(
t
)
,
cos
⁡
(
t
)
)
Τ for t
∈
[
-
2,0
],
where(49)
A
=
[
1
0
0
0
0
0
0
1
0
]
,
B
=
[
1.5
0
0
0
0
-
0.4
0
0.8
0.5
]
,
D
1
=
[
0.2
0.4
0
0
0
0.1
0
0.2
-
0.2
]
,
D
2
=
[
0.1
0
0
0
0
0.1
0
0.15
0.8
]
,
G
1
=
[
0.2
0
0
0
0
-
0.2
0
0.1
0.25
]
,
G
2
=
[
0.1
0
0
0
0
-
0.2
0
0.05
0.25
]
,
and C
1
=
-
0.3
A and C
2
=
-
0.5
A.
Here the matrix coefficients satisfy Theorem9. Hence, the system is asymptotically stable.
We choose the A-stable TSRK methods as follows [22]:
(50)
C
A
b
^
T
s
θ
b
~
T
=
51
32
5151
9760
2601
2440
103
256
-
10609
156160
s
73439
156160
636886846889
1074516737280
s
61448158637
134314592160
16977449
36697976
21872982199
1074516737280
52658918227
134314592160
,
where
(51)
c
i
=
∑
j
=
1
s
a
~
i
j
,
∑
j
=
1
s
(
b
^
j
+
b
~
j
)
=
1
+
θ
.
It can be easily seen that the A-stable TSRK method is asymptotically stable, which illuminates the conclusion of Theorem9.


## 5. Conclusions

This paper develops the asymptotic stability of the two-step Runge-Kutta methods for neutral delay integro differential-algebraic equations with many delays. It studies the asymptotic stability of the analytical solution and introduces two step Runge-Kutta methods and the stability region. It also deals with the asymptotic stability of two-step Runge-Kutta method for neutral delay integro differential-algebraic equations with many delays and proves that the A-stable two-step Runge-Kutta methods are asymptotically stable for neutral delay integro differential-algebraic equations with many delays.



---

*Source: 918371-2014-06-09.xml*