Numerical solution to highly nonlinear neutral-type stochastic differential equation

Applied Numerical Mathematics 140 (2019) 48–75

Contents lists available at ScienceDirect

Applied Numerical Mathematics www.elsevier.com/locate/apnum

Numerical solution to highly nonlinear neutral-type stochastic differential equation Shaobo Zhou a , Hai Jin b,∗ a b

School of Mathematics and Statistics, Huazhong University of Science and Technology, Wuhan, Hubei 430074, PR China School of Mathematics and Statistics, Beijing Institute of Technology, Beijing 100081, PR China

a r t i c l e

i n f o

Article history: Received 15 March 2015 Received in revised form 26 December 2018 Accepted 22 January 2019 Available online 28 January 2019 Keywords: Neutral-type stochastic differential equation Backward Euler–Maruyama method Strong convergence Convergence rate Boundedness

a b s t r a c t In the paper, our main aim is to investigate the strong convergence of the implicit numerical approximations for neutral-type stochastic differential equations with superlinearly growing coefficients. After providing mean-square moment boundedness and mean-square exponential stability for the exact solution, we show that a backward Euler– Maruyama approximation converges strongly to the true solution under polynomial growth conditions for sufficiently small step size. Imposing a few additional conditions, we examine the p-th moment uniform boundedness of the exact and approximate solutions by the stopping time technique, and establish the convergence rate of one half, which is the same as the convergence rate of the classical Euler–Maruyama scheme. Finally, several numerical simulations illustrate our main results. © 2019 IMACS. Published by Elsevier B.V. All rights reserved.

1. Introduction Neutral stochastic differential systems provide an excellent mathematical modeling framework for numerous applications in natural science and technology, such as circuit analysis, computer aided design, real time simulation of mechanical systems, power systems, chemical process simulation, population dynamics and automatic control. For instance, they are frequently used for the study of distributed networks containing lossless transmission lines (see [1,20]). So such models have received more and more attention (see [5,13,25]). Most neutral-type stochastic differential equations cannot be solved explicitly, so the numerical methods have received an increasing attention (see [17,19,22,26,28]). Some research efforts have been devoted to the strong convergence of the explicit Euler–Maruyama (EM) methods with the order 1/2 for the linear or semi-linear neutral-type systems (see [14]). Recently, several authors established the convergence in probability of the explicit EM method for the nonlinear neutraltype stochastic systems. For example, Milo˘sevic´ [9] investigated the convergence in probability of the EM approximations for nonlinear neutral stochastic differential equations with variable delay under the Khasminskii-type conditions. Zhou et al. [23] established the convergence in probability of the EM approximations for neutral stochastic functional differential equations under the polynomial growth conditions. Obradovic´ et al. [11] investigated the convergence in probability of the explicit EM method for neutral stochastic systems with unbounded delay and Markovian switching under nonlinear growth conditions. Adding the linear growth condition on the drift coefficient, they also proved the almost sure exponential stability

*

Corresponding author. E-mail addresses: [email protected] (S. Zhou), [email protected] (H. Jin).

https://doi.org/10.1016/j.apnum.2019.01.014 0168-9274/© 2019 IMACS. Published by Elsevier B.V. All rights reserved.

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

49

of the EM method in the case of the bounded delay. Milo˘sevic´ [10] showed that the appropriate EM equilibrium solution is globally almost surely asymptotically exponentially stable under the linear growth condition. The EM methods have simple algebraic structures, cheap computational costs and acceptable convergence rates, but they are not suited to stochastic systems with super-linearly growing coefficients. Hutzenthaler et al. [2] proved that in the case of super-linearly growing coefficients, the EM approximation may not converge in the strong L p -sense nor in the weak sense to the exact solution. So Mao and Szpruch [6,7,12] proposed an implicit EM scheme and demonstrated that its numerical approximations strongly converges to the true solutions to stochastic differential equations with one-sided linear growing and one-sided Lipschitz drift and super-linearly growing diffusion coefficients. Then Zhou et al. [24,21] proved strong convergence of the backward EM approximation for stochastic differential delay equation under one-sided Lipschitz condition and one-sided polynomial growth condition on the drift coefficient, and polynomial growth condition on the diffusion coefficient. Zong et al. [28] established the strong convergence of split-step scheme for stochastic differential delay equations with non-globally Lipschitz condition coefficients. Recently, stability and convergence of the implicit EM scheme have become an active research field for nonlinear neutral systems. Wang and Chen [16] proved that the semi-implicit Euler method is mean-square stable under the global Lipschitz condition and the linear growth condition. Zhou [22] investigated the exponential stability of the backward EM scheme for neutral system with super-linearly growing coefficients. Zong and Wu [27,29] showed that the implicit method can share the mean square exponential stability of the exact solution for neutral stochastic delay differential equations under the Lyapunov conditions by the Chebyshev inequality and the Borel–Cantelli lemma. Yan et al. [18] proved the strong convergence of the split-step theta method for neutral stochastic delay differential equations with convergence rate of one half. Ji and Yuan [3] investigate the convergence of the tamed EM scheme for neutral stochastic differential delay equations under global and local non-Lipschitz conditions, respectively and provide the convergence rate of tamed EM scheme under the global Lipschitz condition. In the paper, our main aim is to investigate the strong convergence of the implicit numerical approximations for neutral-type stochastic differential equations with polynomially growing coefficients. After providing mean-square moment boundedness of the exact and approximate solutions, we show that a backward EM approximation converges strongly to the true solution under polynomial growth conditions for sufficiently small step size, but the convergence rate is still open, owing to no results on the p-th moment uniform boundedness of the exact and approximate solutions. Imposing a few additional conditions on the drift and diffusion coefficients, we establish the p-th moment uniform boundedness of the exact and approximate solutions by the stopping time technique, and develop the convergence rate of one half, which could be the same as the convergence rate of the classical EM scheme. At the same time, we analyze the key role of the additional conditions in establishing the uniform boundedness. Finally, several numerical experiments are provided to illustrate our results. The structure of the paper is as follows: In the next section, we introduce one-sided polynomial growth conditions under which a neutral-type stochastic differential delay system has a unique global solution, and prove the mean-square boundedness of the solution. Section 3 introduces the backward EM and the forward–backward EM scheme to neutral-type stochastic differential equations and shows the boundedness of numerical solutions. Section 4 proves a strong convergence theorem on the basis of a series of lemmas. Section 5 examines the strong convergence rates of the backward EM scheme under some additional conditions, which is the same as the convergence rate of the classical EM scheme one half. 2. Stability Throughout this paper, unless otherwise specified, let |x| be the Euclidean norm in x∈ Rn . If A is a vector or matrix, its transpose is denoted by A T . If A is a matrix, its trace norm is denoted by | A | = trace( A T A ), while its operator norm is denoted by  A  = sup{| Ax| : |x| = 1}. Let (, F , {Ft }t ≥0 , P) be a complete probability space with a filtration {Ft }t ≥0 , satisfying the usual conditions (i.e., it is increasing and right continuous and F0 contains all P-null sets). Let R+ = [0, +∞). Let τ > 0, denote by C ([−τ , 0]; Rn¯ ) the family of continuous functions from [−τ , 0] to Rn¯ with the norm b ϕ  = sup |ϕ (θ)|. Denote by C F ([−τ , 0]; Rn¯ ) the family of all F0 -measurable bounded C ([−τ , 0]; Rn¯ )-valued random −τ ≤θ≤0

0

variables ξ = {ξ(θ) : −τ ≤ θ ≤ 0}. Let a ∨ b = max{a, b} and a ∧ b = min{a, b}. ¯ Consider an n-dimensional neutral-type stochastic differential delay system

d[x(t ) − u (x(t − τ ))] = f (x(t ), x(t − τ ))dt + g (x(t ), x(t − τ ))dw (t ) n¯

(2.1) n¯

n¯

n¯

n¯

n¯

¯ n¯ ×m

on t ≥ 0 with initial data x0 ≡ ξ = {ξ(θ) : −τ ≤ θ ≤ 0} ∈ C F ([−τ , 0]; R ), and f : R × R −→ R , g : R × R −→ R are 0 Borel-measurable and locally Lipschitz continuous. That is, for each integer R ≥ 1, there is a positive constant K R such that b

| f (x1 , y 1 ) − f (x2 , y 2 )|2 ∨ | g (x1 , y 1 ) − g (x2 , y 2 )|2 ≤ K R (|x1 − x2 |2 + | y 1 − y 2 |2 ) for xi , y i ∈ Rn¯ with |xi | ∨ | y i | ≤ R (i = 1, 2). For the purpose of the stability and convergence, assume that f (0, 0) = g (0, 0) = u (0) = 0. The classical condition of the existence-and-uniqueness and stability imposed on the drift and diffusion coefficients is the linear growth condition. In the paper, we replace the linear growth condition by polynomial growth conditions.

50

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

¯ , λˆ , b0 , b¯ , b, α such that (H1) (polynomial growth conditions) There exist positive constants a0 , a¯ 0 , a, a¯ , λ0 , λ 2 ˜x, f (x, y ) ≤ −a0 |x|2 + a¯ 0 | y |2 + a¯ | y |α +2 − a|x|α +2 ,

¯ y |α +1 + λ| ˆ x|α +1 , | f (x, y )| ≤ λ0 |x| + λ| 2 2 α + 2 + b|x|α +2 | g (x, y )| ≤ b0 |x| + b¯ | y | for x, y ∈ Rn , where  x(t ) = x(t ) − u (x(t − τ )), y (t ) = x(t − τ ). (H2) (contractive condition) For any x1 , x2 ∈ Rn¯ , there exists a constant

κ (0 < κ < 1) such that

|u (x1 ) − u (x2 )| ≤ κ |x1 − x2 |. ¯ Then for any initial data x0 , there almost surely exists a unique Theorem 2.1. Let (H1)–(H3) hold with a0 > b0 + a¯ 0 , a > a¯ + b + b. global solution x(t ) to Eq. (2.1) on t ≥ −τ . Moreover, the solution is mean-square bounded in the sense of

E|x(t )|2 ≤ sup E|x(s)|2 ≤ −τ ≤s≤t

(1 − κ + a¯ 0 τ )Eξ 2 + E|˜x(0)|2 + (¯a + b¯ )τ Eξ α +2 . (1 − κ )2

Proof. Applying the standing truncation technique (see Mao [4], Theorem 3.2.2, P 95 ) to Eq. (2.1) for any initial data x0 , there exists a unique maximal local strong solution to Eq. (2.1) on −τ < t < νe , where νe is the explosion time. To show this solution is global, we only need to show that νe = ∞ a.s. Let k0 > 0 be sufficiently large such that k0 > |x0 |. For each integer k ≥ k0 , define the stopping time

νk = inf{t ∈ [0, νe ) : |˜x(t )| ≥ k}, (k ∈ N ),

(2.2)

where throughout this paper, we set inf ∅ = ∞ (as usual, ∅ = the empty set). By the definition of the stopping time νk , it is obvious that νk is an increasing function with k, so νk → ν∞ ≤ νe (k → ∞)a.s. If we can show that ν∞ = ∞ a.s., then νe = ∞ a.s. which implies that x(t ) is global. In other words, we only prove that P(νk ≤ t ) → 0(k → ∞, t > 0). Define V (x) = |x|2 . Since P(νk ≤ t ) V (˜x(νk )) ≤ E V (˜x(t ∧ νk )) and V (˜x(νk )) = |˜x(νk )|2 = k2 → ∞, then we only need to prove that E V (˜x(t ∧ νk )) < +∞. By the Itô formula, we obtain t∧νk

V (˜x(t ∧ νk )) = V (˜x(0)) +

t∧νk

L V (˜x(s), y (s))ds + 0

V x (˜x(s)) g (x(s), y (s))dw (s),

(2.3)

0

where L V (˜x, y ) = 2x˜ T f (x, y ) + | g (x, y )|2 . According to (H1), it is easy to get

L V (˜x, y ) ≤ a¯ 0 (| y |2 − |x|2 ) + (¯a + b¯ )(| y |α +2 − |x|α +2 )

− (a0 − b0 − a¯ 0 )|x|2 − (a − a¯ − b − b¯ )|x|α +2 .

(2.4)

¯ This, together with (2.3) and (2.4), implies Recalling that a0 > b0 + a¯ 0 , a > a¯ + b + b. t∧νk

V (˜x(t ∧ νk )) ≤ V (˜x(0)) + (¯a + b¯ )

(| y (s)|α +2 − |x(s)|α +2 )ds

0 t∧νk

+ a¯ 0

t∧νk 2

2

(| y (s)| − |x(s)| )ds + 0

V x (x(s)) g (x(s), y (s))dw (s).

(2.5)

0

According to the integral property, we have t∧νk

(| y (s)|α +2 − |x(s)|α +2 )ds =

t ∧νk −τ

|x(s)|α +2 ds −

−τ

0

t∧νk

|x(s)|α +2 ds ≤

0

|x(s)|α +2 ds.

−τ

0

Substituting this into (2.5) and taking expectation, we arrive at

0 E V (˜x(t ∧ νk )) ≤ E V (˜x(0)) + a¯ 0 −τ

|x(s)| ds + (¯a + b¯ )E

0

2

|x(s)|α +2 ds,

−τ

which implies that Eq. (2.1) almost surely admits a unique global solution.

