# Stability Analysis of One-Leg Methods for Nonlinear Neutral Delay Integrodifferential Equations

Authors: Yuexin Yu; Liping Wen
Journal: Discrete Dynamics in Nature and Society
 (2015)

---

## Abstract

This paper is concerned with the numerical solution of nonlinear neutral delay integrodifferential equations (NDIDEs). The adaptation of one-leg methods is considered. It is proved that anA-stable one-leg method can preserve the global stability and a strongly A-stable one-leg method can preserve the asymptotic stability of the analytical solution of nonlinear NDIDEs. Numerical tests are given to confirm the theoretical results.


---

## Body

## 1. Introduction

In this paper, we consider the initial value problem (IVP) of nonlinear neutral delay integrodifferential equations:(1)
y
′
t
=
f
t
,
y
t
,
y
t
-
τ
,
y
′
t
-
τ
,
∫
t
-
τ
t
g
t
,
ξ
,
y
ξ
d
ξ
,
t
≥
0
,
y
t
=
φ
t
,
-
τ
≤
t
≤
0
,where τ
>
0 is a constant delay, φ
:
[
-
τ
,
0
]
→
C
N is a given continuously differential function, and f
:
[
0
,
+
∞
)
×
C
N
×
C
N
×
C
N
×
C
N
→
C
N and g
:
[
0
,
+
∞
)
×
R
×
C
N
→
C
N are given continuous mappings.Neutral delay integrodifferential equations (NDIDEs) arise widely in scientific and engineering fields such as physics, biology, medicine, economics, and control system (see [1–3] and the references therein). Generally speaking, it is difficult to obtain the analytical solutions of such equations. In view of this, people began to study the numerical solutions of the equations. For the special cases of NDIDEs, such as delay differential equations, delay integrodifferential equations, and neutral delay differential equations, the theory of computational methods has been studied by many authors and a great deal of interesting results have been found in recent 30 years. But, for NDIDEs, only a few results have been presented in the literature. In 2005, Zhao et al. [4] discussed the asymptotic stability of analytical solution and numerical solution (obtained by linear θ-methods and BDF methods) of linear neutral Volterra delay integrodifferential system: (2)
A
u
′
t
+
B
u
t
+
C
u
′
t
-
τ
+
D
u
t
-
τ
+
G
∫
t
-
τ
t
u
x
d
x
=
0
,
t
>
0
,
u
t
=
ϕ
t
,
-
τ
≤
t
≤
0
,where A, B, C, D, G
∈
R
d
×
d, τ
>
0 and the matrix A may be singular. Later, Xu and Zhao [5] further considered the asymptotic stability of Runge-Kutta methods for system (2). In 2008, Zhang and Vandewalle [6] dealt with the asymptotic stability of exact and discrete solutions of neutral multidelay integrodifferential equations. Sufficient conditions for the asymptotic stability of the analytical solution have been derived, and the asymptotic stability criteria of Runge-Kutta methods and linear multistep methods were constructed. Wu and Gan [7] investigated a test equation for one-dimension linear NDIDEs and got some delay-dependent stability results.For the nonlinear NDIDEs (1), Yu and Li [8] discussed the global stability and asymptotic stability of (
k
,
l
) algebraically stable Runge-Kutta methods. Recently, Hu and Huang [9] considered the analytical and numerical stability of nonlinear NDIDEs. Sufficient conditions for the analytical stability of nonlinear NDIDEs are derived, and they proved that any A-stable linear multistep method can preserve the asymptotic stability of the analytical solution of nonlinear NDIDEs (1). For another case of nonlinear NDIDEs, namely, the NDIDEs of “Hale’s form” (3)
y
t
-
N
y
t
-
τ
′
=
f
t
,
y
t
,
y
t
-
τ
,
∫
t
-
τ
t
g
t
,
ξ
,
y
ξ
d
ξ
,
t
≥
0
,
y
t
=
φ
t
,
-
τ
≤
t
≤
0
,Yu et al. [10] and Zhang et al. [11, 12] investigated the stability of Runge-Kutta methods and one-leg methods, respectively.In this paper we are interested in the stability of one-leg methods for nonlinear NDIDEs (1). It is proved that an A-stable one-leg method can preserve the global stability and a strongly A-stable one-leg method can preserve the asymptotic stability of the analytical solution of nonlinear NDIDEs. Numerical tests are given to confirm the theoretical results in the end.