(2.6)

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

51

On the other hand, letting k → ∞ in (2.6) then we have E V (˜x(t )) ≤ E˜x(0)2 + a¯ 0 τ Eξ 2 + (¯a + b¯ )τ Eξ α +2 . Applying the inequality (x + y )2 ≤ (1 −  )−1 x2 +  −1 y 2 , x, y > 0, 0 <  < 1, we may obtain

sup E|x(s)|2 ≤ Eξ 2 + sup E|x(s)|2

−τ ≤s≤t

0≤s≤t

≤ Eξ 2 +

1 1 sup E|˜x(s)|2 + sup E|u (x(s − τ ))|2 . 1 −  0≤s≤t  0≤s≤t

Making use of (H2), we may compute

sup E|u (x(s − τ ))|2 ≤ κ 2 sup E|x(s − τ )|2 ≤ κ 2 sup E|x(s)|2 .

0≤s≤t

−τ ≤s≤t

0≤s≤t

Therefore

sup E|x(s)|2 ≤ Eξ 2 +

−τ ≤s≤t

Choosing

1

sup E|˜x(s)|2 +

1 −  0≤s≤t

κ2 sup E|x(s)|2 .  −τ ≤s≤t

 = κ and recalling that κ < 1, we may obtain

E|x(t )|2 ≤ sup E|x(s)|2 ≤ −τ ≤s≤t

(1 − κ + a¯ 0 τ )Eξ 2 + E˜x(0)2 + (¯a + b¯ )τ Eξ α +2 . 2 (1 − κ )2

¯ Then for any initial data x0 , the solution is mean-square exponenTheorem 2.2. Let (H1)–(H2) hold with a0 > b0 + a¯ 0 , a > a¯ + b + b. tially stable in the sense of lim sup

t →∞

1 t

logE|x(t )|2 ≤ −ε ,

where ε = min{ε1 , ε2 , τ1 log κ1 }, ε1 and ε2 are the roots of

a − b − e ετ (¯a + b¯ ) = 0,

a0 − b0 − a¯ eετ − 2ε (1 + κ 2 eετ ) = 0, respectively. Proof. By the Itô formula, we have

eε(t ∧νk ) V (˜x(t ∧ νk )) = V (˜x(0)) +

t∧νk

eε s [L V (˜x(s), y (s)) + ε V (˜x(s))]ds + M (t ∧ νk ),

(2.7)

0

 t ∧ν where M (t ∧ νk ) = 0 k eε s V x (˜x(s)) g (x(s), y (s))dw (s) is a real-valued continuous local martingale with M (0) = 0. Noting 2 that |˜x| = |x − u ( y )|2 ≤ 2|x|2 + 2κ 2 | y |2 , which yields

L V (˜x, y ) + ε V (˜x) ≤ (¯a0 + 2εκ 2 )(| y |2 − e ετ |x|2 ) + (¯a + b¯ )(| y |α +2 − e ετ |x|α +2 )

¯ ετ )|x|α +2 . − (a0 − b0 − a¯ 0 e ετ − 2ε − 2εκ 2 e ετ )|x|2 − (a − a¯ e ετ − b − be ¯ which implies there exists a positive constant Recalling that a0 > b0 + a¯ 0 , a > a¯ + b + b,

(2.8)

ε < ε1 ∧ ε2 such that

a0 − b0 − a¯ 0 eετ − 2ε (1 + κ 2 eετ ) > 0, a − b − eετ (¯a + b¯ ) > 0. This, together with (2.7) and (2.8), yields

eε(t ∧νk ) E V (˜x(t ∧ νk )) ≤ E V (˜x(0)) + (¯a0 + 2εκ 2 )E

t∧νk

eε s (| y (s)|2 − eετ |x(s)|2 )ds 0

+ (¯a + b¯ )E

t∧νk

eε s (| y (s)|α +2 − eετ |x(s)|α +2 )ds.

0

(2.9)

52

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

Making use of the property of the integral, it is easy to estimate t∧νk

eε s (| y (s)|α +2 − eετ |x(s)|α +2 )ds ≤

0

eε(s+τ ) |x(s)|α +2 ds.

−τ

0

Substitution this into (2.9), one may get

eε(t ∧νk ) E V (˜x(t ∧ νk )) ≤ E˜x(0) + (¯a + b¯ )E

0

eε(s+τ ) |x(s)|α +2 ds + (¯a0 + 2κ 2 ε )E

−τ

0

eε(s+τ ) |x(s)|2 .

−τ

Letting k → ∞, we may obtain

eεt E V (˜x(t )) ≤ E˜x(0)2 + (¯a + b¯ )τ eετ Eξ α +2 + (¯a0 + 2εκ 2 )τ eετ Eξ 2 . Applying the inequality (x + y )2 ≤ (1 −  )−1 x2 +  −1 y 2 , x, y > 0, 0 <  < 1, we may obtain

sup eε s E|x(s)|2 ≤ Eξ 2 + sup eε s E|x(s)|2 ≤ Eξ 2 +

−τ ≤s≤t

0≤s≤t

1 1 sup eε s E|˜x(s)|2 + sup eε s E|u (x(s − τ ))|2 . 1 −  0≤s≤t  0≤s≤t

Making use of (H2), we may compute

sup eε s E|u (x(s − τ ))|2 ≤ κ 2 sup eε s E|x(s − τ )|2 ≤ κ 2 eετ

0≤s≤t

0≤s≤t

sup eε s E|x(s)|2 .

−τ ≤s≤t

This implies

sup eε s E|x(s)|2 ≤ Eξ 2 +

−τ ≤s≤t

Letting

1

sup eε s E|˜x(s)|2 +

1 −  0≤s≤t

κ 2 eετ sup eε s E|x(s)|2 .  −τ ≤s≤t

 = κ , and recalling that ε ≤ τ1 log κ1 , so we obtain

eεt E|x(t )|2 ≤ That is, lim sup

t →∞

(1 − κ )Eξ 2 + E˜x(0)2 + (¯a + b¯ )τ eετ Eξ α +2 + (¯a0 + 2εκ 2 )τ eετ Eξ 2 . (1 − κ )(1 − κ eετ )

1 logE|x(t )|2 t

≤ −ε . 2

3. Backward EM scheme In the section, we shall introduce the discrete backward EM scheme and forward backward EM scheme, then prove the boundedness of the discrete backward EM approximation. Let the step size ∈ (0, 1) be a fraction of τ , namely = τ /m for some integer m. The discrete backward EM approximate solution is defined by

X tk+1 = X tk + u ( X tk+1−m ) − u ( X tk−m ) + f ( X tk+1 , X tk+1−m ) + g ( X tk , X tk−m ) w tk

(3.1)

with X tk = ξ(k ) for −m ≤ k ≤ 0, tk = k , w tk = w ((k + 1) ) − w (k ). Moreover, the increments w tk , k = 1, 2, . . . , are independent N (0, )-distributed Gaussian random variables Ftk -measurable at the mesh-point tk . To guarantee that this method is well defined, we impose one-sided Lipschitz condition on the drift coefficient f (x, y ). (H3) (One-sided Lipschitz condition) There exists a positive constant λ > 0 such that for any xi ∈ Rn (i = 1, 2)

x1 − x2 , f (x1 , y ) − f (x2 , y ) ≤ λ|x1 − x2 |2 . Applying a fixed point theorem one can prove that Eq. (3.1) has a unique solution X tk+1 given X tk if λ < 1, then the backward EM scheme (3.1) is well defined (see e.g. [8]). From now on we always assume that < λ−1 . Denote by X˜ tk = X tk − u ( X tk−m ), then the discrete equation (3.1) reduces to

X˜ tk+1 = X˜ tk + f ( X tk+1 , X tk+1−m ) + g ( X tk , X tk−m ) w tk .

(3.2)

Once we compute the discrete values X tk from the backward EM sequence (3.1), then we define the discrete forward backward EM scheme

Xˆ tk+1 = Xˆ tk + u ( Xˆ tk+1−m ) − u ( Xˆ tk−m ) + f ( X tk , X tk−m ) + g ( X tk , X tk−m ) w tk with Xˆ tk = ξ(k ) = X tk for −m ≤ k ≤ 0, tk = k , where Xˆ 0 = X 0 = x0 .

(3.3)

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

53

The forward backward EM scheme is due to Mao (see [7]) who first introduced it for the convergence analysis. Define

η(t ) = tk for t ∈ [tk , tk+1 ), k ≥ 0, the continuous forward backward EM schemes is now defined by X¯ (t ) − u ( X¯ (t − τ )) = X¯ 0 − u ( X¯ (−τ )) +

t

t f ( X η(s) , X η(s−τ ) )ds +

0

g ( X η(s) , X η(s−τ ) )dw (s)

(3.4)

0

with X¯ (t ) = ξ(t ) for −τ ≤ t ≤ 0. Obviously, X η(s−τ ) = X tk−m for s ∈ [tk , tk+1 ), moreover the continuous and discrete forward backward EM scheme coincide at the grid point, i.e. X¯ (tk ) = Xˆ t . For convenience, let R > 0 be arbitrary and define the k

sequence of stopping time

σ R = inf{k ≥ 0 : | Xtk | > R }.

(3.5)

Clearly, when k ∈ [0, σ R ( w )], | X t i ( w )| ≤ R for i = 1, · · ·, k − 1, but we may have | X tk ( w )| > R. So the following lemma is not

˜ trivial. For convenience, denote by X˜¯ (t ) = X¯ (t ) − u ( X¯ (t − τ )) and Xˆ tk = Xˆ tk − u ( Xˆ tk−m ).

¯ 0 + a¯ 0 , a > a¯ + m ¯ (b + b¯ ). Then for any p ≥ 2 and sufficiently large R, there Lemma 3.1. Let (H1), (H2) and (H3) hold with a0 > mb exists a constant C ( p , R ) such that

E[| X tk | p I [0,σ R ] (k)] < C ( p , R ), ∀k ≥ 0, where the indicative function I [0,σ R ] (k) = 1 if k ∈ [0, σ R ], otherwise it is zero. C ( p , R ) is a constant that depends on p and R, and the product of C ( p , R ) and a positive constant is still denoted by C ( p , R ). Proof. According to the discrete sequence (3.2), we have

| X˜ tk+1 |2 = X˜ tk+1 , f ( X tk+1 , X tk+1−m )

+ X˜ tk+1 , X˜ tk + g ( X tk , X tk−m ) w tk . Making use of the inequality x, y ≤ |x|| y | ≤ 12 (|x|2 + | y |2 ), yields

1 1 | X˜ tk+1 |2 ≤ X˜ tk+1 , f ( X tk+1 , X tk+1−m )

+ | X˜ tk+1 |2 + | X˜ tk + g ( X tk , X tk−m ) w tk |2

≤−

a0 2 1

| X t k + 1 |2 +

a¯ 0 2

2

a¯

2

a

| X tk+1−m |2 + | X tk+1−m |α +2 − | X tk+1 |α +2 2

2

1

+ | X˜ tk+1 |2 + | X˜ tk + g ( X tk , X tk−m ) w tk |2 , 2

2

where the second inequality used the condition (H1). Rearranging the inequality, one can get

| X˜ tk+1 |2 ≤ a¯ | X tk+1−m |α +2 + a¯ 0 | X tk+1−m |2 + | X˜ tk + g ( X tk , X tk−m ) w tk |2 ≤ a¯ | X tk+1−m |α +2 + a¯ 0 | X tk+1−m |2 + 2| X˜ tk |2 + 2| g ( X tk , X tk−m ) w tk |2 . Therefore, by (H3), we have

| X tk |2 ≤ 2| X˜ tk |2 + 2|u ( X tk−m )|2 ≤ 2| X˜ tk |2 + 2κ 2 | X tk−m |2 ≤ (2a0 + 2κ 2 )| X t |2 + 4| X˜ t |2 + 2a¯ | X t |α +2 + 4| g ( X t k−m

p

k −1

p

p

k−m

p

p

k −1

, X tk−1−m ) w tk−1 |2 .

p

The inequality (x + y + z + w ) 2 ≤ 4 2 −1 (x 2 + y 2 + z 2 + w 2 ), ∀x, y , z, w > 0, implies that p

p (α +2)

p

| X tk | p ≤ 4 2 −1 [(2a¯ ) 2 | X tk−m | 2 + 2 p | X˜ t | p + 2 p | g ( X t k −1

k −1

p

+ (2a0 + 2κ 2 ) 2 | X tk−m | p , X tk−1−m ) w tk−1 | p ].

Taking expectation on two sides of the above inequality, yields p

p

E[| X tk | p I [0,σ R ] (k)] ≤ 4 2 −1 {(2a¯ ) 2 E[| X tk−m | p 2

p (α +2) 2

I [0,σ R ] (k)] + 2 p E[| X˜ tk−1 | p I [0,σ R ] (k)]

+(2a0 + 2κ 2 ) E[| X tk−m | p I [0,σ R ] (k)] +2 p E[| g ( X tk−1 , X tk−1−m ) w tk−1 | p I [0,σ R ] (k)]}.

(3.6)

By the Hölder inequality, we may compute 1

1

E[| g ( X tk−1 , X tk−1−m ) w tk−1 | p I [0,σ R ] (k)] ≤ (E[| g ( X tk−1 , X tk−1−m )|2p I [0,σ R ] (k)]) 2 (E| w tk−1 |2p ) 2 .

(3.7)

54

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

¯ Noting that w tk−1 is an m-dimensional Brownian motion, and making use of the property of normal distribution, then ¯ ), dependent on p and m, ¯ such that there exists a constant C ( p , m

¯ ) p . E| w tk−1 |2p = C ( p , m This, together with (3.7), yields p

E[| g ( X tk−1 , X tk−1−m ) w tk−1 | p I [0,σ R ] (k)] ≤ C ( p , R ) 2 . Substituting this into (3.6), one can obtain

E[| X tk | p I [0,σ R ] (k)] ≤ C ( p , R ). 2 ¯ 0 − a¯ 0 > 0, a − mb ¯ − a¯ − m ¯ b¯ > 0, then there exists a sufficiently Theorem 3.2. Assume that (H1), (H2) and (H3) hold with a0 − mb small ∗ ∈ (0, 1) such that the approximate solution { X tk } defined by Eq. (3.1) satisfying

E| X tk |2 ≤

1 1 − κ eετ

Eξ 2 +

C0

(1 − κ )(1 − κ eετ )

,

where

C 0 = E| X˜ (0)|2 + [2κ 2 (eε − 1) + a¯ 0 ]

−1 

eε(i +m) Eξ 2 + a0 Eξ 2

i =−m

¯ ε ) + aEξ α +2 + (¯a + be

−1 

eε(i +m) Eξ α +2 ,

(3.8)

i =−m

¯ ετ ) ¯ 0 eετ − a¯ 0 eετ − 2ε eετ (1 + κ 2 eετ ) = 0, and a − mbe ¯ ετ − eετ (¯a + m ¯ be and ε ≤ 1 ∧ ε1 ∧ ε2 ∧ τ1 log κ1 , ε1 , ε2 are the roots of a0 − mb = 0, respectively. Proof. According to (3.1), we compute

| X˜ tk+1 |2 = X˜ tk+1 , f ( X tk+1 , X tk+1−m ) + X˜ tk + g ( X tk , X tk−m ) w tk 1

≤ X˜ tk+1 , f ( X tk+1 , X tk+1−m )

+ (| X˜ tk+1 |2 + | X˜ tk |2 ) 2

1

+ | g ( X tk , X tk−m )|2 | w tk |2 + X˜ tk , g ( X tk , X tk−m )

w tk . 2

Making use of (H1), yields

| X˜ tk+1 |2 ≤ | X˜ tk |2 − a0 | X tk+1 |2 − a| X tk+1 |α +2 + a¯ | X tk+1−m |α +2 ¯ 0 | X tk |2 + mb ¯ | X t k |α + 2 + m ¯ b¯ | X tk−m |α +2 + sk , + a¯ 0 | X tk+1−m |2 + mb where

¯ ) + 2 X˜ tk , g ( X tk , X tk−m )

w tk . sk = | g ( X tk , X tk−m )|2 (| w tk |2 − m According to (3.8), it is easy to get

eε(k+1) | X˜ tk+1 |2 − eεk | X˜ tk |2 ≤ (eε − 1)eεk | X˜ tk |2 − a0 eε(k+1) | X tk+1 |2 − aeε(k+1) | X tk+1 |α +2

+ a¯ eε(k+1) | X tk+1−m |α +2 + a¯ 0 eε(k+1) | X tk+1−m |2 ¯ (b| X tk |α +2 + b¯ | X tk−m |α +2 ) + eε(k+1) sk . ¯ 0 eε(k+1) | X tk |2 + eε(k+1) m + mb Summing up the above inequality from 0 to k ∧ σ R − 1, we have

(3.9)

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

eε(k∧σ R ) | X˜ tk∧σ |2 ≤ | X˜ (0)|2 + (eε − 1)

k∧ σ R −1

R

eε i | X˜ t i |2 − a0

k∧ σ R −1

i =0

−a

k∧ σ R −1

eε(i +1) | X t i+1 |α +2 + a¯

eε(i +1) | X t i+1−m |α +2

i =0

σ R −1 k∧

¯ 0 eε(i +1) | X t i+1−m |2 + mb

σ R −1 k∧

i =0

¯ +m

e ε ( i + 1 ) | X t i |2

i =0

k∧ σ R −1

eε(i +1) (b| X t i |α +2 + b¯ | X t i−m |α +2 ) +

i =0

k∧ σ R −1

eε(k∧σ R ) | X˜ tk∧σ |2 ≤ | X˜ (0)|2 + 2(eε − 1)[ R

eε i | X t i |2 + a¯ 0

i =1

−a

k ∧σ R

e ε i | X t i |2 + κ 2

k∧σ R −1−m

eε(i +m) | X t i |2 ]

i =−m

σ R −m k∧

eε(i +m) | X t i |2

eε i | X t i |α +2 + a¯

k∧ σ R −m

eε(i +m) | X t i |α +2

i =−m+1

k∧ σ R −1

¯ b¯ e ε ( i + 1 ) | X t i |2 + m

k∧σ R −1−m

i =0

¯ + mb

(3.10)

i =−m+1

i =1

¯ 0 + mb

eε ( i + 1 ) s i .

k∧σ −1−m ε(i +m) k∧σ R −1 ε i e | X t i−m |α +2 = i =−Rm e | X t i |α +2 . i =0

i =0

∧σ R k

k∧ σ R −1 i =0

Noting that | X˜ t i |2 ≤ 2| X t i |2 + 2|u ( X t i−m )|2 ≤ 2(| X t i |2 + κ 2 | X t i−m |2 ), These, together with (3.9), yield

− a0

e ε ( i + 1 ) | X t i + 1 |2

i =0 k∧ σ R −1

i =0

+ a¯ 0

55

eε(i +1+m) | X t i |α +2

i =−m

k∧ σ R −1

e ε ( i + 1 ) | X t i |α + 2 + S k ,

(3.11)

i =0

where S k =

k∧σ R −1 ε(i +1) e si . Rearranging the inequality yields i =0

¯ 0 eε − a¯ 0 eετ + 2eε (1 + κ 2 eετ ) eε(k∧σ R ) | X˜ tk∧σ |2 ≤ | X˜ (0)|2 − [a0 − mb

e− ε − 1

R

¯ ε )] ¯ ε − eετ (¯a + m ¯ be − [a − mbe

k∧ σ R −1

]

k∧ σ R −1

e ε i | X t i |2

i =0

e ε i | X t i |α + 2

i =0

−1 

+ 2κ 2 (eε − 1)

eε(i +m) | X t i |2 − a0 eε(k∧σ R ) | X tk∧σ |2 + a0 | X (0)|2 R

i =−m

− aeε(k∧σ R ) | X tk∧σ R |α +2 + a| X (0)|α +2 + a¯ ¯ b¯ +m

−1 

−1 

eε(i +1+m) | X t i |α +2 + a¯ 0

i =−m

i =−m+1

Denote by

¯ 0 eε − a¯ 0 eετ + 2eε (1 + κ 2 eετ ) A (ε ) = a0 − mb

e− ε − 1

¯ ε ). ¯ ε − eετ (¯a + m ¯ be B (ε ) = a − mbe Noting that

ε < 1, and using the Taylor expansion, we have

e− ε − 1

> −ε +

ε2 2!

−

ε 3 2 3!

+ · · · > −ε ,

,

−1 

eε(i +m) | X t i |α +2

i =−m+1

eε(i +m) | X t i |2 + S k .

(3.12)

56

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

then

¯ 0 eε − a¯ 0 eετ − 2ε eε (1 + κ 2 eετ ) > a0 − mb ¯ 0 eετ − a¯ 0 eετ − 2ε eετ (1 + κ 2 eετ ). A (ε ) > a0 − mb ¯ 0 − a¯ 0 > 0, there exists a Since a0 − mb

ε1 > 0 such that

¯ 0 eετ − a¯ 0 eετ − 2ε eετ (1 + κ 2 eετ ) = 0, a0 − mb ¯ ετ ), since a − mb ¯ ετ − eετ (¯a + m ¯ be ¯ − a¯ − m ¯ b¯ > 0, then there which implies that for ε ≤ ε1 , A (ε ) ≥ 0. Noting B (ε ) > a − mbe exists a ε2 > 0 such that

¯ ετ ) = 0, ¯ ετ − eετ (¯a + m ¯ be a − mbe which implies that for

ε ≤ ε2 , B (ε) ≥ 0. This, together with (3.11), yields

eε(k∧σ R ) E| X˜ tk∧σ |2 ≤ C 0 , R

which is independent of R. Letting R → ∞, yields

eεk E| X˜ tk |2 ≤ C 0 . Recalling X tk = X˜ tk + u ( X tk−m ), for 0 <  < 1, we have

eεk E| X tk |2 ≤

1 1−

eεk E| X˜ tk |2 +

1 εk e E|u ( X tk−m )|2 .



Therefore

max eε i E| X t i |2 ≤ Eξ 2 + max eε i E| X t i |2

−m≤i ≤k

0≤i ≤k

≤ Eξ 2 +

1 1 max eε i E| X˜ t i |2 + max eε i E|u ( X t i−m )|2 . 1 −  0≤i ≤k  0≤i≤k

(3.13)

Noting

max eε i E|u ( X t i−m )|2 ≤ κ 2 max eε i E| X t i−m |2 = κ 2

0≤i ≤k

0≤i ≤k

max

−m≤i ≤k−m

eε(i +m) E| X t i |2 ≤ κ 2 eετ max eε i E| X t i |2 . −m≤i ≤k

This, together with (3.13), yields

max eε i E| X t i |2 ≤ Eξ 2 +

−m≤i ≤k

Letting

1

C0 1−

+

κ 2 ετ e max eε i E| X t i |2 .  −m≤i ≤k

1

 = κ , since ε < τ log κ , then 1 − κ eετ > 0, so we have (1 − κ eετ ) max eε i E| X t i |2 ≤ Eξ 2 + −m≤i ≤k

C0 1−κ

.

This implies that

E| X tk |2 ≤ eεk E| X tk |2 ≤ sup eε i E| X t i |2 ≤ −m≤i ≤k

1 1 − κ eετ

Eξ 2 +

C0

(1 − κ )(1 − κ eετ )

. 2

From the process of the proof, we may see that one-sided polynomial growth condition (H2) plays an important role to guarantee the boundedness of the approximation solution. 4. Strong convergence In this section, our main aim is to establish the convergence theorem on the basis of a series of lemmas. Lemma 4.1. Suppose that (H1), (H2) and (H3) hold. Then for any integer p ≥ 2, there exists a positive constant C ( p , R ) such that

E[| Xˆ tk − X tk | p I [0,σ R ] (k)] ≤ C ( p , R ) p ,

∀k ∈ N,

E[| X˜ˆ tk − X˜ tk | p I [0,σ R ] (k)] ≤ C ( p , R ) p ,

∀k ∈ N.

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

57

Proof. According to (3.1) and (3.2), we have

Xˆ tk+1 − X tk+1 = u ( Xˆ tk+1−m ) − u ( X tk+1−m ) − [u ( Xˆ tk−m ) − u ( X tk−m )]

+ Xˆ tk − X tk + [ f ( X tk , X tk−m ) − f ( X tk+1 , X tk+1−m )] ,

k ≥ 0.

Summing up the above equalities from 0 to k − 1, yields

Xˆ tk − X tk = u ( Xˆ tk−m ) − u ( X tk−m ) + [ f ( X t0 , X t−m ) − f ( X tk , X tk−m )] , where X t−m = x(−τ ) = ξ . For any x, y ∈ Rn , |x| ≤ R , | y | ≤ R, it is easy to see

| f (x, y )|2 ≤ C ( R )(1 + |x|2 + | y |2 ),

(4.1)

where C ( R ) = 2( K R ∨ | f (0, 0)|2 ). Assumption (H2) gives

max E[| Xˆ t i − X t i | p I [0,σ R ] (i )]

1≤i ≤k

≤ 2 p −1 max E[|u ( Xˆ t i−m ) − u ( X t i−m )| p I [0,σ R ] (i )] 1≤i ≤k

+2

p −1

max E[| f ( X t0 , X t−m ) − f ( X t i , X t i−m )| p p I [0,σ R ] (i )]

1≤i ≤k

≤ 2 p −1 κ p max E[| Xˆ t i − X t i | p I [0,σ R ] (i )] 1≤i ≤k

+ 2 p −1 max E[| f ( X t0 , X t−m ) − f ( X t i , X t i−m )| p p I [0,σ R ] (i )] 1≤i ≤k

≤ For any

2

p −1

max E[| f ( X t0 , X t−m ) − f ( X t i , X t i−m )| p p I [0,σ R ] (i )]. 1 − 2 p −1 κ p 1≤i ≤k

κ < 0.5, by (4.1) and Lemma 3.1, there exists a positive constant C ( p , R ), dependent on p and R, such that

E[| Xˆ tk − X tk | p I [0,σ R ] (k)] ≤ max E[| Xˆ t i − X t i | p I [0,σ R ] (i )] ≤ C ( p , R ) p . 1≤i ≤k

The other estimate is easy to obtain similarly.

2

¯ Then there exists a positive integer R 0 ¯ 0 + a¯ 0 , a > a¯ + mb ¯ +m ¯ b. Lemma 4.2. Assume that (H1), (H2) and (H3) hold with a0 > mb sufficiently large and a ∗ sufficiently small such that for every R > R 0 and any < ∗ T∧ρ R

E

q | X˜¯ (s) − X˜ η(s) |q ds ≤ C ( R , T ) 2 , E

0

T∧ρ R

q

| X¯ (s) − X η(s) |q ds ≤ C ( R , T ) 2 ,

q = 1, 2,

0

where ρ R = inf{t > 0, | X¯ (t )| ≥ R or | X η(t ) | > R }. Proof. For s ∈ [tk , tk+1 ), by (3.4),

| X˜¯ (s) − X˜ˆ η(s) | I [t

s  s    ( s ∧ ρ ) = f ( X , X ) dh + g ( X , X ) dw ( h )   I [tk ,tk+1 ) (s ∧ ρ R ) ≤ C ( R ) 1/2 . R t t t t , t ) k k−m k k−m k k +1 tk

tk

58

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

This, together with Lemma 4.1, yields T∧ρ R

E

| X˜¯ (s) − X˜ η(s) |ds ≤ E

0

T∧ρ R

| X˜¯ (s) − X˜ˆ η(s) |ds + E

0

T∧ρ R

| X˜ˆ η(s) − X˜ η(s) |ds ≤ C ( R , T ) 2 . 1

0

Again, by (3.4), for s ∈ [tk , tk+1 ), it is easy to see that

| X¯ (s) − Xˆ η(s) | I [tk ,tk+1 ) (s ∧ ρ R ) = |u ( X¯ (s − τ )) − u ( Xˆ η(s−τ ) )| I [tk ,tk+1 ) (s ∧ ρ R ) s   s   +  f ( X tk , X tk−m )dh + g ( X tk , X tk−m )dw (h) I [tk ,tk+1 ) (s ∧ ρ R ). tk 1

1

tk 1

¯ 2 2 , by Lemma 3.1, we have Note that E[| w (h)|] ≤ [E( w (h))2 ] 2 ≤ m max E[| X¯ (s) − Xˆ η(s) | I [tk ,tk+1 ) (s ∧ ρ R )] ≤ C ( R ) 2 + max κ E[| X¯ (s) − Xˆ η(s) | I [tk ,tk+1 ) (s ∧ ρ R )]. 1

k ≥0

k ≥0

Therefore we have

E| X¯ (s) − Xˆ η(s) | I [tk ,tk+1 ) (s ∧ ρ R ) ≤ (1 − κ )−1 C ( R ) 2 . 1

This implies T∧ρ R

E

| X¯ (s) − Xˆ η(s) |ds ≤ C ( R , T ) 2 . 1

(4.2)

0

By Lemma 4.1 and the Lyapunov inequality, one may get T∧ρ R

E

| Xˆ η(s) − X η(s) |ds ≤ C ( R , T ) .

(4.3)

0

These give T∧ρ R

E

| X¯ (s) − X η(s) |ds ≤ C ( R , T ) 2 . 1

0

The other two estimates are easy to obtain similarly.

2

¯ Then for any given ¯ 0 + a¯ 0 , a > a¯ + mb ¯ +m ¯ b. Lemma 4.3. Assume that (H1), (H2) and (H3) hold with a0 > mb positive integer R 0 sufficiently large and a ∗ sufficiently small such that for every R > R 0 and any < ∗

P(ρ R ≤ T ) ≤ ε ,

ε > 0, there exists a

T > 0.

Proof. The Itô formula implies

E| X˜¯ ( T ∧ ρ R )|2 ≤ E| X˜ 0 |2 − (a0 − b0 )E

T∧ρ R

| X η(s) | ds + (¯a + b¯ )E

0 T∧ρ R

+ a¯ 0 E

T∧ρ R

2

0 T∧ρ R

| X η(s−τ ) )|2 ds − (a − b)E 0