## 2. Problem ClassD
(
α
,
β
1
,
β
2
,
β
3
,
γ
,
δ
) and Its Stability

Let〈
·
,
·
〉 denote the inner product and · the corresponding norm in space C
N. Assume that the continuous mappings f and g in problem (1) satisfy the following conditions: (4)
R
e
u
1
-
u
2
,
f
t
,
u
1
,
v
,
w
,
x
-
f
t
,
u
2
,
v
,
w
,
x
≤
α
u
1
-
u
2
2
,
∀
t
≥
0
,
u
1
,
u
2
,
v
,
w
,
x
∈
C
N
,
(5)
f
t
,
u
,
v
1
,
w
1
,
x
1
-
f
t
,
u
,
v
2
,
w
2
,
x
2
≤
β
1
v
1
-
v
2
+
β
2
w
1
-
w
2
+
β
3
x
1
-
x
2
,
∀
t
≥
0
,
u
,
v
1
,
v
2
,
w
1
,
w
2
,
x
1
,
x
2
∈
C
N
,
(6)
g
t
,
θ
,
s
1
-
g
t
,
θ
,
s
2
≤
γ
s
1
-
s
2
,
∀
t
≥
0
,
t
-
τ
≤
θ
≤
t
,
s
1
,
s
2
∈
C
N
,where α, β
1, β
2, β
3, and γ are real constants and β
2
<
1. Furthermore, we also consider the function F
(
t
,
u
,
v
,
w
,
x
,
r
,
s
)
:
=
f
(
t
,
u
,
v
,
f
(
t
-
τ
,
v
,
w
,
x
,
r
)
,
s
) and assume that it is continuous and satisfies (7)
F
t
,
u
,
v
1
,
w
,
x
,
r
,
s
-
F
t
,
u
,
v
2
,
w
,
x
,
r
,
s
≤
δ
v
1
-
v
2
,
∀
t
≥
τ
,
u
,
v
1
,
v
2
,
w
,
x
,
r
,
s
∈
C
N
,where δ
>
0 is a real constant.Throughout this paper, we assume that problem (1) has a unique exact solution y
(
t
), and we use the symbol D
(
α
,
β
1
,
β
2
,
β
3
,
γ
,
δ
) to denote the problem class consisting of all of problem (1) satisfying conditions (4)-(7).Remark 1.
When the right-hand side function of problem (1) does not possess the term y
′
(
t
-
τ
), problem (1) degenerates into an IVP of delay integrodifferential equations (DIDEs): (8)
y
′
t
=
f
t
,
y
t
,
y
t
-
τ
,
∫
t
-
τ
t
g
t
,
ξ
,
y
ξ
d
ξ
,
t
≥
0
,
y
t
=
φ
t
,
-
τ
≤
t
≤
0
.The stability of numerical methods for DIDEs has been investigated in [13–17].Remark 2.
When the right-hand side function of problem (1) does not possess the integral term, problem (1) degenerates into an IVP of neutral delay differential equations (NDDEs): (9)
y
′
t
=
f
t
,
y
t
,
y
t
-
τ
,
y
′
t
-
τ
,
t
≥
0
,
y
t
=
φ
t
,
-
τ
≤
t
≤
0
.The stability of numerical methods for NDDEs has been researched in [18–22].Remark 3.
When the right-hand side function of problem (1) does not possess the term y
′
(
t
-
τ
) and the integral term, problem (1) degenerates into an IVP of delay differential equations (DDEs): (10)
y
′
t
=
f
t
,
y
t
,
y
t
-
τ
,
t
≥
0
,
y
t
=
φ
t
,
-
τ
≤
t
≤
0
.The stability results of numerical methods for DDEs can be found in [23–28] and so forth.For problems of the classD
(
α
,
β
1
,
β
2
,
β
3
,
γ
,
δ
), Hu and Huang derived the following stability results (see [9]).Theorem 4.
Suppose problem (1) belongs to the class D
(
α
,
β
1
,
β
2
,
β
3
,
γ
,
δ
) satisfying α
+
(
δ
+
β
3
γ
τ
)
/
(
1
-
β
2
)
≤
0. Then one has (11)
y
t
-
z
t
≤
max
⁡
φ
0
-
ψ
0
,
κ
,
∀
t
>
0
,where(12)
κ
=
sup
-
τ
≤
t
≤
0
⁡
β
1
+
β
3
γ
τ
φ
t
-
ψ
t
+
β
2
φ
′
t
-
ψ
′
t
-
αand z
(
t
) denotes the solution of any given perturbed problem of (1): (13)
z
′
t
=
f
t
,
z
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 continuously differential function.Theorem 5.
Suppose problem (1) belongs to the class D
(
α
,
β
1
,
β
2
,
β
3
,
γ
,
δ
) satisfying α
+
(
δ
+
β
3
γ
τ
)
/
(
1
-
β
2
)
<
0. Then one has (14)
lim
t
→
+
∞
⁡
y
t
-
z
t
=
0
.Inequality (11) characterizes the stability property and relation (14) characterizes the asymptotic stability property of problem (1), respectively.