| X η(s−τ ) )|α +2 ds

| X η(s) |α +2 ds

0 T∧ρ R

+ 2C ( R )E 0

| X˜¯ (s) − X˜ η(s) |ds.

(4.4)

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

59

Lemma 4.2 leads to

E| X˜¯ ( T ∧ ρ R )|2 ≤ E| X˜ 0 |2 − (a0 − b0 − a¯ 0 )E

T∧ρ R

| X η(s) | ds + (¯a + b¯ )E

+ a¯ 0 E

| X η(s) |α +2 ds

−τ

0

0

0

2

| X η(s) )|2 ds − (a − b − a¯ − b¯ )E

T∧ρ R

−τ

1

| X η(s) |α +2 ds + C ( R , T ) 2 .

0

¯ yields ¯ 0 + a¯ 0 , a > a¯ + mb ¯ +m ¯ b, Noting that a0 > mb 1 E| X˜¯ ( T ∧ ρ R )|2 ≤ E¯x(0)2 + (¯a + b¯ )τ Eξ α +2 + a¯ 0 τ Eξ 2 + C ( R , T ) 2 .

Case 1: inf{t > 0 : | X¯ (t )| ≥ R } ≤ inf{t > 0 : | X η(t ) | > R }. If | X¯ (ρ R )| = R and | X¯ (ρ R − τ )| < R. Therefore

(4.5)

ρ R = ∞, then the conclusion is obvious. So for ρ R < ∞, we have

1 1 1 | X˜¯ (ρ R )|2 ≥ | X¯ (ρ R )|2 − |u ( X¯ (ρ R − τ ))|2 ≥ R 2 − κ 2 R 2 = ( − κ 2 ) R 2 . 2

2

Case 2: inf{t > 0 : | X¯ (t )| ≥ R } > inf{t > 0 : | X η(t ) | > R }. Clearly,

1

1

2

2

| X˜ ρ R |2 ≥ | X ρ R |2 − |u ( X ρ R −τ )|2 ≥

(4.6)

2

ρ R is the grid point, and | X ρ R | > R , | X ρ R −τ | ≤ R. Moreover,

R2 − κ 2 R2.

Applying the same technique as in Lemma 4.1, it is easy to get

X˜¯ ρ R = X˜ ρ R − f ( X ρ R , X ρ R −τ ) + f ( X t0 , X t−m ) . Using the inequality (a + b)2 ≤ 2a2 + 2b2 and (H2), we have

1

| X˜¯ (ρ R )|2 ≥ | X˜ ρ R − f ( X ρ R , X ρ R −τ ) |2 − | f ( X t0 , X t−m ) |2 2 1

≥ | X˜ ρ R |2 − X˜ ρ R , f ( X ρ R , X ρ R −τ )

− | f ( X 0 , X t−m ) |2 2 1

κ2

4

2

≥( −

)R2 +

a0 2

| X ρ R |2 −

a¯ 2

| X ρ R − τ |α + 2 +

a 2

| X ρ R |α + 2 −

a¯ 0 2

| X ρ R −τ |2 − | f ( X t0 , X t−m )|2 2 .

This implies

| X˜¯ (ρ R )|2 ≥

1 − 2κ 2 + 2a0 − 2a¯ 0 4

R2 +

a − a¯ α +2

− | f ( X t0 , X t−m )|2 2 . R 2

Recalling a0 > a¯ 0 + b0 , a > a¯ and using (4.6) and (4.7), there exist positive constants c 1 and c 2 such that

E| X˜¯ (ρ R )|2 ≥ c 1 R 2 − c 2 . Note that

E| X˜¯ ( T ∧ ρ R )|2 ≥ E[| X˜¯ (ρ R )|2 I {ρ R ≤T } ] ≥ P(ρ R ≤ T )(c 1 R 2 − c 2 ), which implies

P(ρ R ≤ T ) ≤ For any given

1 E|¯x(0)|2 + a¯ 0 τ Eξ 2 + (¯a + b¯ )τ Eξ α +2 + C ( R , T ) 2 . c1 R 2 − c2

ε > 0, choose an N sufficiently large such that for any R ≥ N

E|¯x(0)| + a¯ 0 τ Eξ 2 + (¯a + b¯ )τ Eξ α +2 ε ≤ . 2 c1 R 2 − c2 2

Choose a ∗ < 1/λ sufficiently small such that for any < ∗ 1

C ( R , T ) 2 c1

R2

− c2

ε

≤ , 2

which implies that P(ρ R ≤ T ) ≤ ε as required.

2

(4.7)

60

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

¯ for sufficiently large R, there exists a positive constant ¯ 0 , a > a¯ + mb ¯ +m ¯ b, Lemma 4.4. Let (H1), (H2) and (H3) hold with a0 > a¯ 0 + mb C ( R , T ) such that

E[ sup | X¯ (t ∧ γ R ) − x(t ∧ γ R )|2 ] ≤ C ( R , T ) , 0≤t ≤ T

where γ R = ν R ∧ ρ R . Proof. The local Lipschitz condition and Lemma 4.2 leads to t ∧γ R

t ∧γ R 2

E

| X η(s) − x(s)|2 + | X η(s−τ ) − x(s − τ )|2 )ds

| f ( X η(s) , X η(s−τ ) ) − f (x(s), x(s − τ ))| ds ≤ K R E 0

0 t ∧γ R

[| X η(s) − X¯ (s)|2 + | X¯ (s) − x(s)|2 ]ds

≤ 2K R E 0

t ∧γ R

≤ 2K R C ( R , T ) + 2K R E

| X¯ (s) − x(s)|2 ds.

0

By (2.1) and (3.4), for any T 1 ∈ [0, T ], the Hölder’s and the Burkhölder–Davis–Gundy’s inequalities imply

E[ sup | X¯ (t ∧ γ R ) − x(t ∧ γ R )|2 ] 0≤t ≤ T 1

≤ 3E[ sup |u ( X¯ (t ∧ γ R − τ ) − u (x(t ∧ γ R − τ )|2 + 3T E 0≤t ≤ T 1

T 1 ∧γ R

| f ( X η(s) , X η(s−τ ) ) − f (x(s), x(s − τ ))|2 ds 0

T 1 ∧γ R

| g ( X η(s) , X η(s−τ ) ) − g (x(s), x(s − τ ))|2 ds

+ 12E 0

≤ 3κ E[ sup | X¯ (t ∧ γ R ) − x(t ∧ γ R )|2 + 6( T + 4)(2K R C ( R , T ) + 2K R E

t ∧γ R

2

0≤t ≤ T 1

≤

| X¯ (s) − x(s)|2 ds)

0

12K R ( T + 4)C ( R , T ) 1 − 3κ

2

+

12K R ( T + 4) 1 − 3κ 2

T 1 0

E[ sup | X¯ (s ∧ γ R ) − x(s ∧ γ R )|2 ]ds. 0≤t ≤s

The Gronwall inequality implies that the desired conclusion holds.

2

¯ Then for any given T n = n > 0, n ∈ N, s ∈ [1, 2), ¯ 0 + a¯ 0 , a > a¯ + mb ¯ +m ¯ b. Theorem 4.5. Let (H1), (H2) and (H3) hold with a0 > mb the backward EM scheme (3.1) has the property lim E| X T n − x( T n )|s = 0.

→0

Proof. Let e ( T n ) = X T n − x( T n ), by the Young’s inequality

xs y ≤

δs 2

x2 +

2−s 2

2δ 2−s

2

y 2−s ,

∀x, y , δ > 0,

gives that

E|e ( T n )|s = E[|e ( T n )|s I {ν R >T n ,ρ R >T n } ] + E[|e ( T n )|s I {ν R ≤T n or ρ R ≤T n } ]

≤ 2s−1 E[| X¯ ( T n ) − x( T n )|s I {ν R >T n ,ρ R >T n } ] + E[| X¯ ( T n ) − X T n )|s I {ν R >T n ,ρ R >T n } ] +

δs 2

E[|e ( T n )|2 ] +

2−s 2

2δ 2−s

P(ν R ≤ T n or ρ R ≤ T n ).

(4.8)

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

61

The Hölder inequality and Lemma 4.4 imply that

E[| X¯ ( T n ) − x( T n )|s I {ν R >T n ,ρ R >T n } ] ≤ (E[| X¯ ( T n ) − x( T n )|2 I {ν R >T n ,ρ R >T n } ]) 2 ≤ C ( R , T n ) 2 . s

s

(4.9)

For 0 < s < p < +∞, making use of the Hölder inequality, one may get s

E[| X¯ ( T n ) − X T n |s I {ν R >T n ,ρ R >T n } ] ≤ (E[| X¯ ( T n ) − X T n | p I {ν R >T n ,ρ R >T n } ]) p s

≤ (E[| X¯ ( T n ) − X T n | p I {ρ R >T n } ]) p . Note that to

(4.10)

ρ R > T n means | X¯ ( T n )| < R , | X T n | ≤ R by the definition of ρ R , which implies that n < σ R , then Lemma 4.1 leads

E[| X¯ ( T n ) − X T n | p I {ρ R >T n } ] ≤ C ( p , R ) p . This, together with (4.10), implies s

E[| X¯ ( T n ) − X T n |s I {ν R >T n ,ρ R >T n } ] ≤ (C ( p , R ) p ) p = C ( p , R ) s . This, together with (4.8) and (4.9), yields

s δs 2−s E[|e ( T n )|s ] ≤ 2s−1 C ( R , T n ) 2 + C ( p , R ) s ] + E[|e ( T n )|2 ] + [P(ν R ≤ T n ) + P(ρ R ≤ T n )]. 2 2 2δ 2−s For any given

δs 2

 > 0, by the Hölder inequality and Theorem 2.1 and Lemma 3.2, we choose δ > 0 such that

E[|e ( T n )|2 ] ≤ δ sE[|x( T n )|2 + | X T n |2 ] ≤

 3

.

According to Theorem 2.1, there exists an N 0 sufficiently large such that for any R ≥ N 0

2−s 2δ

s 2−s

P(ν R ≤ T n ) ≤

 3

.

According to Lemma 4.3, choose a < ∗ sufficiently small such that



s

2s−1 C ( R , T n ) 2 + C ( p , R ) s ] + as required.

2−s 2δ

2 2−s

P(ρ R ≤ T n ) ≤

 3

2

In the following, we consider a two-dimensional example to test the above result. Example 4.6. Consider a two dimensional neutral-type stochastic differential delay equation



d(x1 (t ) − a1 x1 (t − τ )) = (−x1 (t ) + a2 x2 (t ) − 2x31 (t ))dt + a3 x2 (t − τ )dw (t ), d(x2 (t ) − a1 x2 (t − τ )) = (−2x2 (t ) + a4 x1 (t ) − 2x32 (t ))dt + a5 x1 (t − τ )dw (t ),

(4.11)

where all ai > 0 (i = 1, 2, 3, 4, 5), τ > 0, w (t ) is a scalar Brownian motion. Let x = (x1 , x2 ), f = ( f 1 , f 2 ), g = ( g 1 , g 2 ), x˜ (t ) = x(t ) − a1 x(t − τ ), y (t ) = x(t − τ ), f 1 (x, y ) = −x1 (t ) + a2 x2 (t ) − 2x31 (t ), f 2 (x, y ) = −2x2 (t ) + a4 x1 (t ) − 2x32 (t ), g 1 (x, y ) = a3 y 2 (t ), g 2 (x, y ) = a5 y 1 (t ). Noting that x21 + x1 x2 + x22 ≥ 2|x1 x2 | + x1 x2 ≥ 0, which implies

x − x¯ , f (x, y ) − f (¯x, y ) = (x1 − x¯ 1 )(−x1 + a2 x2 − 2x31 + x¯ 1 − a2 x¯ 2 + 2x¯ 31 ) + (x2 − x¯ 2 )(−2x2 + a4 x1 − 2x32 + 2x¯ 2 − a4 x¯ 1 + 2x¯ 32 ) = −[(x1 − x¯ 1 )2 + 2(x2 − x¯ 2 )2 ] + (a2 + a4 )(x1 − x¯ 1 )(x2 − x¯ 2 ) − 2(x1 − x¯ 1 )2 (x21 + x1 x¯ 1 + x¯ 21 ) − 2(x2 − x¯ 2 )2 (x22 + x2 x¯ 2 + x¯ 22 ) ≤ −|x − x¯ |2 + 0.5(a2 + a4 )[(x1 − x¯ 1 )2 + (x2 − x¯ 2 )2 ] = −|x − x¯ |2 + 0.5(a2 + a4 )|x − x¯ |2 = (−1 + 0.5(a2 + a4 ))|x − x¯ |2 .

62

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

This implies that f (x, y ) satisfies the one-sided Lipschitz condition. By the Young inequality, it is easy to obtain

x˜ T f (x, y ) = (x1 − a1 y 1 )(−x1 + a2 x2 − 2x31 ) + (x2 − a1 y 2 )(−2x2 + a4 x1 − 2x32 )

= −x21 + a2 x1 x2 − 2x41 + a1 x1 y 1 − a1 a2 x2 y 1 + 2a1 x31 y 1 − 2x22 + a4 x1 x2 − 2x42 + 2a1 x2 y 2 − a1 a4 x1 y 2 + 2a1 x32 y 2 ≤ −x21 + 0.5(a2 + a4 )(x21 + x22 ) − 2x41 + 0.5a1 (x21 + y 21 ) + 0.5a1a2 (x22 + y 21 ) + 0.5a1 (3x41 + y 41 ) − 2x22 − 2x42 + a1 (x22 + y 22 ) + 0.5a1 a4 (x21 + y 22 ) + 0.5a1 (3x42 + y 42 ) = −(1 − 0.5a2 − 0.5a4 − 0.5a1 − 0.5a1a4 )x21 − (2 − 0.5a2 − 0.5a4 − a1 − 0.5a1a2 )x22 + a1 (0.5 + 0.5a2 ) y 21 + a1 (1 + 0.5a4 ) y 22 + 0.5a1 ( y 41 + y 42 ) − (2 − 1.5a1 )(x41 + x42 ). Denote by amin = min{1 − 0.5a2 − 0.5a4 − 0.5a1 − 0.5a1 a4 , 2 − 0.5a2 − 0.5a4 − a1 − 0.5a1 a2 }, amax = max{a1 (0.5 + 0.5a2 ), a1 (1 + 0.5a4 )}, if amin > 0, 2 − 1.5a1 > 0, then

x˜ T f (x, y ) ≤ −amin (x21 + x22 ) + amax ( y 21 + y 22 ) + 0.5a1 ( y 41 + y 42 ) − (2 − 1.5a1 )(x41 + x42 ). Note that 0.5(x21 + x22 )2 ≤ x41 + x42 ≤ (x21 + x22 )2 , so we have

x˜ T f (x, y ) ≤ −amin |x|2 + amax | y |2 + 0.5a1 ( y 21 + y 22 )2 − 0.5(2 − 1.5a1 )(x21 + x22 )2

≤ −amin |x|2 + amax | y |2 + 0.5a1 | y |4 − (1 − 0.75a1 )|x|4 . Moreover, we also have

| g (x, y )|2 = a23 y 22 + a25 y 21 ≤ max{a23 , a25 }( y 21 + y 22 ) = max{a23 , a25 }| y |2 . Clearly, if amin > amax + 0.5 max{a23 , a25 }, 2 − 1.5a1 > a1 , so all conditions in Theorems 2.1 and 4.5 hold, then Eq. (4.11) has a unique global solution with mean-square boundedness and mean-square exponential stability. Moreover, on the basis of Theorem 4.5, the corresponding backward EM approximate solution converges strongly to the true solution. Especially, choose a1 = 0.2, a2 = 0.2, a3 = 0.3, a4 = 0.5, a5 = 0.4, it is easy to see the above conditions satisfy. Let τ = 1, Eq. (4.11) may be written as



d(x1 (t ) − 0.2x1 (t − 1)) = (−x1 (t ) + 0.2x2 (t ) − 2x31 (t ))dt + 0.3x2 (t − 1)dw (t ), d(x2 (t ) − 0.2x2 (t − 1)) = (−2x2 (t ) + 0.5x1 (t ) − 2x32 (t ))dt + 0.4x1 (t − 1)dw (t ).

(4.12)

The corresponding backward EM approximate solutions are defined by

⎧ X 1 (tk+1 ) = X 1 (tk ) + 0.2 X 1 (tk+1−m ) − 0.2 X 1 (tk−m ) + [− X 1 (tk+1 ) ⎪ ⎪ ⎪ ⎨ + 0.2 X 2 (tk+1 ) − 2 X 13 (tk+1 )] + 0.3 X 2 (tk−m ) w k , ⎪ X 2 (tk+1 ) = X 2 (tk ) + 0.2 X 2 (tk+1−m ) − 0.2 X 2 (tk−m ) + [−2 X 2 (tk+1 ) ⎪ ⎪ ⎩ + 0.5 X 1 (tk+1 ) − 2 X 23 (tk+1 )] + 0.4 X 1 (tk−m ) w k ,

(4.13)

where < ∗ , ∗ < λ−1 = 1/0.65, m = −1 , k = 1, 2, ..., N, T = N . In Figs. 1, we plot the backward EM method of Eq. (4.12) for initial value x1 (0) = x2 (0) = 1 (the left figure) and x1 (0) = 1, x2 (0) = 0.5 (the right figure) with τ = 1 for step size = 0.1. In Figs. 2, we plot the backward EM paths of Eq. (4.11) with the initial values x1 (0) = 1, x2 (0) = 0.5, time delay τ = 1 and step size = 0.1 for two different sets of parameters a1 = 0.2, a2 = 0.2, a3 = 0.3, a4 = 0.5, a5 = 0.4 and a1 = 0.5, a2 = 0.4, a3 = 0.2, a4 = 0.6, a5 = 0.4, respectively. These figures illustrate that the numerical solutions have stability and convergence property. 5. Strong convergence rate The previous section establishes the strong convergence in the sense of lim E| X T n − x( T n )|s = 0 for T n = n > 0, n ∈ N,

→0

s ∈ [1, 2). We may observe that the strong convergence rate is still open under the general polynomial growth condition (H1), owing to no results on the p-th moment uniform boundedness of the exact and approximate solutions. In the section, we shall first investigate the uniform boundedness of the true and numerical solutions, then prove the strong convergence with order 1/2. To establish the convergence rate of the backward EM scheme, it is necessary to impose some stronger conditions on the drift and diffusion coefficients. (H4) For x, y , x1 , x2 , y 1 , y 2 ∈ R, there exist a constant

x1 − x2 , f (x1 , y ) − f (x2 , y ) ≤ μ|x1 − x2 |2 ,

μ ∈ R and positive constants μ¯ , l, li > 0 (i = 1, 2) such that (5.1)

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

63

Fig. 1. Computer simulation of the paths x1 (t ) and x2 (t ) of Eq. (4.11) with different parameters.

Fig. 2. Computer simulation of the paths x1 (t ) and x2 (t ) for Eq. (4.11) for different initial datas.

| f (x, y 1 ) − f (x, y 2 )| ≤ v 1 ( y 1 , y 2 )| y 1 − y 2 |, q

q

(5.2)

| f (x1 , y ) − f (x2 , y )| ≤ l(1 + |x1 | + |x2 | )|x1 − x2 |,

(5.3)

¯ |x1 − x2 | + v 2 ( y 1 , y 2 )| y 1 − y 2 |, | g (x1 , y 1 ) − g (x2 , y 2 )| ≤ μ

(5.4)

where v i ( y 1 , y 2 ) = li (1 + | y 1 | + | y 2 || ), i = 1, 2. By (H5), it is easy to see that qi

qi

2 x, f (x, y ) ≤ (1 + 2μ)|x|2 + v 21 ( y , 0)| y |2 ,

(5.5)

¯ |x| + v 2 ( y , 0)| y |. | g (x, y )| ≤ μ

(5.6)

By this and (H1), it is not difficult to observe that show the main result.

α = 2q2 . Moreover, we also need to impose a condition on initial data to

b (H5) For p ≥ 2, ξ ∈ C F ([−τ , 0]; Rn ), there exists positive constant ζ such that 0

E(

sup

−τ ≤s≤t ≤0

|ξ(t ) − ξ(s)|2 ) ≤ ζ (t − s).

Theorem 5.1. Let conditions (H1)–(H5) hold, for any given T > 0 and p ≥ 2, there exists a positive constant c such that E[ sup |x(s)| p ] ≤ c, where c represents a generic positive constant, whose value varies with each appearance throughout the section. 0≤s≤ T

64

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

Proof. Define

νk = inf{t ∈ [0, νe ) : |x(t )| ≥ k}, (k ∈ N ). The Itô formula gives t

|˜x(t )|2 = |˜x(0)|2 +

t x˜ T (s) g (x(s), y (s))dw (s).

(2 ˜x(s), f (x(s), y (s)) + | g (x(s), y (s))|2 )ds + 2 0

0

Applying the same procedure as in Theorem 2.1, we may obtain

t

|˜x(t )| ≤ |˜x(0)| + a¯ 0 ξ  τ + (¯a + b¯ )ξ α +2 τ + 2 2

2

x˜ T (s) g (x(s), y (s))dw (s).

2

0

By this, it is not difficult to get

p

E[ sup |˜x(s ∧ νk )| ] = 2

p /2−1

0≤s≤t

s∧νk

p

(|˜x(0)| + a¯ 0 ξ  τ + (¯a + b¯ )|ξ α +2 τ ) 2 + 2 p −1 E[ sup | 2

2

0≤s≤t

p

x˜ T (ν ) g (x(ν ), y (ν ))dw (ν )| 2 ]. 0

(5.7) By the Burkhölder–Davis–Gundy and the Hölder inequality, we may compute s∧νk

x˜ (ν ) g (x(ν ), y (ν ))dw (ν )| ] ≤ c p E[

E[ sup | 0≤s≤t

t∧νk

p 2

T

0

p

|˜xT (s) g (x(s), y (s))|2 ds] 4 0

≤ c p E[t

t∧νk

p 2 −1

1

|˜xT (s) g (x(s), y (s))| p ds] 2 0 p 2 −1

p

≤ E[ sup |˜x(s)| (t 0≤s≤t

≤2 where c p =

p ( 64 )4 p

for 0 <

−p

t∧νk

c 2p 0

p −2

p

E[ sup |˜x(s ∧ νk )| ] + 2 0≤s≤t

p

p

E[ sup |˜x(s ∧ νk )| ] ≤ c + 0.5E[ sup |˜x(s ∧ νk )| ] + 2

2p −3

0≤s≤t

T

p 2 −1

t∧νk

c 2p E

| g (x(s), y (s))| p ds, 0

p ( p /2) p /2+1 p < 4; c p = 4 for p = 4; c p = ( 2( p /2−1) p/2−1 ) 4 for p

0≤s≤t

1

| g (x(s), y (s))| p ds)] 2

T

p −1 2

> 4. This, together with (5.7), implies

t∧νk

c 2p E

| g (x(s), y (s))| p ds. 0

Since t∧νk

t ∧νk −τ p

t∧νk p

|x(s − τ )| ds =

p

−τ

0

|x(s)| p ds,

|x(s)| ds ≤ ξ  τ + 0

by this and (5.6), it is not difficult to get t∧νk

E[ sup |˜x(s ∧ νk )| p ] ≤ c + c E 0≤s≤t

t∧νk

|x(s)| p ds + c E 0

0

t∧νk

t∧νk

0

|x(s − τ )| p (q2 +1) ds

0

|x(s − τ )| p (q2 +1) ds.

|x(s)| p ds + c E

≤ c + cE

t∧νk

|x(s − τ )| p ds + c E

(5.8)

0

Making use of the Young inequality

(a + b) p ≤ (1 −  )1− p a p +  1− p b p , a, b > 0, 0 <  < 1, we may obtain

(5.9)

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

65

E[ sup |x(s ∧ νk )| p ] = Eξ  p + E[ sup |˜x(s ∧ νk ) + u (x(s ∧ νk − τ ))| p ] −τ ≤s≤t

0≤s≤t

≤ Eξ  + (1 −  )1− p E[ sup |˜x(s ∧ νk )| p ] +  1− p κ p E[ sup |x(s ∧ νk )| p ]. p

−τ ≤s≤t

0≤s≤t

Letting

 = κ , then we have

E[ sup |x(s ∧ νk )| p ] ≤ (1 − κ )−1 Eξ  p + (1 − κ )− p E[ sup |˜x(s ∧ νk )| p ]. −τ ≤s≤t

0≤s≤t

This, together with (5.8), yields t∧νk p

E[ sup |x(s ∧ νk )| ] ≤ a p + b p E 0≤s≤t

t∧νk

|x(s − τ )| p (q2 +1) ds.

p

|x(s)| ds + c p E 0

(5.10)

0

In the following, we adopt the same technique as in Wu et al. [15, Theorem 2.1] and Zong et al. [28, Lemma 4.2]. For any p ≥ 2, define

p i = ([ T /τ ] + 2 − i ) p (q2 + 1)[ T /τ ]+1−i ,

i = 1 , 2 , · · · , [ T /τ ] + 1 ,

it is easy to see that

(q2 + 1) p i +1 < p i , p [T /τ ]+1 = p ≥ 2, For t 1 ∈ [0, τ ], E|x(s − τ )|

p 1 (q2 +1)

≤ Eξ 

i = 1, 2, · · · , [ T /τ ].

p 1 (q2 +1)

p1

E[ sup |x(s ∧ νk )| ] ≤ a p 1 + c p 1 τ Eξ 

for s ∈ [0, t 1 ]. By (5.10), we have p 1 (q2 +1)

0≤s≤t 1

t1 E[ sup |x(u ∧ νk )| p 1 ]ds,

+ b p1 0

0≤ u ≤ s

(5.11)

where a p 1 , b p 1 , c p 1 is similar to a p , b p , c p . The Gronwall inequality yields

E[ sup |x(s ∧ νk )| p 1 ] ≤ (a p 1 + c p 1 T Eξ  p 1 (q2 +1) )eb p1 τ . 0≤ s ≤ τ

Letting k → ∞, the Fatou theorem implies

E[ sup |x(s)| p 1 ] ≤ (a p 1 + c p 1 T Eξ  p 1 (q2 +1) )eb p1 τ ≡ c 1 . 0≤ s ≤ τ

For t 1 ∈ [0, 2τ ], it follows from (5.10) that

t1 p2

E[ sup |x(s ∧ νk )| ] ≤ a p 2 + b p 2 0≤s≤t 1

2τ p2

E|x(s ∧ νk )| ds + c p 2 E 0

|x(s ∧ νk − τ )| p 2 (q2 +1) ds.

(5.12)

0

The Lyapunov inequality gives

τ p2

p1

E[ sup |x(t ∧ νk )| ] ≤ a p 2 + c p 2 0≤t ≤t 1

(E|x(s ∧ νk )| )

p 2 (q2 +1) p1

t1 E|x(s ∧ νk )| p 2 ds

ds + b p 2

−τ

0 p 2 (q2 +1) p1

≤ a p2 + c p2 T c1

t1 E[ sup |x(ν ∧ νk )| p 2 ]ds.

+ b p2 0

0≤ ν ≤ s

The Gronwall inequality and the Fatou Lemma imply p 2 (q2 +1) p1

E[ sup |x(s)| p 2 ] ≤ (a p 2 + c p 2 T c 1 0≤ s ≤ 2 τ

)e 2b p2 τ ≡ c 2 .

Repeating this procedure, we may finally obtain for p [ T /τ ]+1 = p > 2, there exists a constant c [ T /τ ]+1 such that

E[ sup |x(s)| p ] ≤ c [T /τ ]+1 . 0≤ s ≤ T

Choosing c = c [ T /τ ]+1 , yields the desired result.

2

(5.13)

66

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

From the procedure of the proof in Theorem 5.1, we may see that (5.4) plays an important role in establishing the estimation (5.8), where the order of the delay item |x(s − τ )| is higher p, and the order of non-delay item |x(s)| is exactly p, which makes the Gronwall inequality be able to use in estimating the boundedness later. Theorem 5.2. Let conditions (H1)–(H5) hold, for any p ≥ 2, k ≥ 1, there exists a constant c independent on such that E[ sup | X tk | p ] ≤ c. k∈[0, N ]

Proof. By Lemma 3.1, we may obtain immediately

| X˜ tk+1 |2 ≤ a¯ | X tk+1−m |α +2 + a¯ 0 | X tk+1−m |2 + 2| X˜ tk |2 + 2| g ( X tk , X tk−m ) w k |2 . Summing up the above inequalities from 0 to k ∧ σ R , the result is

| X˜ tk∧σ R +1 |2 ≤ | X˜ 0 |2 +

The inequality (

n  i =1

k ∧σ R i =0

xi ) p ≤ n p −1

51− p E[ sup | X˜ tk∧σ 0≤k≤l

| X˜ t i |2 + a¯

n  i =1

k ∧σ R

| X t i+1−m |α +2 + a¯ 0

i =0

k ∧σ R

0≤k≤l

i =0

n  i =1

xi ) p ≤ n p −1

n  i =1

k ∧σ R

| X˜ t i |2 ) p ] + E[ sup a¯ p p ( 0≤k≤l

i =0

k ∧σ R

0≤k≤l

0≤k≤l

| g ( X t i , X t i−m ) w i |2 .

p

p

k ∧σ R

k ∧σ R

xi , p ≥ 1, n ∈ N, yields

+ a¯ 0 p E[ sup (

E[ sup (

| X t i+1−m |2 + 2

i =0

|2p ] ≤ | X˜ 0 |2p + E[ sup ( R +1

Again, applying the inequality (

k ∧σ R

i =0 k ∧σ R

| X t i+1−m |2 ) p ] + 2 p E[ sup ( 0≤k≤l

i =0

| X t i+1−m |α +2 ) p ]

| g ( X t i , X t i−m ) w i |2 ) p ].

(5.14)

i =0

x2i , p ≥ 1, xi > 0 and the Doob martingale inequality, we may compute

| g ( X tk , X t i−m ) w i |2 ) p ] ≤ N p −1 E[ sup

∧σ R k

0≤k≤l

i =0

2p

≤ N p −1 (

2p − 1 2p

≤ N p −1 (

2p − 1

| g ( X t i , X t i−m ) w i |2p ]

i =0

)2p

l ∧σ R

E| g ( X t i , X t i−m ) w i |2p

i =0

)2p

l ∧σ R

E| g ( X t i , X t i−m ) w i |2p .

i =0

2p

Note that E| w i |2p ≤ ( 2p −1 )2p (2p − 1)!! p , by (5.6), we may compute k ∧σ R

E[ sup ( 0≤k≤l

| g ( X t i , X t i−m ) w i |2 ) p ] ≤ N p −1 (

i =0

≤ N p −1 (

≤ c

2p 2p − 1 2p 2p − 1

l ∧σ R i =0

)4p (2p − 1)!! p

l ∧σ R

E| g ( X t i , X t i−m )|2p

i =0

)4p (2p − 1)!! p

l ∧σ R

¯ | X t i | + (1 + | X t i−m |q2 )| X t i−m |)2p E(μ

i =0

E(| X t i |2p + | X t i−m |2p (q2 +1) + | X t i−m |2p ).

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

67

Note that | X˜ t i |2p ≤ | X t i − u ( X t i−m )|2p ≤ 22p −1 | X t i |2p + 22p −1 κ 2p | X t i−m |2p , this, together with (5.14), yields

51− p E[ sup | X˜ tk∧σ 0≤k≤l

k ∧σ R

|2p ] ≤ | X˜ 0 |2p + N p −1 22p −1 E R +1 + N p −1a¯ p p E

k ∧σ R

| X t i |2p + N p −1 22p −1 κ 2p E

i =0

+ c

| X t i−m |2p

i =0

| X t i+1−m | p (α +2) + N p −1a¯ 0 p E p

k ∧σ R

i =0 l ∧σ R

k ∧σ R

| X t i+1−m |2p

i =0

E(| X t i |2p + | X t i−m |2p (q2 +1) + | X t i−m |2p ).

(5.15)

i =0

Note that l ∧σ R

| X t i−m |2p ≤

i =0

−1 

Eξ 2p +

l ∧σ R

i =−m

| X t i |2p ≤ mEξ 2p +

l ∧σ R

i =0

| X t i |2p .

i =0

This implies

E[ sup | X˜ tk∧σ R |2p ] ≤ c + c

l 

0≤k≤l+1

E| X t i∧σ R |2p + c

i =0

+ c

l  i =0

l  i =0

E| X t i∧σ R +1−m | p (α +2)

E| X t i∧σ R +1−m |2p + c

l  i =0

E| X t i∧σ R −m |2p (q2 +1) .

(5.16)

According to x˜ (t ) = x(t ) − u (x(t − τ )) and (5.9), we may compute

E[

sup

−m≤k≤l+1

| X tk∧σ R |2p ] ≤ Eξ 2p + E[ sup | X˜ tk∧σ R + u ( X tk∧σ R −m )|2p ] 0≤k≤l+1

+ (1 −  )1−2p E[ sup | X˜ tk∧σ R |2p ] +  1−2p κ 2p E[ sup | X tk∧σ R −m |2p ]

2p

≤ Eξ 

0≤k≤l+1

1−2p

2p

+ (1 −  )

≤ Eξ  Let

0≤k≤l+1

E[ sup | X˜ tk∧σ R | ] +  2p

1−2p

0≤k≤l+1

2p

κ E[

sup

−m≤k≤l+1

| X tk∧σ R |2p ].

 = κ , so we have E[

sup

−m≤k≤l+1

| X tk∧σ R |2p ] ≤ c + c E[ sup | X˜ tk∧σ R |2p ].

(5.17)

0≤k≤l+1

By (5.16) and (5.17), then we have

E[ sup | X tk∧σ R |2p ] ≤ c + c

l 

0≤k≤l+1

E| X t i∧σ R |2p + c

i =0

+ c

l  i =0

l  i =0

E| X t i∧σ R +1−m | p (α +2)

E| X t i∧σ R +1−m |2p + c

l  i =0

E| X t i∧σ R −m |2p (q2 +1) .

(5.18)

Let p i = ([ N /m] + 2 − i ) p (α + 2)[ N /m]+1−i , i = 1, 2, ..., [ N /m] + 1, it is easy to see that p i +1 (α + 2) < p i and p [ N /m]+1 = p ≥ 1, i = 1, 1, ..., [N/m]. Recalling α = 2q2 , so we have

E[ sup | X tk∧σ R |2p 1 ] ≤ c + c 0≤k≤m

m −1  i =0

+ c

E| X t i∧σ R |2p 1 + c

m −1  i =0

m −1  i =0

E| X t i∧σ R +1−m | p 1 (α +2)

E| X t i∧σ R +1−m |2p 1 + c

m −1  i =0

E| X t i∧σ R −m |2p 1 (q2 +1)

68

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

≤c+c

m −1 

m −1 

0≤k≤i

i =0

+ c

E[ sup | X tk∧σ R |2p 1 ] + c

m −1 

Eξ 2p 1 + c

m −1 

i =0

Eξ  p 1 (α +2)

i =0

Eξ 2p 1 (q2 +1) .

(5.19)

i =0

The discrete Gronwall inequality implies

E[ sup | X tk∧σ R |2p 1 ] ≤ c + c Eξ  p 1 (α +2) + c Eξ 2p 1 + c Eξ 2p 1 (q2 +1) . 0≤k≤m

Letting k → ∞, the Fatou lemma gives

E[ sup | X tk |2p 1 ] ≤ c + c Eξ  p 1 (α +2) + c Eξ 2p 1 + c Eξ 2p 1 (q2 +1) ≡ c 1 . 0≤k≤m

For k ∈ [0, 2m], by (5.18), we have

E[ sup | X tk∧σ R |2p 2 ] ≤ c + c

2m −1 

0≤k≤2m

i =0

+ c

E| X t i∧σ R |2p 2 + c

2m −1  i =0

(E| X t i∧σ R +1−m |2p 1 ) p 2 (α +2)/2p 1

2m −1 

2m −1 

i =0

i =0

(E| X t i∧σ R +1−m |2p 1 ) p 2 / p 1 + c

(E| X t i∧σ R −m |2p 1 ) p 2 (q2 +1)/ p 1 .

Again, by the discrete Gronwall inequality, we have

E[ sup | X tk∧σ R |2p 2 ] ≤ c + c

2m −1 

0≤k≤2m

p (α +2)/2p 1

c1 2

+ c

2m −1 

i =0

p / p1

c1 2

+ c

2m −1 

i =0

p (q2 +1)/ p 1

c1 2

.

i =0

Letting k → ∞, the Fatou lemma gives

E[ sup | X tk |2p 2 ] ≤ c + c 0≤k≤2m

2m −1 

p (α +2)/2p 1

c1 2

+ c

2m −1 

i =0

p / p1

c1 2

+ c

i =0

2m −1 

p (q2 +1)/ p 1

c1 2

≡ c2 .

i =0

Repeating the previous procedures, for any p ≥ 1, we may get for p [ T /τ ]+1 = p ≥ 1, there exists a c [ T /τ ]+1 = c such that

E[ sup | X tk |2p ] ≤ c . 2 0≤k≤ N

Theorem 5.3. Under condition (H1)–(H5), for p ≥ 2, E[ sup | X¯ (t )| p ] ≤ c, where c represents a generic positive constant, whose −τ ≤t ≤ T

value varies with each appearance. Proof. Recalling that (3.4), we may obtain

E[ sup | X˜¯ (t )| p ] ≤ 3 p −1 |˜e (0)|2 + 3 p −1 E[ sup | 0≤u ≤t

0≤u ≤t

+ 3 p −1 E[ sup | 0≤u ≤t

t f ( X η(s) , X η(s−τ ) )ds| p ] 0

u | g ( X η(s) , X η(s−τ ) )dw (s)| p ]. 0

By the Hölder inequality and Theorem 5.2, it is easy to obtain

u p

E[ sup | 0≤u ≤t

f ( X η(s) , X η(s−τ ) )ds| ] ≤ t 0

p −1

t | f ( X η(s) , X η(s−τ ) )| p ds ≤ c .

E 0

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

69

By the Burkhölder–Davis–Gundy and Hölder inequality, we may compute

u E[ sup | 0≤u ≤t

t t p 2 2 g ( X η(s) , X η(s−τ ) )dw (s)| ] ≤ c p E( | g ( X η(s) , X η(s−τ ) )| ds) ≤ c E| g ( X η(s) , X η(s−τ ) )| p ds ≤ c . p

0

0

0

So we may obtain E[ sup | X˜¯ (s)| p ] ≤ c. Again, by (5.9), We may compute −τ ≤s≤t

E[ sup | X¯ (s)| p ] ≤ Eξ  p + E[ sup | X¯˜ (s) + u ( X¯ (s − τ ))| p ] −τ ≤s≤t

0≤s≤t

≤ Eξ  + (1 −  )1− p E[ sup | X˜¯ (s)| p +  1− p κ p | X¯ (s − τ )| p ]. p

0≤s≤t

Let

 = κ , we may arrive at E[ sup | X¯ (s)| p ] ≤ c . −τ ≤s≤t

¯ (t ) − X η(t ) | p ≤ c p /2 for any t ∈ [−τ , T ]. Theorem 5.4. Let condition (H1)–(H5) hold, E| X Proof. The inequality (

n  i =1

xi ) p ≤ n p −1

n  i =1

x2i , p ≥ 1, xi > 0 yields

E| X¯ (s) − X η(s) | p ≤ 2 p −1 E| X¯ (s) − Xˆ η(s) | p + 2 p −1 E| Xˆ (s) − X η(s) | p . According to (3.4) and (5.9) with

 = κ and Theorem 5.2, we may compute

E| X¯ (s) − Xˆ η(s) | p ≤  1− p E|u ( X¯ (s − τ )) − u ( Xˆ η(s−τ ) )| p + (1 −  )1− p 2 p −1 p −1

s E| f ( X tk , X tk−m )| p dh tk

s

+ (1 −  )1− p 2 p −1 c p E(

| g ( X tk , X tk−m )|2 dh) p /2

tk

≤ κ E| X¯ (s − τ ) − Xˆ η(s−τ ) | p + c p + c p /2 . p

p

For s ∈ [−τ , 0], E| X¯ (s) − Xˆ η(s) | p = E|ξ(s) − ξ(tk )| p ≤ ζ 2 2 by (H5), so we have

sup

s∈[−τ , T ]

E| X¯ (s) − Xˆ η(s) | p ≤ sup E| X¯ (s) − Xˆ η(s) | p + sup E| X¯ (s) − Xˆ η(s) | p s∈[−τ ,0] p 2

s∈[0, T ]

p 2

≤ ζ + κ sup E| X¯ (s − τ ) − Xˆ η(s−τ ) | p + c p /2 s∈[0, T ]

≤κ

sup

s∈[−τ , T ]

E| X¯ (s) − Xˆ η(s) | p + c p /2 ≤ c p /2 .

By (5.9) and Lemma 4.1, we have

E| Xˆ η(s) − X η(s) | p = E|u ( Xˆ (s − τ )) − u ( X η(s−τ ) ) + ( f ( X 0 , X −m ) − f ( X tk , X tk−m )) | p ≤  1− p κ p E| Xˆ (s − τ ) − X η(s−τ ) | + (1 −  )1− p c p . Letting

 = κ , then we have sup

s∈[−τ , T ]

E| Xˆ η(s) − X η(s) | p ≤ κ

sup

s∈[−τ , T ]

E| Xˆ (s) − X η(s) | + (1 − κ )1− p c p ≤ (1 − κ )− p c p .

(5.20)

70

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

This, together with (5.20), yields

sup

s∈[−τ , T ]

E| X¯ η(s) − X η(s) | p ≤ c p /2 . 2 p

Theorem 5.5. Assume (H1)–(H5) hold, for any t ∈ [0, T ], E[ sup |x(s) − X¯ (s)| p ] ≤ c 2 , p ≥ 2. 0≤s≤t

Proof. Let e (t ) = x(t ) − X¯ (t ), by (2.1) and (3.4), it is easy to see that

t e˜ (t ) =

t ( f (x(s), x(s − τ )) − f ( X η(s) , X η(s−τ ) ))ds +

0

( g (x(s), x(s − τ )) − g ( X η(s) , X η(s−τ ) ))dw (s). 0

The Itô formula and (5.5) give

t [2e˜ (s)T ( f (x(s), x(s − τ )) − f ( X η(s) , X η(s−τ ) )) + | g (x(s), x(s − τ )) − g ( X η(s) , X η(s−τ ) )|2 ]ds

|˜e (t )|2 ≤ |˜e (0)|2 + 0

t 2e˜ (s) T ( g (x(s), x(s − τ )) − g ( X η(s) , X η(s−τ ) ))dw (s)

+ 0

≡ J 1 + J 2 + J 3, where e˜ (0) = 0. The Hölder inequality implies p 2

E[ sup | J 1 | ] ≤ ( T

t

p −1

2 p |˜e (s)| p | f (x(s), x(s − τ )) − f ( X η(s) , X η(s−τ ) )| p ds)1/2

E

0≤s≤t

0

≤3

− 2p

p

p 2

E[ sup |˜e (s)| ] + 3 T

p −1 p −2

2

0≤s≤t

t | f (x(s), x(s − τ )) − f ( X η(s) , X η(s−τ ) )| p ds.

E 0

Note that

| f (x(s), x(s − τ )) − f ( X η(s) , X η(s−τ ) )| p ≤ 3 p −1 | f (x(s), x(s − τ )) − f ( X¯ (s), x(s − τ ))| p + 3 p −1 | f ( X¯ (s), x(s − τ )) − f ( X¯ (s), X η(s−τ ) )| p + 3 p −1 | f ( X¯ (s), X η(s−τ ) ) − f ( X η(s) , X η(s−τ ) )| p . This, together with (5.2) and (5.3), implies p 2

E[ sup | J 1 | ] ≤ 3 0≤s≤t

− 2p

t p

E[ sup |˜e (s)| ] + c E 0≤s≤t

(1 + |x(s)|q + | X¯ (s)|q ) p |x(s) − X¯ (s)| p ds

0

t (1 + |x(s − τ )|q1 + | X η(s−τ ) |q1 ) p |x(s − τ ) − X η(s−τ ) | p ds

+ cE 0

t +c E

(1 + | X¯ (s)|q + | X η(s) |q ) p | X¯ (s) − X η(s) | p ds.

0

Note that

E|x(s) − X η(s) | p ≤ 2 p −1 (E|x(s) − X¯ (s)| p + E| X¯ (s) − X η(s) | p ) ≤ 2 p −1 (E|x(s) − X¯ (s)| p + c p /2 ) and

E|x(s) − X η(s) | p = E|ξ(s) − ξ(tk )| p ≤ ζ p /2 p /2 for s ∈ [tk , tk+1 ), k = −m, · · ·, −1, by (H3) and Theorems 5.1, 5.2 and 5.4, we may compute

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

p

p

E[ sup | J 1 | 2 ] ≤ 3− 2 E|˜e (s)| p + c E 0≤s≤t

p

≤ 3− 2 E|˜e (s)| p + c E

t

|x(s) − X¯ (s)| p ds + c E

t |x(s) − X η(s) | p ds + c E

0

−τ

t

0

t

t

t

|x(s) − X¯ (s)| p ds + c E

| X¯ (s) − X η(s) | p ds

0

| X¯ (s) − X η(s) | p ds

0

(|x(s) − X¯ (s)| p + | X¯ (s) − X η(s) | p )ds + c E

+ cE

71

0 |x(s) − X η(s) | p ds

−τ

0 p

≤ 3− 2 E|˜e (s)| p + +c E

t

|x(s) − X¯ (s)| p ds + c p /2 .

0

By (5.4), we may compute

u p p E[ sup | J 2 | 2 ] = E[ sup ( | g (x(s), x(s − τ )) − g ( X η(s) , X η(s−τ ) )|2 ds) 2 ] 0≤u ≤t

0≤u ≤t

≤ T p /2−1 E

0

t | g (x(s), x(s − τ )) − g ( X η(s) , X η(s−τ ) )| p ds 0

≤ T p /2−1 E

t ¯ |x(s) − X η(s) | + v 2 (x(s − τ ), X η(s−τ ) )|x(s − τ ) − X η(s−τ ) |] p ds [μ 0

≤2

p −1

T

p −1 2

t ¯ p |x(s) − X η(s) | p + | v 2 (x(s − τ ), X η(s−τ ) )(x(s − τ ) − X η(s−τ ) )| p ]ds [μ

E 0

≤2

p −1

T

p −1 2

t

t

μ¯ |x(s) − X η(s) | ds + E p

(E

p

p

v 2 (x(s), X η(s) )|x(s) − X η(s) | p ds).

−τ

0

According to Theorems 5.1 and 5.2, we may compute p

p

E[ sup | J 2 | 2 ] ≤ 2 p −1 T 2 −1 [E 0≤u ≤t

t

p

¯ p + v 2 (x(s), X η(s) ))|x(s) − X η(s) | p ds (μ 0

0

p

v 2 (x(s), X η(s) )|x(s) − X η(s) | p ds]

+E −τ

t ≤ cE

|x(s) − X¯ (s)| p ds + c E

0

t

| X¯ (s) − X η(s) | p ds + c p /2

0

t ≤ cE

|x(s) − X¯ (s)| p ds + c p /2 .

0

By the Burkhölder–Davis–Gundy and Hölder inequality, we may compute

u

p

0≤s≤t

0≤u ≤t

0

t p ≤ c p E( 4|˜e (s)|2 | g (x(s), x(s − τ )) − g ( X η(s) , X η(s−τ ) )|2 ds) 4 0

p

e˜ (s) T ( g (x(s), x(s − τ )) − g ( X η(s) , X η(s−τ ) ))dw (s)| 2 ]

E[ sup | J 3 | 2 ] ≤ E[ sup |

72

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

t 1 ≤ E( T p /2−1 2 p c 2p |˜e (s)| p | g (x(s), x(s − τ )) − g ( X η(s) , X η(s−τ ) )| p ds) 2 0 p

≤ 3− 2 E sup |˜e (s)| p + c E 0≤s≤t

≤3

≤3

− 2p

(|x(s) − X η(s) | p + |x(s − τ ) − X η(s−τ ) | p )ds 0

t 0 p E sup |˜e (s)| + c E( |x(s) − X η(s) | ds + |x(s) − X η(s) | p ds) p

0≤s≤t

− 2p

t

−τ

0

t p

E sup |˜e (s)| + c E 0≤s≤t

|x(s) − X¯ (s)| p ds + c

0

p

≤ 3− 2 E sup |˜e (s)| p + c E 0≤s≤t

t

t

| X¯ (s) − X η(s) | p ds + c p /2

0

|x(s) − X¯ (s)| p ds + c p /2 .

0

By J 1 , J 2 and J 3 , we have

p

t

2

p

E[ sup |˜e (s)| ] ≤ E[ sup |˜e (s)| ] + c E 3

0≤s≤t

0≤s≤t

|x(s) − X¯ (s)| p ds + c p /2 ≤ c E

0

t

|x(s) − X¯ (s)| p ds + c p /2 .

(5.21)

0

Since x(s) = X¯ (s) = ξ(s) for s ∈ [−τ , 0], by (5.9), we have p

p

p

p

E[ sup |e (s)| p ] ≤ 3 2 −1 (E[ sup | J 1 | 2 ] + E[ sup | J 2 | 2 ] + E[ sup | J 3 | 2 ]) −τ ≤s≤t

0≤s≤t

1− p

≤ (1 −  )

0≤s≤t

1− p

E[ sup |u (x(s − τ )) − u ( X¯ (s − τ ))| p ]

p

1− p

κ E[ sup |x(s) − X¯ (s)| p ]

E[ sup |˜e (s)| ] +  0≤s≤t

1− p

≤ (1 −  )

E[ sup |˜e (s)| ] + 

0≤s≤t

0≤s≤t

−p

≤ (1 − κ )

0≤s≤t

p

p

−τ ≤s≤t

p

E[ sup |˜e (s)| ]. −τ ≤s≤t

This, together with (5.21), implies

t p

E[ sup |e (u )| ] ≤ c E 0≤u ≤t

|x(s) − X¯ (s)| p ds + c p /2 ≤ c

0

t E[ sup |e (u )| p ]ds + c p /2 . 0

0≤ u ≤ s

The Gronwall inequality gives

E[ sup |e (s)| p ] ≤ c p /2 . 2 0≤s≤t

Remark 5.6. Section 4 discusses the strong convergence under very general polynomial growth condition (H1), but the convergence rate is still open, owing to no result on the uniform boundedness of the true and numerical solutions at present. Imposing some stronger conditions on the drift and diffusion coefficients in the fifth section, we obtain the strong convergence rate of one half after the p-th moment uniform boundedness is established. It is easy to see that the p-th moment uniform boundedness is the key result in examining the convergence rate, which greatly depends on conditions (5.2) and (5.4). We may observe that the non-delay item x is linearly growing and the delay item x(t − τ ) may be polynomially growing. Example 5.7. Consider a one-dimensional neutral differential system

d(x(t ) − 0.1x(t − τ )) = (−x(t ) + 0.2x5 (t − τ ) − 2x5 (t ))dt + (0.2x(t ) + 0.1x3 (t − τ ))dw (t ),

(5.22)

where w (t ) is a scalar Brownian motion. Let x˜ = x − 0.1 y , f (x, y ) = −x + 0.2 y 5 − 2x5 , y = x(t − τ ). We may compute

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

73

x˜ T f (x, y ) = (x − 0.1 y )(−x + 0.2 y 5 − 2x5 ) = −x2 + 0.2xy 5 − 2x6 + 0.1xy − 0.02 y 6 + 0.2x5 y

≤ −x2 + 0.2(x6 + 5 y 6 )/6 − 2x6 + 0.1(x2 + y 2 )/2 − 0.02 y 6 + 0.2(5x6 + y 6 )/6 = −0.95x2 + 0.05 y 2 + 0.18 y 6 − 1.8x6 . It is easy to get

| g (x, y )|2 = (0.2x + 0.1 y 3 )2 ≤ 0.4x2 + 0.2 y 6 . Therefore condition (H1) is satisfied. In the following, we consider one-sided Lipschitz condition

x − x¯ , f (x, y ) − f (¯x, y ) = (x − x¯ )(−x − 2x5 + x¯ + 2x¯ 5 ) = −(x − x¯ )2 [1 + 2(x4 + x3 x¯ + x2 x¯ 2 + xx¯ 3 + x¯ 4 )]. Note that (a2 + b2 )2 /2 ≤ a4 + b4 ≤ (a2 + b2 )2 , a, b > 0, then

x4 + x3 x¯ + x2 x¯ 2 + xx¯ 3 + x¯ 4 = (x4 + x¯ 4 ) + xx¯ (x2 + x¯ 2 ) + x2 x¯ 2

≥ (x2 + x¯ 2 )2 /2 + xx¯ (x2 + x¯ 2 ) + x2 x¯ 2 /2 ≥ (x2 + x¯ 2 + xx¯ )2 /2 ≥ 0. That is, the one-sided Lipschitz condition with λ = −1

x − x¯ , f (x, y ) − f (¯x, y ) ≤ −(x − x¯ )2 holds. Similarly,

| f (x1 , y ) − f (x2 , y )| = |(−x1 − 2x51 ) − (−x2 − 2x52 )| = |x1 − x2 + 2(x1 − x2 )(x41 + x31 x2 + x21 x22 + +x1 x32 + x42 )| = [1 + 2(x41 + x31 x2 + x21 x22 + +x1 x32 + x42 )]|x1 − x2 | ≤ (1 + 5x41 + 5x42 )|x1 − x2 |, where the third equality uses x3 y ≤ (3x4 + y 4 )/4 and x2 y 2 ≤ (x4 + y 4 )/2. It is not difficult to obtain

| f (x, y 1 ) − f (x, y 2 )| ≤ 0.5( y 41 + y 42 )| y 1 − y 2 | and

| g (x1 , y 1 ) − g (x2 , y 2 )| ≤ 0.2|x1 − x2 | + 0.1| y 31 − y 32 | = 0.2|x1 − x2 | + 0.3( y 21 + y 22 )| y 1 − y 2 |. Therefore all conditions (H1)–(H5) hold in Theorems 2.1, 4.5 and 5.5 hold, then Eq. (5.22) has a unique global solution with uniform boundedness and mean-square exponential stability. On the basis of Theorems 4.5 and 5.5, the corresponding backward EM approximate solutions

X (tk+1 ) = X (tk ) + 0.1 X (tk+1−m ) − 0.1 X (tk−m ) + (− X (tk+1 ) + 0.2 X 5 (tk+1−m )

− 2 X 5 (tk+1 )) + (0.2 X (tk ) + 0.1 X 3 (tk−m )) w (tk )

(5.23)

with τ = 0.1 converges strongly to the true solution of one half, where < ∗ , ∗ < λ−1 = 1, m = (10 )−1 , k = 1, 2, ..., N, T = N . In Figs. 3, we plot the backward EM method of Eq. (5.22) with T = 10, τ = 0.1, = 0.1 (the left figure) and T = 5,

= 0.05, τ = 1 (the right figure) for initial value x(0) = 1. These figures illustrate that the numerical solutions have stability and convergence property. 6. Conclusion In this paper, we have discussed the strong convergence of the backward Euler methods for neutral-type stochastic differential delay equations with super-linearly growing coefficients. The significant contributions of this paper is the convergence rate of one half for backward Euler method. The strong convergence result under very general polynomial growth condition (H1) is established in Section 4, but the convergence rate is still open under (H1), owing to no results on the p-th moment uniform boundedness of the exact and approximate solutions. To obtain the convergence rate of one half, some stronger conditions are imposed on the drift and diffusion coefficients in the fifth section, the p-th moment uniform boundedness of the exact and approximate solutions are established by making use of Gronwall inequality, which are the key results in examining the convergence rate.

74

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

Fig. 3. Computer simulation of the paths x(t ) with T = 10, τ = 0.1, = 0.1 (the left) and T = 5,

τ = 0.05, = 1 (the right).

Acknowledgements The authors would like to thank the referees for their detailed comments and helpful suggestions. This research was supported in part by the National Natural Science Foundation of China (Grant Nos. 11422110; 11571126). References [1] H. Chen, Y. Zhang, P. Hu, Novel delay-dependent robust stability criteria for neutral stochastic delayed neural networks, Neurocomputing 73 (2010) 2554–2561. [2] M. Hutzenthaler, A. Jentzen, P.E. Kloeden, Strong and weak divergence in finite time of Euler’s method for stochastic differential equations with nonglobally Lipschitz continuous coefficients, Proc. R. Soc. A, Math. Phys. Eng. Sci. 467 (2130) (2011) 1563. [3] Y. Ji, C. Yuan, Tamed EM scheme of neutral stochastic differential delay equations, J. Comput. Appl. Math. 326 (2017) 337–357. [4] X. Mao, Exponential Stability of Stochastic Differential Equation, Dekker, 1994. [5] X. Mao, Razumikhin-type theorems on exponential stability of neutral stochastic functional differential equation, SIAM J. Math. Anal. 28 (1997) 389–401. [6] X. Mao, L. Szpruch, Strong convergence and stability of implicit numerical methods for stochastic differential equations with non-globally Lipschitz continuous coeffcients, J. Comput. Appl. Math. 238 (2012) 14–28. [7] X. Mao, L. Szpruch, Strong convergence rates for backward Euler–Maruyama method for nonlinear dissipative-type stochastic differential equations with super-linear diffusion coefficients, stochastic, an inter, Stoch. Int. J. Probab. Stoch. Process. 123 (2012) 1–28. [8] X. Mao, Y. Shen, A. Gray, Almost sure exponential stability of backward Euler–Maruyama discretization for hybrid stochastic differential equation, J. Comput. Appl. Math. 235 (2011) 1213–1226. ´ Highly nonlinear neutral stochastic differential equations with time-dependent delay and the Euler–Maruyama method, Math. Comput. [9] M. Milo˘sevic, Model. 54 (2011) 2235–2251. [10] M. Milosevic, Almost sure exponential stability of solutions to highly nonlinear neutral stochastic differential equations with time-dependent delay and the Euler–Maruyama approximation, Math. Comput. Model. 57 (2013) 887–899. ´ M. Milo˘sevic, ´ Stability of a class of neutral stochastic differential equations with unbounded delay and Markovian switching and the [11] M. Obradovic, Euler–Maruyama method, J. Comput. Appl. Math. 309 (2017) 244–266. [12] L. Szpruch, X. Mao, D.J. Higham, J. Pan, Numerical simulation of a strongly nonlinear Ait-Sahalia-type interest rate model, BIT Numer. Math. 51 (2011) 405–425. [13] W. Wang, Y. Chen, Mean-square stability of semi-implicit Euler method for nonlinear neutral stochastic delay differential equations, Appl. Numer. Math. 61 (2011) 696–701. [14] F. Wu, X. Mao, Numerical solutions of neutral stochastic functional differential equations, SIAM J. Numer. Anal. 46 (2008) 1821–1841. [15] F. Wu, X. Mao, K. Chen, The Cox–Ingersoll–Ross model with delay and strong convergence of its Euler–Maruyama approximate solutions, Appl. Numer. Math. 59 (2009) 2641–2658. [16] F. Wu, S. Hu, C. Huang, Robustness of general decay stability of nonlinear neutral stochastic functional differential equations with infinite delay, Syst. Control Lett. 59 (2011) 195–202. [17] S. Xu, Y. Chu, J. Lu, Y. Zou, Exponential dynamic output feedback controller design for stochastic neutral systems with distributed delays, IEEE Trans. Syst. Man Cybern., Part A, Syst. Hum. 36 (2006) 540–548. [18] Z. Yan, A. Xiao, X. Tang, Strong convergence of the split-step theta method for neutral stochastic delay differential equations, Appl. Numer. Math. 120 (2017) 215–232. [19] D. Zhang, L. Yu, H ∞ filtering for linear neutral systems with mixed time-varying delays and nonlinear perturbations, J. Franklin Inst. 347 (2010) 1374–1390. [20] L. Zhou, G. Hu, Almost sure exponential stability of neutral stochastic delayed cellular neural networks, J. Control Theory Appl. 6 (2008) 195–200. [21] S. Zhou, Strong convergence and stability of backward Euler Maruyama scheme for highly nonlinear hybrid stochastic differential delay equation, Calcolo 52 (2015) 445–473. [22] S. Zhou, Exponential stability of numerical solution to neutral stochastic functional differential equation, Appl. Math. Comput. 266 (2015) 441–461. [23] S. Zhou, Z. Fang, Numerical approximation of nonlinear neutral stochastic functional differential equation, J. Appl. Math. Comput. 41 (2013) 427–444. [24] S. Zhou, H. Jin, Strong convergence of implicit numerical methods for nonlinear stochastic functional differential equations, J. Comput. Appl. Math. 324 (2017) 241–257. [25] S. Zhou, Z. Wang, D. Feng, Stochastic functional differential equations with infinite delay, J. Math. Anal. Appl. 357 (2009) 416–426. [26] S. Zhou, M. Xue, F. Wu, Robustness of hybrid neutral differential systems perturbed by noise, J. Syst. Sci. Complex. 27 (2014) 1138–1157.

S. Zhou, H. Jin / Applied Numerical Mathematics 140 (2019) 48–75

75

[27] X. Zong, F. Wu, Exponential stability of the exact and numerical solutions for neutral stochastic delay differential equations, Appl. Math. Model. 40 (2016) 19–30. [28] X. Zong, F. Wu, C. Huang, Theta schemes for SDDEs with non-globally Lipschitz continuous coefficients, J. Comput. Appl. Math. 278 (2015) 258–277. [29] X. Zong, F. Wu, C. Huang, Exponential mean square stability of the theta approximations for neutral stochastic differential delay equations, J. Comput. Appl. Math. 286 (2017) 172–185.