## 3. Stability Analysis of One-Leg Methods for NDIDEs

Consider using a one-legk-step method (for ordinary differential equations) (15)
ρ
E
y
n
=
h
f
σ
E
t
n
,
σ
E
y
n
,
n
=
0,1
,
2
,
…
,to solve problem (1); we have (16)
ρ
E
y
n
=
h
f
σ
E
t
n
,
σ
E
y
n
,
σ
E
y
n
-
m
,
y
~
n
,
G
n
,
n
=
0,1
,
2
,
…
,where h
=
τ
/
m, m is an arbitrarily given positive integer, t
n
=
n
h, E is the translation operator, E
y
n
=
y
n
+
1, y
n is an approximation to y
(
t
n
), y
n
=
φ
(
t
n
) for n
≤
0, and G
n is an approximation to ∫
σ
(
E
)
t
n
-
τ
σ
(
E
)
t
n
g
(
σ
(
E
)
t
n
,
ξ
,
y
(
ξ
)
)
d
ξ, which can be computed by the repeated trapezoidal rule: (17)
G
n
=
h
∑
k
=
1
m
-
1
g
σ
E
t
n
,
σ
E
t
n
-
k
,
σ
E
y
n
-
k
+
1
2
g
σ
E
t
n
,
σ
E
t
n
,
σ
E
y
n
+
1
2
g
σ
E
t
n
,
σ
E
t
n
-
m
,
σ
E
y
n
-
m
.
y
~
n is an approximation to y
′
(
σ
(
E
)
t
n
-
τ
), which is obtained by using the following formula: (18)
y
~
n
=
f
σ
E
t
n
-
τ
,
σ
E
y
n
-
m
,
σ
E
y
n
-
2
m
,
y
~
n
-
m
,
G
n
-
m
,where y
~
n
=
φ
′
(
σ
(
E
)
t
n
-
τ
) for σ
(
E
)
t
n
-
τ
<
0 and ρ
(
x
)
=
∑
j
=
0
k
α
j
x
j and σ
(
x
)
=
∑
j
=
0
k
β
j
x
j are generating polynomials which are assumed to have real coefficients and no common divisor. We also assume ρ
(
1
)
=
0, ρ
′
(
1
)
=
σ
(
1
)
=
1, and σ
′
(
1
)
≥
0.Similarly, applying the same method to perturbed problem (13), we have (19)
ρ
E
z
n
=
h
f
σ
E
t
n
,
σ
E
z
n
,
σ
E
z
n
-
m
,
z
~
n
,
H
n
,
n
=
0,1
,
2
,
…
,
H
n
=
h
∑
k
=
1
m
-
1
g
σ
E
t
n
,
σ
E
t
n
-
k
,
σ
E
z
n
-
k
+
1
2
g
σ
E
t
n
,
σ
E
t
n
,
σ
E
z
n
+
1
2
g
σ
E
t
n
,
σ
E
t
n
-
m
,
σ
E
z
n
-
m
,where z
n and H
n are approximations to z
(
t
n
) and ∫
σ
(
E
)
t
n
-
τ
σ
(
E
)
t
n
g
(
σ
(
E
)
t
n
,
ξ
,
z
(
ξ
)
)
d
ξ, respectively, z
n
=
ψ
(
t
n
) for n
≤
0, and z
~
n is an approximation to z
′
(
σ
(
E
)
t
n
-
τ
), which can be computed by (20)
z
~
n
=
f
σ
E
t
n
-
τ
,
σ
E
z
n
-
m
,
σ
E
z
n
-
2
m
,
z
~
n
-
m
,
H
n
-
m
,where z
~
n
=
ψ
′
(
σ
(
E
)
t
n
-
τ
) for σ
(
E
)
t
n
-
τ
<
0.For a real symmetric positivek
×
k matrix G
=
[
g
i
j
], the norm ·
G is defined by (21)
U
G
=
∑
i
,
j
=
1
k
g
i
j
u
i
,
u
j
1
/
2
,
∀
U
=
u
1
T
,
u
2
T
,
…
,
u
k
T
T
∈
C
k
N
.Theorem 6.
Assume that one-leg method (15) is A-stable. Then the numerical solutions {
y
n
} and {
z
n
}, obtained by using corresponding method (16) to problems (1) and (13) which belong to the class D
(
α
,
β
1
,
β
2
,
β
3
,
γ
,
δ
) with α
+
(
δ
+
β
3
γ
τ
)
/
(
1
-
β
2
)
≤
0, respectively, satisfy the global stability inequality (22)
y
n
-
z
n
≤
C
max
⁡
max
0
≤
j
≤
k
-
1
⁡
y
j
-
z
j
,
Γ
,
n
=
1,2
,
3
,
…for all h
>
0, where C depends only on the method, β
1, β
2, β
3, γ, δ, and τ, and (23)
Γ
=
max
-
τ
≤
t
≤
0
⁡
φ
t
-
ψ
t
,
φ
′
t
-
ψ
′
t
.Proof.
Let(24)
w
n
=
y
n
-
z
n
,
W
n
=
w
n
T
,
w
n
+
1
T
,
…
,
w
n
+
k
-
1
T
T
,
q
n
=
n
+
σ
′
1
m
,where · denotes the integer part; then q
n
τ
≤
σ
(
E
)
t
n
<
(
q
n
+
1
)
τ.
SinceA-stability is equivalent to G-stability (cf. [29]), there is a k
×
k real symmetric positive definite matrix G such that, for any real sequence {
a
i
}
i
=
0
k, the following inequality holds: (25)
A
1
T
G
A
1
-
A
0
T
G
A
0
≤
2
σ
E
a
0
ρ
E
a
0
,where A
i
=
(
a
i
,
a
i
+
1
,
…
,
a
i
+
k
-
1
)
T
(
i
=
0,1
). Therefore, we can easily obtain (cf. [29]) (26)
W
n
+
1
G
2
≤
W
n
G
2
+
2
R
e
σ
E
w
n
,
ρ
E
w
n
.Using condition (4), we have (27)
2
R
e
σ
E
w
n
,
ρ
E
w
n
=
2
h
R
e
σ
E
·
w
n
,
f
σ
E
t
n
,
σ
E
y
n
,
σ
E
y
n
-
m
,
y
~
n
,
G
n
-
f
σ
E
t
n
,
σ
E
z
n
,
σ
E
z
n
-
m
,
z
~
n
,
H
n
≤
2
h
α
σ
E
w
n
2
+
2
h
σ
E
w
n
·
f
σ
E
t
n
,
σ
E
z
n
,
σ
E
y
n
-
m
,
y
~
n
,
G
n
-
f
σ
E
t
n
,
σ
E
z
n
,
σ
E
z
n
-
m
,
z
~
n
,
H
n
.When σ
(
E
)
t
n
-
τ
<
0, that is, q
n
=
0, (27) leads to (28)
2
R
e
σ
E
w
n
,
ρ
E
w
n
≤
2
h
α
σ
E
w
n
2
+
2
h
σ
E
w
n
β
1
+
β
2
Γ
+
β
3
G
n
-
H
n
.When σ
(
E
)
t
n
-
τ
≥
0, that is, q
n
≥
1, using conditions (4)–(7), (27) leads to (29)
2
R
e
σ
E
w
n
,
ρ
E
w
n
≤
2
h
α
σ
E
w
n
2
+
2
h
σ
E
·
w
n
f
σ
E
t
n
,
σ
E
z
n
,
σ
E
y
n
-
m
,
f
σ
E
t
n
-
τ
,
σ
E
y
n
-
m
,
σ
E
y
n
-
2
m
,
y
~
n
-
m
,
G
n
-
m
,
G
n
-
f
σ
E
t
n
,
σ
E
z
n
,
σ
E
z
n
-
m
,
f
σ
E
t
n
-
τ
,
σ
E
z
n
-
m
,
σ
E
z
n
-
2
m
,
z
~
n
-
m
,
H
n
-
m
,
H
n
≤
2
h
α
σ
E
w
n
2
+
2
h
σ
E
w
n
δ
σ
E
·
w
n
-
m
+
β
3
G
n
-
H
n
+
β
2
f
σ
E
t
n
-
τ
,
σ
E
·
z
n
-
m
,
σ
E
y
n
-
2
m
,
y
~
n
-
m
,
G
n
-
m
-
f
σ
E
t
n
-
τ
,
σ
E
z
n
-
m
,
σ
E
z
n
-
2
m
,
z
~
n
-
m
,
H
n
-
m
≤
2
h
α
σ
E
·
w
n
2
+
2
h
σ
E
w
n
δ
∑
j
=
1
q
n
β
2
j
-
1
σ
E
w
n
-
j
m
+
β
3
∑
j
=
0
q
n
-
1
β
2
j
G
n
-
j
m
-
H
n
-
j
m
+
β
2
q
n
β
1
+
β
2
Γ
+
β
3
G
n
-
q
n
m
-
H
n
-
q
n
m
=
2
h
α
σ
E
w
n
2
+
2
h
σ
E
w
n
δ
∑
j
=
1
q
n
β
2
j
-
1
σ
E
w
n
-
j
m
+
β
3
∑
j
=
0
q
n
β
2
j
G
n
-
j
m
-
H
n
-
j
m
+
β
2
q
n
β
1
+
β
2
Γ
,where, here and below, we define ∑
k
=
s
t equal to 0 for t
<
s. Combining (28) and (29) yields (30)
2
R
e
σ
E
w
n
,
ρ
E
w
n
≤
2
h
α
σ
E
w
n
2
+
2
h
σ
E
w
n
δ
∑
j
=
1
q
n
β
2
j
-
1
σ
E
w
n
-
j
m
+
β
3
∑
j
=
0
q
n
β
2
j
G
n
-
j
m
-
H
n
-
j
m
+
β
2
q
n
β
1
+
β
2
Γ
.Substituting (30) into (26) and using condition α
+
(
δ
+
β
3
γ
τ
)
/
(
1
-
β
2
)
≤
0, we obtain (31)
W
n
+
1
G
2
≤
W
n
G
2
+
2
h
α
σ
E
w
n
2
+
2
h
σ
E
w
n
·
δ
∑
j
=
1
q
n
β
2
j
-
1
σ
E
w
n
-
j
m
+
β
3
∑
j
=
0
q
n
β
2
j
G
n
-
j
m
-
H
n
-
j
m
+
β
2
q
n
β
1
+
β
2
Γ
≤
W
n
G
2
+
2
h
α
σ
E
w
n
2
+
h
δ
∑
j
=
0
q
n
β
2
j
σ
E
w
n
2
+
h
δ
∑
j
=
1
q
n
β
2
j
-
1
σ
E
w
n
-
j
m
2
+
h
β
2
q
n
β
1
+
β
2
2
δ
Γ
2
+
2
h
2
β
3
γ
σ
E
w
n
∑
j
=
0
q
n
β
2
j
∑
k
=
1
m
-
1
σ
E
w
n
-
j
m
-
k
+
1
2
σ
E
w
n
-
j
m
+
1
2
σ
E
w
n
-
j
+
1
m
≤
W
n
G
2
+
h
2
α
+
δ
1
-
β
2
+
β
3
γ
τ
1
-
β
2
σ
E
w
n
2
+
h
δ
∑
j
=
1
q
n
β
2
j
-
1
σ
E
w
n
-
j
m
2
+
h
β
2
q
n
β
1
+
β
2
2
δ
Γ
2
+
h
2
β
3
γ
∑
j
=
0
q
n
β
2
j
∑
k
=
1
m
-
1
σ
E
w
n
-
j
m
-
k
2
+
1
2
σ
E
w
n
-
j
m
2
+
1
2
σ
E
w
n
-
j
+
1
m
2
≤
W
n
G
2
+
h
α
σ
E
w
n
2
+
h
δ
∑
j
=
1
q
n
β
2
j
-
1
σ
E
w
n
-
j
m
2
+
h
β
2
q
n
β
1
+
β
2
2
δ
Γ
2
+
h
2
β
3
γ
∑
j
=
0
q
n
β
2
j
∑
k
=
1
m
-
1
σ
E
w
n
-
j
m
-
k
2
+
1
2
σ
E
w
n
-
j
m
2
+
1
2
σ
E
w
n
-
j
+
1
m
2
.By induction, (31) gives (32)
W
n
+
1
G
2
≤
W
0
G
2
+
∑
i
=
0
n
h
α
σ
E
w
i
2
+
∑
i
=
0
n
h
δ
∑
j
=
1
i
+
σ
′
1
/
m
β
2
j
-
1
σ
E
w
i
-
j
m
2
+
h
·
β
1
+
β
2
2
δ
Γ
2
∑
i
=
0
n
β
2
i
+
σ
′
1
/
m
+
h
2
β
3
γ
∑
i
=
0
n
∑
j
=
0
i
+
σ
′
1
/
m
β
2
j
∑
k
=
1
m
-
1
σ
E
w
i
-
j
m
-
k
2
+
1
2
σ
E
w
i
-
j
m
2
+
1
2
σ
E
w
i
-
j
+
1
m
2
≤
W
0
G
2
+
∑
i
=
0
n
h
α
σ
E
w
i
2
+
h
δ
1
-
β
2
∑
i
=
0
n
σ
E
·
w
i
-
m
2
+
h
m
δ
1
-
β
2
β
1
+
β
2
2
Γ
2
+
h
2
β
3
γ
m
1
-
β
2
·
∑
i
=
-
m
n
σ
E
w
i
2
≤
W
0
G
2
+
∑
i
=
0
n
h
α
+
δ
1
-
β
2
+
β
3
γ
τ
1
-
β
2
σ
E
w
i
2
+
τ
β
1
+
β
2
2
δ
1
-
β
2
Γ
2
+
δ
+
β
3
γ
τ
1
-
β
2
h
m
max
-
m
≤
i
≤
-
1
⁡
σ
E
w
i
2
≤
W
0
G
2
+
δ
+
β
3
γ
τ
1
-
β
2
τ
m
a
x
-
m
≤
i
≤
-
1
σ
E
w
i
2
+
τ
β
1
+
β
2
2
δ
1
-
β
2
Γ
2
.Let λ
1 and λ
2 denote the maximum and minimum eigenvalues of the matrix G, respectively. Then, we have (33)
λ
2
w
n
+
k
2
≤
λ
1
∑
j
=
0
k
-
1
w
j
2
+
δ
+
β
3
γ
τ
1
-
β
2
τ
max
-
m
≤
i
≤
-
1
⁡
σ
E
w
i
2
+
τ
β
1
+
β
2
2
δ
1
-
β
2
Γ
2
.Hence (34)
w
n
+
k
2
≤
k
λ
1
λ
2
max
0
≤
j
≤
k
-
1
⁡
w
j
2
+
δ
+
β
3
γ
τ
τ
1
-
β
2
λ
2
max
-
m
≤
i
≤
-
1
⁡
σ
E
w
i
2
+
τ
β
1
+
β
2
2
δ
1
-
β
2
λ
2
Γ
2
.Therefore, there is a real constant C depending only on the method, β
1, β
2, β
3, γ, δ, and τ such that inequality (22) holds and this completes the proof of Theorem 6.Remark 7.
It is well known that many one-leg methods, such as implicit Euler method, the second-order BDF formula method, and one-legθ-methods (1
/
2
≤
θ
≤
1), are all A-stable. Therefore, in terms of Theorem 6, the corresponding methods are globally stable for solving the nonlinear NDIDEs of the class D
(
α
,
β
1
,
β
2
,
β
3
,
γ
,
δ
) which satisfies the condition α
+
(
δ
+
β
3
γ
τ
)
/
(
1
-
β
2
)
≤
0.Next, we further discuss the asymptotic stability of the one-leg method. One-leg method (15) is called strongly A-stable if it is A-stable and the modulus of any root of σ
(
x
) is strictly less than 1.Theorem 8.
Assume that one-leg method (15) is strongly A-stable. Then the numerical solutions {
y
n
} and {
z
n
}, obtained by using corresponding method (16) to problems (1) and (13) which belong to the class D
(
α
,
β
1
,
β
2
,
β
3
,
γ
,
δ
) with α
+
(
δ
+
β
3
γ
τ
)
/
(
1
-
β
2
)
<
0, respectively, satisfy (35)
lim
n
→
+
∞
⁡
y
n
-
z
n
=
0for all h
>
0. Relation (35) characterizes the asymptotic stability property of method (16).Proof.
In terms of the proof of Theorem6, we have (36)
W
n
+
1
G
2
≤
W
0
G
2
+
∑
i
=
0
n
h
α
+
δ
+
β
3
γ
τ
1
-
β
2
σ
E
w
i
2
+
τ
β
1
+
β
2
2
δ
1
-
β
2
Γ
2
+
δ
+
β
3
γ
τ
1
-
β
2
τ
max
-
m
≤
i
≤
-
1
⁡
σ
E
w
i
2
.Since α
+
(
δ
+
β
3
γ
τ
)
/
(
1
-
β
2
)
<
0, it is easily obtained from (36) that (37)
lim
n
→
+
∞
⁡
σ
E
y
n
-
z
n
=
0
.By analogy with the proof of Theorem  4.3 in [27], we have (38)
lim
n
→
+
∞
⁡
y
n
-
z
n
=
0
,and this completes the proof of Theorem 8.Remark 9.
It is well known that many one-leg methods, such as implicit Euler method, the second-order BDF formula method, and one-legθ-methods (1
/
2
<
θ
≤
1), are all strongly A-stable. Therefore, in terms of Theorem 8, the corresponding methods are asymptotically stable for solving the nonlinear NDIDEs of the class D
(
α
,
β
1
,
β
2
,
β
3
,
γ
,
δ
) which satisfies the condition α
+
(
δ
+
β
3
γ
τ
)
/
(
1
-
β
2
)
<
0.


## 4. Numerical Experiments

Example 1.
Consider the one-dimensional parabolic problem with neutral type(39)
∂
u
x
,
t
∂
t
=
∂
2
u
x
,
t
∂
x
2
+
0.01
sin
⁡
u
x
,
t
-
1
+
∂
u
x
,
t
-
1
∂
t
+
∫
t
-
1
t
u
x
,
s
d
s
,
0
<
x
<
1
,
t
≥
0
,
u
x
,
t
=
x
-
x
2
e
-
t
,
0
<
x
<
1
,
-
1
≤
t
≤
0
,
u
0
,
t
=
u
1
,
t
=
0
.After application of the numerical method of lines, we obtain the following NDIDEs: (40)
u
i
′
t
=
1
Δ
x
2
u
i
-
1
t
-
2
u
i
t
+
u
i
+
1
t
+
0.01
sin
⁡
u
i
t
-
1
+
u
i
′
t
-
1
+
∫
t
-
1
t
u
i
s
d
s
,
t
≥
0
,
u
0
t
=
u
N
t
≡
0
,
t
≥
0
,
u
i
t
=
i
Δ
x
1
-
i
Δ
x
e
-
t
,
t
≤
0
,
i
=
1,2
,
…
,
N
-
1
,where Δ
x is the spatial step, N is a natural number such that N
Δ
x
=
1, x
i
=
i
Δ
x, i
=
1,2
,
…
,
N
-
1, and u
i
(
t
)
=
u
(
x
i
,
t
). Then, problem (40) belongs to the class D
(
α
,
β
1
,
β
2
,
β
3
,
γ
,
δ
) with (41)
α
=
-
4
N
2
s
i
n
2
π
2
N
,
β
1
=
β
2
=
0.01
,
β
3
=
γ
=
1
,
δ
=
0.01
1
+
4
N
2
,where the inner product is standard inner product. We take Δ
x
=
0.1 for the numerical method of lines; thus the condition α
+
(
δ
+
β
3
γ
τ
)
/
(
1
-
β
2
)
≤
0(<0) is satisfied, which means the analytical solution of problem (40) is stable and asymptotically stable.
We use the 2-step one-leg methods of order 2:(42)
3
2
y
n
+
2
-
2
y
n
+
1
+
1
2
y
n
=
h
f
t
n
+
2
,
5
4
y
n
+
2
-
1
2
y
n
+
1
+
1
4
y
n
,which is A-stable and strongly A-stable, for solving problem (40) and its perturbed problem, where the initial function of the perturbed problem is v
i
(
t
)
=
i
Δ
x
(
1
-
i
Δ
x
)
e
-
t
+
5
sin
⁡
t
,
-
1
≤
t
≤
0
,
i
=
1,2
,
…
,
9
.
As a comparison, we also use the 2-step one-leg method of order 3:(43)
y
n
+
2
-
y
n
+
1
=
h
f
t
n
+
3
2
h
,
5
12
y
n
+
2
+
8
12
y
n
+
1
-
1
12
y
n
,which is not A-stable, for solving problem (40) and its perturbed problem. We denote the numerical solutions of problem (40) and its perturbed problems u
n and v
n, where u
n and v
n are approximations to [
u
1
(
t
n
)
,
u
2
(
t
n
)
,
…
,
u
9
(
t
n
)
]
T and [
v
1
(
t
n
)
,
v
2
(
t
n
)
,
…
,
v
9
(
t
n
)
]
T, respectively. The values u
n
-
v
n obtained by different methods are listed in Figure 1.Figure 1
Valuesu
n
-
v
n obtained by method (42) (a) or by method (43) (b) with h
=
0.1.
(a)
(b)Example 2.
Consider the following neutral delay integrodifferential equations, which are the models for population dynamics (see [1]): (44)
U
′
t
=
b
0
U
t
+
b
1
V
t
+
b
2
U
t
-
τ
+
b
3
V
t
-
τ
+
b
4
V
′
t
-
τ
+
∫
t
-
τ
t
k
1
t
,
s
,
U
s
,
V
s
d
s
,
V
′
t
=
c
1
V
t
+
c
2
U
t
-
τ
+
c
3
V
t
-
τ
+
c
4
V
′
t
-
τ
+
∫
t
-
τ
t
k
2
t
,
s
,
U
s
,
V
s
d
s
,
t
>
0
,
U
t
=
φ
t
,
V
t
=
ψ
t
,
-
τ
≤
t
≤
0
.We choose b
0
=
-
4
,
b
1
=
1
,
b
2
=
1
,
b
3
=
-
0.5
,
b
4
=
0.1
,
c
1
=
-
10
,
c
2
=
0.5
,
c
3
=
-
1
,
c
4
=
0.05, k
1
(
t
,
s
,
U
,
V
)
=
0.5
U
s
sin
⁡
s, k
2
(
t
,
s
,
U
,
V
)
=
0.5
V
s
cos
⁡
s, and τ
=
1. Then, problem (44) belongs to the class D
(
α
,
β
1
,
β
2
,
β
3
,
γ
,
δ
) with (45)
α
≈
-
3.84
,
β
1
=
1.5
,
β
2
≈
0.112
,
β
3
=
0.5
,
γ
=
1
,
δ
≈
2.622
.Thus the condition α
+
(
δ
+
β
3
γ
τ
)
/
(
1
-
β
2
)
≤
0(<0) is satisfied, which means the analytical solution of problem (44) is stable and asymptotically stable. We use the 2-step one-leg method (42) for solving problem (44) and its perturbed problem, where the initial functions of problem (44) and its perturbed problem are φ
(
t
)
=
1.1, ψ
(
t
)
=
0.9 and φ
~
(
t
)
=
1.1
-
sin
⁡
π
/
2
t, ψ
~
(
t
)
=
0.9
-
sin
⁡
(
π
/
2
)
t for -
1
≤
t
≤
0, respectively. We denote the numerical solutions of problem (44) and its perturbed problems W
n and W
~
n, where W
n and W
~
n are approximations to [
U
(
t
n
)
,
V
(
t
n
)
]
T and [
U
~
(
t
n
)
,
V
~
(
t
n
)
]
T, respectively. The values W
n
-
W
~
n obtained by method (42) with different stepsize h are listed in Figure 2.
From Figures1-2, one can see that the values u
n
-
v
n (or W
n
-
W
~
n) obtained by method (42) are bounded and tend to zero. These coincide with Theorems 6 and 8. However, for method (43), the situation is inverse as one can see that the values u
n
-
v
n are divergent as n
→
+
∞.Figure 2
ValuesW
n
-
W
~
n obtained by method (42) with h
=
0.05 (a) and h
=
0.02 (b).
(a)
(b)


## 5. Conclusion

In the present paper, the adaptation of one-leg methods is applied for solving nonlinear NDIDEs. It is proved that anA-stable one-leg method can preserve the global stability and a strongly A-stable one-leg method can preserve the asymptotic stability of the analytical solution of nonlinear NDIDEs. Investigating the stability of other numerical methods, for example, multistep Runge-Kutta methods, for solving nonlinear NDIDEs will be our future work.



---

*Source: 325364-2015-07-14.xml*