Calcolo DOI 10.1007/s10092-016-0173-4
Numerical approximation of stochastic differential delay equation with coefficients of polynomial growth Shaobo Zhou1 · Chaozhu Hu2
Received: 23 November 2014 / Accepted: 20 January 2016 © Springer-Verlag Italia 2016
Abstract Although numerical methods of nonlinear stochastic differential delay equations (SDDEs) have been discussed by several authors, there is so far little work on the numerical approximation of SDDE with coefficients of polynomial growth. The main aim of the paper is to investigate convergence in probability of the EulerMaruyama (EM) approximate solution for SDDE with one-sided polynomial growing drift coefficient and polynomial growing diffusion coefficient. Moreover, we prove the existence-and-uniqueness of almost surely exponentially stable global solution for this nonlinear stochastic delay system. Finally, a computer simulation confirms the efficiency of our numerical method. Keywords Stochastic differential delay equation · Convergence in probability · One-sided polynomial growth conditions · Euler Maruyama method · Stopping time Mathematics Subject Classification
65C20
1 Introduction Stochastic delay models play an important role in science and engineering fields. Most of real world models are nonlinear, such as, population models, neural net-
B
Shaobo Zhou
[email protected] Chaozhu Hu
[email protected]
1
School of Mathematics and Statistics, Huazhong University of Science and Technology, Wuhan 430074, China
2
School of Science, Hubei University of Technology, Wuhan 430068, China
123
S. Zhou, C. Hu
works. Explicit solutions can rarely be obtained for highly nonlinear SDDEs. The numerical approximation for nonlinear stochastic delay system has recently received more and more attention (see [1–6]). The convergence of numerical solution is one of key problems in numerical analysis. Many research efforts were devoted to the convergence of the EM numerical methods for stochastic delay systems with linear growing coefficients several years ago (see [7–12]). Recently, several authors have devoted to strong convergence of implicit EM schemes for nonlinear SDDEs under some specific conditions. Wang et al. [13] developed a split step backward Euler scheme and proved mean-square convergence of approximate solution under conditions which the drift coefficient f (x, y) satisfies onesided Lipschitz condition in x and globally Lipschitz in y, but the diffusion coefficient g(x, y) is globally Lipschitz. Tretyakov et al. [14] established strong convergence of the drift-implicit Euler, balanced and fully implicit mean-square approximate solutions for nonlinear SDDEs. Mao et al. [15] proposed the backward and forward-backward EM method and proved that the numerical solution converges strongly to the true solution to stochastic differential equation ( for short SDE ) with one-sided linear growth condition. Motivated by Mao’s works [15–17], Zhou [18] showed that the implicit EM scheme can preserve boundedness of moment, and the numerical approximation converges strongly to the true solution to SDDE with Markovian switching under onesided polynomial growing and one-sided Lipschitz drift and super-linearly growing diffusion coefficients. The implicit EM scheme is effective on establishing strong convergence of the approximate solution for nonlinear SDEs, but the methods are computationally expensive. Compare with the implicit method, the explicit EM method has its simple algebraic structure, cheap computational cost and may guarantee the acceptable convergence rate under the global Lipschitz condition (see [7,8]). Higham et al. [19] proved strong convergence of the classical EM method under two very strong assumptions that pth moments of both true solution and numerical solution are bounded for SDE without global Lipschitz condition. However, Hutzenthaler et al. [20] proved that in the case of super-linearly growing coefficients, the EM approximation may not converge in the strong L p -sense to the exact solution. So Hutzenthaler et al. [21] developed a tamed EM method and proved the approximate solution converges strongly to SDEs with one-sided Lipschitz drift coefficient and the linear growth diffusion coefficient. Recently, Kumar et al. [22] proved strong convergence of the EM approximation to SDDEs with linear growing non-delay state argument and nonlinear growing delay state argument, which relaxed the classical linear growth condition. Moreover, several authors have devoted to convergence in probability of the EM approximate solutions for nonlinear delay systems. Mao et al. [23] established convergence in probability of the explicit EM approximate solution to SDDE under the Khasminskii-type conditions several years ago. Milo˘s evic´ [24] studied convergence in probability of the EM approximation for neutral SDE with variable delay under the Khasminskii-type conditions. In the paper, our main aim is to establish convergence in probability of the EM approximate solution for SDDE with one-sided polynomial growing drift coefficient and polynomial growing diffusion coefficient. We shall remove the global Lipschitz and linear growth condition and establish a new criterion on convergence of numerical solution under the local Lipschitz and one-sided
123
Numerical approximation of stochastic differential...
polynomial growth condition, which are very general so that many highly nonlinear SDDEs may obey these conditions. We shall also prove the existence-and-uniqueness and almost surely exponential stability of the global solution for this nonlinear stochastic delay system. Finally, a highly nonlinear example is given, and the computer simulation confirms the efficiency of our numerical method. The structure of the paper is as follows: In the next section, we prove that existenceand-uniqueness and almost surely exponential stability to the global solutions for SDDEs with coefficients of the polynomial growth. Section 3 establishes convergence in probability of numerical solution for highly nonlinear SDDE based on a series of lemmas. Finally, we consider a highly nonlinear example to illustrate our theory, which implies the results in the paper are very general and can cover a wider class of nonlinear SDDEs.
2 Global solution 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 p > 0, τ > 0, denote by p L Ft ([−τ, 0]; Rn ) the family of all Ft -measurable and C([−τ, 0]; Rn )-valued random variables ξ such that Eξ p < +∞, where ξ = sup−τ ≤θ≤0 |ξ(θ )|. Let w(t) = (w1 (t), w2 (t), . . . , wm (t)) be m-dimensional Brownian motion. Consider an n-dimensional stochastic differential delay equation dx(t) = f (x(t), x(t − τ ))dt + g(x(t), x(t − τ ))dw(t)
(2.1)
p
on t ≥ 0 with the initial data x0 = ξ ∈ L Ft ([−τ, 0]; Rn ), f : Rn × Rn −→ Rn , g : Rn × Rn −→ Rn×m are locally Lipschitz continuous. Moreover, for convenience, denote by y(t) = x(t − τ ). (H1) (Local Lipschitz Condition) For each integer R ≥ 1, there exists a positive constant C R such that | f (x1 , y1 ) − f (x2 , y2 )|2 ∨ |g(x1 , y1 )−g(x2 , y2 )|2 ≤ C R (|x1 −x2 |2 + |y1 − y2 |2 ) for any xk , yk ∈ Rn with |xk | ∨ |yk | ≤ R (k = 1, 2). The classical condition of the existence-and-uniqueness and stability imposed on the drift and diffusion coefficients is the linear growth condition or one-sided linear growth condition. In the paper, we replace the these classical conditions by one-sided polynomial growth condition.
123
S. Zhou, C. Hu
(H2) For some positive integer L , there exist positive constants a0 , al , a¯ l , a, b0 , bl , b¯l , α, αl such that x T f (x, y) ≤ −a0 |x|2 +
L al |x|αl +2 + a¯ l |y|αl +2 − a|x|α+2 , l=1
|g(x, y)|2 ≤ b0 |x|2 +
L
bl |x|αl +2 + b¯l |y|αl +2
l=1
for any x, y ∈ Rn . Lemma 2.1 (see [25,26]) Let α, b > 0, κ(x) ∈ C(Rn+ ; R), if κ(x) = o(|x|α )(|x| → ∞), then sup [κ(x) − b|x|α ] < ∞.
x∈Rn+
Lemma 2.2 Suppose that x ≥ 0, a > 0, b > 0, ci ≥ 0(0 ≤ i ≤ l), p >0, α > 0, αi ≥ 0(0 ≤ i ≤ l) satisfying 0 < α1 ≤ α2 ≤ · · · ≤ αl ≤ α. If b ≥ c = li=1 ci , then there exists a positive constant p0 satisfying p0 =
0, (α
α
α1 = α,
1 α 1 − α1 )( α1α ) α−α1 ,
α1 = α.
such that a ≥ p0 c and ax p + bx α+ p −
l
ci x αi + p ≥ (a − p0 c)x p .
(2.2)
i=1
Proof If α1 = α, noting b ≥ c, then (2.2) reduces to ax p + bx α+ p − cx α+ p ≥ (a − p0 c)x p . It is easy to see that ax p + bx α+ p − cx α+ p = ax p + (b − c)x α+ p ≥ ax p . Clearly, we only need to choose p0 = 0, so the conclusion holds. Now, we only need to prove that for α1 < α, there exists a positive constant p0 such that μ(x) ≡ a + bx α −
l
ci x αi ≥ (a − p0 c) ≥ 0.
(2.3)
i=1
In the following we shall prove this claim in accordance with two cases.
123
Numerical approximation of stochastic differential...
case 1 l = 1. Obviously, c = c1 , then (2.3) reduces to μ(x) = a + bx α − cx α1 .
(2.4) 1
1 α−α1 on Since α1 < α, then μ (x) = 0 has a unique positive solution x0 = ( cα bα ) (0, +∞), then
α1α1 μ(x0 ) = a − c (α − α1 ) αα α
1 α−α α1 1 c α−α 1 . b
α1
1
Noting that b ≥ c and letting p0 = (α − α1 )( α1α ) α−α1 , then we may obtain μ(x0 ) = a − p0 c
c
α1 α−α1
b
≥ a − p0 c ≥ 0.
On the other hand, μ(0) = a > 0, and by Lemma 2.1, we have μ(∞) ≥ 0. Therefore μ(x) ≥ μ(0) ∧ μ(x0 ) ∧ μ(∞) ≥ 0. Case 2 l > 1. For 0 ≤ x ≤ 1, then we have μ(x) ≥ a + bx α − cx α1 .
(2.5)
Making use of the same technique as one used in case 1, we also obtain μ(x) ≥ 0(x ≥ 0). For x > 1, then we have μ(x) ≥ a + bx α − cx αl . Repeating the same process in case 1, we may also obtain a + bx α − cx αl ≥ a − ν(αl )c.
(2.6)
where ν(y) = (α − y)
yy αα
1 α−y
,
0 < y < α.
It is easy to obtain (ln ν(y)) = α(α − y)−2 ln(y/α) < 0, 0 < y < α. So ν(y) is decreasing on (0, α), then ν(α1 ) ≥ ν(αl ), and a > p0 c = ν(α1 )c ≥ ν(αl )c. That is, μ(x) ≥ a − ν(αl )c ≥ a − ν(α1 )c = a − p0 c ≥ 0(∀x > 1). So the conclusion holds for l > 1 and x ≥ 0.
123
S. Zhou, C. Hu
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 ), where νe is the explosion time, and throughout this paper, we set inf ∅ = ∞(as usual, ∅=the empty set). Theorem 2.3 Let (H1) and (H2) hold, and assume that α > αl (1 ≤ l ≤ L), T > 0, p ≥ 1 and L
p−1 p−1 ¯ al + a¯ l + b0 , a > (bl + bl ) . a0 > 2 2
(2.7)
l=1
Then (i) for any initial data ξ, there almost surely exists a unique global solution x(t) to Eq. (2.1) on t ≥ −τ . Moreover, L p − 1 ¯ p(αl + 2) a¯ l + τ Eξ αl + p . E|x| ≤ Eξ + bl 2 αl + p p
p
l=1
(ii) for any ∈ (0, 1), there exists a sufficiently large integer k ∗ = k ∗ ( , T ) such that P{νk ≤ T } ≤ , k ≥ k ∗ . (iii) the solution x(t) to Eq. (2.1) is almost surely exponentially stable, i.e. lim sup t→∞
ε 1 log|x(t)| ≤ − a.s. t p
where ⎫ ⎤ ⎪ bl + b¯l ⎥ ⎪ ⎪ ⎬ ⎢ p−1 ⎢ ⎥ l=1 . , p a log ⎢1 + b ε = min − ⎥ 0 0 L ⎪ ⎪ ⎣ τ 2 ⎪ ⎪ p−1 ¯ ⎦ −1 ⎪ ⎪ a ¯ (α +2)(α + p) + b ⎭ ⎩ l l l l 2 ⎧ ⎪ ⎪ ⎪ ⎨1
⎡
L a− al + a¯ l +
p−1 2
l=1
(2.8) Proof The proof is rather technical, so we divide the whole proof into two steps.
Step 1 Under(H1), applying the standard truncation technique [see Mao[27], Theorem 3.2.2, P95 ] to Eq. (2.1) for any initial data ξ , 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. 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).
123
Numerical approximation of stochastic differential...
Define V (x) = |x| p , since P(νk ≤ t)V (x(νk )) ≤ EV (x(t ∧νk )), that is, we only need to prove that EV (x(t ∧ νk )) < ∞ since V (x(νk )) = |x(νk )| p = k p → ∞(k → ∞). By the Itô formula, we obtain t∧νk EV (x(t ∧ νk )) = EV (ξ(0)) + E LV (x(s))ds, (2.9) 0
where LV (x) = p|x| p−2 x T f (x, y) +
p p−2 p( p − 2) p−4 T |x| |x| |g(x, y)|2 + |x g(x, y)|2 . 2 2
According to (H2), we may obtain
L 2 − a|x|α+2 + al + LV (x)) ≤ p|x| p−2 −(a0 − p−1 b )|x| 0 2 +
L
p a¯ l +
l=1
p−1 ¯ 2 bl
The inequality |x| p−2 |y|αl +2 ≤
l=1
p−1 2 bl
|x|αl +2
|x| p−2 |y|αl +2 .
p−2 αl + p αl + p |x|
+
αl +2 αl + p , αl + p |y|
(2.10) implies
p − 2 αl + p αl + 2 αl + p |x| |y| + αl + p αl + p αl + 2 αl + p |y| = |x|αl + p + − |x|αl + p . αl + p
|x| p−2 |y|αl +2 ≤
(2.11)
Substituting (2.11) into (2.10), yields
p−1 b0 |x|2 − a|x|α+2 LV (x) ≤ p|x| p−2 − a0 − 2 L p−1 α +2 al + a¯ l + bl + b¯l |x| l + 2 l=1
+
L p − 1 ¯ αl + 2 αl + p |y| p a¯ l + − |x|αl + p bl 2 αl + p l=1
L p − 1 ¯ αl + 2 αl + p |y| = p a¯ l + − |x|αl + p − I0 (x), bl 2 αl + p
(2.12)
l=1
where p−1 b0 )|x| p + pa|x|α+ p 2 L p−1 (bl + b¯l ) |x|αl + p . − p al + a¯ l + 2
I0 (x) = p(a0 −
l=1
123
S. Zhou, C. Hu
L ¯ Recalling that α > αl , a0 > p−1 ¯ l + p−1 l=1 (al + a 2 b0 , a > 2 (bl + bl )), Lemma 2.2 p implies that there exists a¯ 0 such that I0 (x) ≥ a¯ 0 |x| ≥ 0. This, together with (2.9) and (2.12), yields EV (x(t ∧ νk ))
t∧νk L p − 1 ¯ αl + 2 E |y(s)|αl + p −|x(s)|αl + p ds. p a¯ l + ≤ EV (ξ(0)) + bl 2 αl + p 0 l=1
(2.13) According to the integral property, we have
t∧νk 0
≤
(|y(s)|αl + p − |x(s)|αl + p )ds ≤
0
−τ
t∧νk −τ
|x(s)|αl + p ds −
t∧νk
|x(s)|αl + p ds
0
|x(s)|αl + p ds.
(2.14)
This, together with (2.13), yields L p − 1 ¯ αl + 2 τ Eξ αl + p < ∞. p a¯ l + bl 2 αl + p l=1 (2.15) This implies that Eq. (2.1) almost surely admits a unique global solution. Noting that P(νk ≤ t)V (x(νk )) ≤ EV (x(t ∧ νk )) and |x(νk )| = k by the definition of stopping time νk , letting k → ∞ yields EV (x(t ∧ νk )) ≤ Eξ p +
lim sup P(νk ≤ t) = 0. k→∞
That is, for any given ∈ (0, 1) and T > 0, there exists a sufficiently large integer k ∗ = k ∗ ( , T ) such that P{νk ≤ T } ≤ , k ≥ k ∗ . Moreover, letting k → ∞ in (2.15), yields L p − 1 ¯ p(αl + 2) a¯ l + τ Eξ αl + p . EV (x(t)) ≤ Eξ + bl 2 αl + p p
l=1
Step 2 In the following, we shall prove that Eq. (2.1) is almost surely exponentially stable. By the Itô formula, we have eεt V (x(t))
t
= V (x(0)) +
¯ eεs [LV (x(s)) + εV (x(s))]ds + M(t),
(2.16)
0
t ¯ where M(t) = 0 eεs Vx (x(s))g(x(s), y(s))dw(s) is a real-valued continuous local ¯ martingale with M(0) = 0. Applying the same technique as one used in (2.11), it is easy to obtain
123
Numerical approximation of stochastic differential...
|x| p−2 |y|αl +2 ≤
p − 2 + (αl + 2)eετ αl + p αl + 2 + |x| (|y|αl + p − eετ |x|αl + p )]. αl + p αl + p
So we may also estimate L p − 1 ¯ αl + 2 αl + p |y| p a¯ l + − eετ |x|αl + p − I¯0 (x), LV (x) + εV (x) ≤ bl 2 αl + p l=1 (2.17) where I¯0 (x) = p a0 − + a¯ l +
p−1 2 b0 p−1 ¯ 2 bl
L |x| p + pa|x|α+ p − p al + l=1 p−2+(αl +2)eετ αl + p . |x| αl + p
−
ε p
p−1 2 bl
(2.18)
By conditions (2.7) and (2.8), it is easy to see that p−1 ε b0 − > 0, 2 p L p−1 p − 1 ¯ p − 2 + (αl + 2)eετ al + bl + a¯ l + > 0. a− bl 2 2 αl + p a0 −
l=1
By Lemma 2.2, there exists a positive aˆ 0 such that I¯0 (x) ≥ aˆ 0 |x| p . This, together with (2.17), yields eεt V (x(t)) ≤ V (x(0)) +
L p − 1 ¯ αl + 2 p a¯ l + bl 2 αl + p l=1
t
×
¯ eεs |y(s)|αl + p − eετ |x(s)|αl + p ds + M(t).
(2.19)
0
Making use of the property of the integral, we may estimate
t
eεs (|y(s)|αl + p − eετ |x(s)|αl + p )ds t = eεs [|x(s − τ )|αl + p − eετ |x(s)|αl + p ]ds 0 t t ≤ eε(s+τ ) |x(s)|αl + p ds − eε(s+τ ) |x(s)|αl + p ds
0
−τ 0
=
−τ
(2.20)
0
eε(s+τ ) |x(s)|αl + p ds.
123
S. Zhou, C. Hu
Substitution this into (2.19) yields L p − 1 ¯ αl + 2 p a¯ l + e V (x(t)) ≤ V (x(0)) + bl 2 αl + p l=1 0 ¯ eε(s+τ ) |x(s)|αl + p ds + M(t). × εt
−τ
Applying the nonnegative semi-martingale convergence theorem (see [28]), we have lim sup eεt V (x(t)) < ∞ a.s. t→∞
That is, there is a finite positive random variable C0 such that lim sup0≤t<∞ eεt V (x(t)) ≤ C0 a.s. This implies lim supt→∞ 1t log|x(t)| ≤ − εp a.s.
3 Numerical method In the section, we shall establish convergence in probability of the EM approximation to Eq. (2.1). We shall first establish several lemmas, then give the main theorem. Now we define the EM approximate solution to Eq. (2.1) as follows z k+1 = z k + f (z k , z k−m¯ ) + g(z k , z k−m¯ )wk ,
(3.1)
where ≥ 0 is a step size which satisfies τ = m ¯ for some integer m. ¯ Let tk = k(k ≥ 0). z k is an approximation to x(tk ) if tk ≤ 0, we have z k = ξ(tk ). Moreover, the increments wk = w(tk+1 ) − w(tk ), k = 1, 2, ..., are independent N (0, )distributed Gaussian random variables Ftk -measurable at the mesh-point tk . For each k > 0, define z¯ (t) = z k , z¯ (t − τ ) = z k−m¯ for t ∈ [tk , tk+1 ), with the initial value z¯ 0 = ξ0 on [−τ, 0]. That is, z¯ (t) =
∞
z k I[k,(k+1)) (t).
k=−m¯
While the continuous EM approximate solution z(t) is to be interpreted as the stochastic integral ⎧ ⎨ ξ, t ∈ [−τ, 0], t t z(t) = f (¯z (s), z¯ (s − τ ))ds + g(¯z (s), z¯ (s − τ ))dw(s), t ∈ [0, T ]. ⎩ ξ(0) + 0
0
(3.2)
Therefore, we have for any t ≥ 0 z(t) = z(k) +
t k
123
f (¯z (s), z¯ (s − τ ))ds +
t
k
g(¯z (s), z¯ (s − τ ))dw(s).
Numerical approximation of stochastic differential...
Clearly, z(tk ) = z k = z¯ (tk ), that is, z(t) and z¯ (t) coincide with the discrete approximate solution at the grid points. Let C be a positive constant independent of , and the product of C and other constants is still denoted by C in the section. Lemma 3.1 Let (H1) hold, then there exists a positive constant C such that sup |z(t ∧ ρ R )|2 ≤ C,
E
−τ ≤t≤T
where ρ R = σ R ∧ ν R , σ R = inf{t ≥ 0, |z(t)| ≥ R}, ν R = inf{t ≥ 0, |x(t)| ≥ R}. Proof It is easy to obtain that E
sup |z(t ∧ ρ R )|
≤ Eξ + E
2
sup |z(t ∧ ρ R )|
2
−τ ≤t≤T
2
.
(3.3)
0≤t≤T
Applying the inequality |a + b + c|2 ≤ 3(|a|2 + |b|2 + |c|2 ), yields E
sup |z(t ∧ ρ R
)|2
≤
3Eξ 2
+ 3E
0≤t≤T
sup |
t∧ρ R 0
0≤t≤T
sup |
+3E
0≤t≤T
t∧ρ R 0
f (¯z (s), z¯ (s
− τ ))ds|2
g(¯z (s), z¯ (s − τ ))dw(s)|2 .
(3.4) According to (H1) and the inequality |a + b|2 ≤ 2(|a|2 + |b|2 ), for |x1 | ∨ |x2 | ≤ R, it is easy to obtain | f (x1 , x2 )|2 ≤ 2| f (x1 , x2 ) − f (0, 0)|2 + 2| f (0, 0)|2 ≤ K R (1 + |x1 |2 + |x2 |2 ). (3.5) Similarly, |g(x1 , x2 )|2 ≤ K R (1 + |x1 |2 + |x2 |2 ), where K R = 2(C R ∨ | f (0, 0)|2 ∨ |g(0, 0)|2 ). Using the Hölder inequality (E|X T Y |)2 ≤ E|X |2 E|Y |2 , we may compute E
t∧ρ R
sup | 0≤t≤T
f (¯z (s), z¯ (s − τ ))ds|2
0 T ∧ρ R
≤ TE
| f (¯z (s), z¯ (s − τ ))| ds 2
0
≤ T K RE ! ≤ T KR
T ∧ρ R
(3.6)
(1 + |¯z (s)| + |¯z (s − τ )| )ds " T E sup |z(ν)|2 ds . T +2 2
2
0
0
−τ ≤ν≤s∧ρ R
123
S. Zhou, C. Hu
Applying the Burkhölder-Davis-Gundy inequality, yields
t∧ρ R
sup |
E
g(¯z (s), z¯ (s − τ ))dw(s)|
2
0
0≤t≤T
T ∧ρ R
≤ 4E
|g(¯z (s), z¯ (s − τ ))|2 ds
0
!
T
≤ 4K R T + 2
E
sup
−τ ≤ν≤s∧ρ R
0
|z(ν)|
2
" ds .
(3.7)
Substituting (3.6) and (3.7) into (3.4), then we have sup |z(t ∧ ρ R )|
E
2
0≤t≤T
!
T
≤ 3Eξ 2 + 3(T + 4)K R T + 2
E
sup
−τ ≤ν≤s∧ρ R
0
"
|z(ν)|2 ds .
This, together with (3.3), yields E
sup |z(t ∧ ρ R )|
2
−τ ≤t≤T
!
T
≤ 4Eξ 2 + 3(T + 4)K R T + 2
" E
0
sup
−τ ≤ν≤s∧ρ R
≤ 4Eξ + 3(T + 4)K R T + 6(T + 4)K R 2
0
T
E|z(ν)|2 ds
E
sup |z(ν ∧ ρ R )|
2
−τ ≤ν≤s
ds.
The Gronwall inequality implies E[ sup |z(t ∧ ρ R )|2 ] ≤ (3(T + 4)K R T + 4Eξ 2 )e6T (T +4)K R ≡ C. −τ ≤t≤T
Lemma 3.2 Under (H1), for any t ∈ [0, T ], there exist a positive constant C, independent of such that
t∧ρ R
E|z(s) − z¯ (s)|2 ds ≤ C.
0
Proof For any t ∈ [tk , tk+1 ), recalling that the definition of z(t) by z(t) = ξ(0) + 0
123
t
t
f (¯z (s), z¯ (s − τ ))ds + 0
g(¯z (s), z¯ (s − τ ))dw(s).
Numerical approximation of stochastic differential...
Since z¯ (t) = z k for t ∈ [tk , tk+1 ), it is obvious that
k
z¯ (t) = ξ(0) +
k
f (¯z (s), z¯ (s − τ ))ds +
0
g(¯z (s), z¯ (s − τ ))dw(s).
0
According to Lemma 3.1, it is not difficult to obtain E[sup−τ ≤t≤T ∧ρ R |¯z (t)|2 ] ≤ C. Since z¯ (t) = z k , z¯ (t − τ ) = z k−m¯ for t ∈ [tk , tk+1 ), therefore z(t) − z¯ (t) =
t k
f (¯z (s), z¯ (s − τ ))ds +
t
g(¯z (s), z¯ (s − τ ))dw(s)
(3.8)
k
= f (¯z (t), z¯ (t − τ ))(t − tk ) + g(¯z (t), z¯ (t − τ ))(w(t) − w(tk )). Applying the inequality |a + b|2 ≤ 2(|a|2 + |b|2 ) and (3.5), yields |z(t) − z¯ (t)|2 ≤ 2| f (¯z (t), z¯ (t − τ ))|2 2 + 2|g(¯z (t), z¯ (t − τ ))|2 |w(t) − w(tk )|2 ≤ 2K R (1+|¯z (t)|2 +|¯z (t −τ )|2 )(2 +|w(t)−w(tk )|2 ), t ∈ [tk , tk+1 ). (3.9) Noting that E|w(s) − w(tk )|2 = m(s − tk ), s ∈ (tk , tk+1 ), Lemma 3.1 implies
t∧ρ R
E|z(s) − z¯ (s)|2 ds 0 t∧ρ R ≤ 2K R 1 + E|¯z (s)|2 + E|¯z (s − τ )|2 2 + E|w(s) − w(tk )|2 ds 0 t∧ρ R ≤ 2K R 1 + E|¯z (s)|2 + E|¯z (s − τ )|2 2 + m ds 0 ≤ 2K R T (1 + 2C) 2 + m = C. Lemma 3.3 Under the conditions of Theorem 2.3, for any ∈ (0, 1) and T > 0 there exists a sufficiently large R¯ = R( , T ) and sufficiently small ∗ such that ¯ < ∗ . P{σ R ≤ T } ≤ , R ≥ R, Proof Recalling that (3.2) and applying the Itô formula, we may compute dV (z(t)) = Vx (z(t)) f (¯z (t), z¯ (t − τ ))dt + Vx (z(t))g(¯z (t), z¯ (t − τ ))dw(t) 1 + trace[g T (¯z (t), z¯ (t − τ ))Vx x (z(t))g(¯z (t), z¯ (t − τ ))]dt 2 ≤ {LV (z(t)) + p|z(t)| p−2 z T (t)[ f (¯z (t), z¯ (t − τ )) − f (z(t), z(t − τ ))] p( p−1) + |z(t)| p−2 |g(¯z (t), z¯ (t −τ ))−g(z(t), z(t −τ ))||g(¯z (t), z¯ (t −τ ))| 2
123
S. Zhou, C. Hu
p( p−1) |z(t)| p−2 |g(z(t), z(t −τ ))||g(¯z (t), z¯ (t −τ ))−g(z(t), z(t −τ ))|}dt 2 +Vx (z(t))g(¯z (t), z¯ (t − τ ))dw(t)
+
≡ [LV (z(t)) + F(z(t), z(t − τ ), z¯ (t), z¯ (t − τ ))]dt +Vx (z(t))g(¯z (t), z¯ (t − τ ))dw(t),
(3.10)
where |¯z (t)| ∨ |¯z (t − τ )| ∨ |z(t)| ∨ |z(t − τ )| ≤ R z, z¯ ∈ Rn . Lemma 3.1 and the local Lipschitz conditions implies that there exists a constant C such that F(z(t), z(t − τ ), z¯ (t), z¯ (t − τ )) ≤ C[|z(t) − z¯ (t)| + |z(t − τ ) − z¯ (t − τ )|]. Therefore we obtain immediately dV (z(t)) ≤ LV (z(t))dt + C[|z(t) − z¯ (t)| + |z(t − τ ) − z¯ (t − τ )|]dt +Vx (z(t))g(¯z (t), z¯ (t − τ ))dw(t).
(3.11)
Integrating the above inequality with respect to t, yields V (z(t ∧ σ R ))
t∧σ R t∧σ R ≤ V (ξ(0)) + LV (z(s))ds + Vx (z(s))g(¯z (s), z¯ (s − τ ))dw(s) 0 0 t∧σ R +C [|z(s) − z¯ (s)| + |z(s − τ ) − z¯ (s − τ )|]ds. (3.12) 0
Recalling that (3.8), for s ∈ (0, t ∧ σ R ), we obtain E[|z(s) − z¯ (s)|] = E[| f (¯z (s), z¯ (s − τ ))(s − tk ) + g(¯z (s), z¯ (s −τ ))(w(s)−w(tk ))|] ≤ E[| f (¯z (s), z¯ (s − τ ))| + |g(¯z (s), z¯ (s − τ ))||(w(s) − w(tk ))|]. (3.13) 1 1 Applying the Hölder inequality (E|x|r ) r ≤ (E|x| p ) p , 0 < r < p < ∞, we have 1
1
E| f (¯z (s), z¯ (s − τ ))| ≤ (E| f (¯z (s), z¯ (s − τ ))|2 ) 2 ≤ [K R E(1 + 2C)] 2 . E|g(¯z (s), z¯ (s − τ ))|2 ≤ K R E(1 + 2C). m 2 |w j (s) − w j (tk )|2 = m(s − tk ) ≤ m. E|w(s) − w(tk )| = E j=1 1
1
The Hölder inequality E|X T Y | ≤ (E|X | p ) p (E|Y |q ) q for 1/ p + 1/q = 1, p, q > 1 implies 1 1 2 2 E|w(s)−w(tk )|2 E(|g(¯z (s), z¯ (s −τ ))||w(s)−w(tk ))|) ≤ E|g(¯z (s), z¯ (s −τ ))|2 1 1 2 ≤ K R E(1 + 2C)) 2 (m .
123
Numerical approximation of stochastic differential...
Substituting these into (3.13), the result is 1
1
1
E|z(s) − z¯ (s)| ≤ [K R (1 + 2C)] 2 ( + (m) 2 )] = C 2 . Therefore
t∧σ R
1
E[|z(s) − z¯ (s)| + |z(s − τ ) − z¯ (s − τ )|]ds ≤ 2C T 2
(3.14)
(3.15)
0
Substituting (3.15) into (3.12), yields 1
EV (z(t ∧ σ R )) ≤ V (ξ(0)) + 2C T 2 + E
t∧σ R 0
LV (z(s))ds.
(3.16)
Repeating the procedure from Theorem 2.3, we can prove that L p−1 ¯ αl +2 a¯ l + 2 bl αl + p E EV (z(t ∧ σ R )) ≤ EV (ξ(0))+2C T + 1 2
1 ≤ M¯ 0 + C T 2 ,
l=1
0 −τ
|z(s)|αl + p ds
(3.17) where M¯ 0 be positive constant. Since P(σ R ≤ t)V (z(σ R )) ≤ EV (z(t ∧ σ R )), then 1 P(σ R ≤ t)V (z(σ R )) ≤ M¯ 0 + C T 2 . Noting that |z(σ R )| = R by the definition of stopping time σ R , therefore for arbitrary ∈ (0, 1), there exists a sufficiently large R¯ = R( , T ) such that M¯ 0 /R p ≤ /2. Moreover, we may choose sufficiently small ∗ such that
1
CT 2 Rp
≤ 2 . Therefore P{σ R ≤ T } ≤ ,
¯ R ≥ R,
≤ ∗ .
Lemma 3.4 Under the conditions of Theorem 2.3, the EM approximate solution converges to the exact solution to Eq. (2.1) in the sense lim E
sup |x(t ∧ ρ R ) − z(t ∧ ρ R )|2 = 0.
→0
0≤t≤T
Proof By the definitions of x(t) and z(t), for any t1 ∈ [0, T ], the inequality (a +b)2 ≤ 2(a 2 + b2 ) implies E
sup |x(t ∧ ρ R ) − z(t ∧ ρ R )|2 0≤t≤t1
# t∧ρ #2 R # # # ≤ 2E sup # ( f (x(s), x(s − τ )) − f (¯z (s), z¯ (s − τ )))ds ## 0≤t≤t1 0 # t∧ρ #2 R # # # +2E sup # (g(x(s), x(s − τ )) − g(¯z (s), z¯ (s − τ )))dw(s)## . 0≤t≤t1
0
123
S. Zhou, C. Hu
The Hölder inequality, (H1) and Lemma 3.2, compute
t∧ρ R
E sup | 0≤t≤t1
0 t1 ∧ρ R
≤ TE 0
| f (x(s), x(s − τ )) − f (¯z (s), z¯ (s − τ ))|2 ds
t1 ∧ρ R
≤ T K RE
( f (x(s), x(s − τ )) − f (¯z (s), z¯ (s − τ )))ds|2
0
≤ T K RE
t1 ∧ρ R
|x(s) − z¯ (s)|2 + |x(s − τ ) − z¯ (s − τ )|2 ds
(2|x(s) − z(s)|2 + 2|z(s) − z¯ (s)|2
0
+2|x(s − τ ) − z(s − τ )|2 + 2|z(s − τ ) − z¯ (s − τ )|2 )ds t1 ∧ρ R (|x(s) − z(s)|2 + |z(s) − z¯ (s)|2 )ds ≤ 4T K R E 0 T
≤ 4T K R
E
sup |x(ν ∧ ρ R ) − z(ν ∧ ρ R )|2 ds + 4T K R C. 0≤ν≤s
0
Applying the Burkhölder-Davis-Gundy inequality implies
t∧ρ R
E sup | 0≤t≤t1
(g(x(s), x(s − τ )) − g(¯z (s), z¯ (s − τ )))dw(s)|2
0
T
≤ 16K R
E
sup |x(ν ∧ ρ R ) − z(ν ∧ ρ R )|
ds + 16K R C.
2
0≤ν≤s
0
Therefore E
sup |x(t ∧ ρ R ) − z(t ∧ ρ R )|
2
0≤t≤T
T
≤ 8(T + 4)K R 0
E
sup |x(ν ∧ ρ R ) − z(ν ∧ ρ R )|2 ds + 8(T + 4)K R C. 0≤ν≤s
The Gronwall inequality implies E
sup |x(t ∧ ρ R ) − z(t ∧ ρ R )|
2
≤ 8(T + 4)K R Ce8K R T (T +4) .
0≤t≤T
Thus lim→0 E[sup0≤t≤T |x(t ∧ ρ R ) − z(t ∧ ρ R )|2 ] = 0.
Lemma 3.5 Under the condition of Theorem 2.3, the EM approximate solution converges to the exact solution of Eq. (2.1) in the sense
123
Numerical approximation of stochastic differential...
lim E
→0
sup |x(t ∧ ρ R ) − z¯ (t ∧ ρ R )|
2
= 0.
0≤t≤T
Proof Recalling that (3.9), for p > 1, t ∈ [tk , tk+1 ), Let N = [T /], we may compute E
sup |z(t ∧ ρ R ) − z¯ (t ∧ ρ R )|
2
0≤t≤T
≤ C + CE 2
≤ C + C 2
sup
m
E
j=1
≤ C + C 2
m
sup
≤ C + C
sup
|w j (t) − w j (tk )|
2
! E
m N j=1
|w(t) − w(tk )|
k=0,1,2,...,N tk ≤t≤tk+1 ∧T
j=1 2
sup
k=0,1,2,...,N tk ≤t≤tk+1 ∧T
2
"1/ p sup
sup
k=0,1,2,...,N tk ≤t≤tk+1 ∧T
|w j (t) − w j (tk )|
! E
k=0
2p
"1/ p sup
tk ≤t≤tk+1 ∧T
|w j (t) − w j (tk )|
2p
.
E|w j (t) − w j (s)|2 p = (2 p − 1)!!(t − s) p and the Doob martingale inequality implies ! E
" sup
tk ≤t≤tk+1 ∧T
|w j (t) − w j (tk )|
2p
≤ =
2p 2p − 1 2p 2p − 1
2 p E|w j (tk+1 ∧ T ) − w j (tk )|2 p 2 p (2 p − 1)!! p .
Noting that the inequality ab ≤ (a + b)2 /4 implies [(2 p − 1)!!]1/ p ≤ p, therefore we have E
sup |z(t ∧ ρ R ) − z¯ (t ∧ ρ R )|2 0≤t≤T
1/ p 2 p 2 p (2 p − 1)!! p ≤ C2 + Cm (N + 1) 2p − 1 2 2p 2 1/ p = C + Cm(T + 1) [(2 p − 1)!!]1/ p ( p−1)/ p 2p − 1 2 2p = C2 + Cmp(T + 1)1/ p ( p−1)/ p . 2p − 1
123
S. Zhou, C. Hu
This, together with Lemma 3.4, yields E
sup |x(t ∧ ρ R ) − z¯ (t ∧ ρ R )|
2
≤ 2E
0≤t≤T
sup |x(t ∧ ρ R ) − z(t ∧ ρ R )|
2
0≤t≤T
+2E
sup |z(t ∧ ρ R ) − z¯ (t ∧ ρ R )|
2
.
0≤t≤T
Taking the limit with respect to → 0, yields lim E
→0
sup |x(t ∧ ρ R ) − z¯ (t ∧ ρ R )|
2
= 0.
0≤t≤T
Theorem 3.6 Under the conditions of Theorem 2.3, for arbitrary T > 0, lim sup |x(t) − z(t)| = 0 in pr obabilit y.
→0 t∈[0,T ]
lim sup |x(t) − z¯ (t)| = 0 in pr obabilit y.
→0 t∈[0,T ]
Proof Let ∈ (0, 1). Denote by = {ω : sup |x(t) − z(t)| ≥ δ}. According to t∈[0,T ]
Theorem 2.3 and Lemma 3.3, for arbitrary ∈ (0, 1), there exists a sufficiently large ¯ such that R¯ = R( , T ) and sufficiently small P{σ R¯ ≤ T } ≤ /3, P{ν R¯ ≤ T } ≤ /3. We may compute P() ≤ P{ ∩ {ρ R¯ > T }} + P{ρ R¯ ≤ T } ≤ P{ ∩ {ρ R¯ > T }} + P{σ R¯ ≤ T } + P{ν R¯ ≤ T }.
(3.18)
Again, we may compute ! δ 2 P{ ∩ {ρ R > T }} ≤ E I{ρ R >T }
" sup |x(t ∧ ρ R ) − z(t ∧ ρ R )|2 ]
−τ ≤t≤T
≤ E sup |x(t ∧ ρ R ) − z(t ∧ ρ R )|2 ].
(3.19)
−τ ≤t≤T
By Lemma 3.4, for sufficiently small , we may obtain P{ ∩ {ρ R¯ > T }} ≤ /3. Therefore we have P() ≤ , which results in the first limit. The proof of the second limit is similar by Lemma 3.5.
123
Numerical approximation of stochastic differential...
4 Numerical simulation In the section, we shall discuss a highly nonlinear SDDE and its numerical simulations to illustrate convergence in probability established in the previous sections. The computer simulation confirms the efficiency of our numerical method. Example 4.1 Consider a nonlinear scalar SDDE dx(t) = [−ax(t) + bx 2 (t)x 3 (t − τ ) + cx 3 (t) − d x 5 (t)]dt +(ex(t) + f x 2 (t)x(t − τ ))dw(t),
(4.1)
where w(t) is a scalar Brownian motion. For a given step size h, the corresponding EM approximation z k is defined by 3 z k+1 = z k + (−az k + bz k2 z k−m + cz k3 − dz k5 )h + (ez k + f z k2 z k−m ))wk .
(4.2)
For each k > 0, define z¯ (t) = z k , z¯ (t − τ ) = z k−m¯ for t ∈ [kh, (k + 1)h), with the initial value z¯ 0 = ξ0 on [−τ, 0]. That is, z¯ (t) =
∞
z k I[kh,(k+1)h) (t).
k=−m¯
We choose Eq. (4.1) with highly nonlinear coefficients in order to test the results established in the previous sections to be very general and cover a wider class of nonlinear SDDEs. Denote by f (x, y) = −ax + bx 2 y 3 + cx 3 − d x 5 , g(x, y) = ex + f x 2 y. By using p q a p+q + p+q b p+q , it is easy to compute the inequality a p bq ≤ p+q x, f (x, y) ≤ −ax 2 + cx 4 + bx 3 y 3 − d x 6 |b| 6 (x + y 6 ) − d x 6 . ≤ −ax 2 + cx 4 + 2 and |g(x, y)|2 = e2 |x|2 + 2e f x 3 y + f 2 x 4 y 2 f2 3 4 1 4 2 2 x + y + (4x 6 + 2y 6 ). ≤ e x + 2|e f | 4 4 6 p−1 2 2 If a > p−1 2 e , d > c + |b| + 2 ( f + 2|e f |) for p ≥ 2, then the assumptions in Theorem 2.3 and Theorem 3.6 are satisfied, so Eq. (4.1) has a unique almost surely exponentially stable global solution and the corresponding EM numerical solution z¯ (t) to Eq. (4.2) converges in probability to the true solution to Eq. (4.1) (Fig. 1). We simulate the solutions with parameters a = 0.5, b = 0.2, c = 0.3, d = 1.5, e = 0.1, f = 0.2 and a = 0.6, b = 0.8, c = 0.3, d = 6.5, e = 1.1, f = 1.2 for different
123
123
t
6
8
10
x(t)
0
1
2
3 t
4
5
6
h=1/1024 h=1/512 h=1/256 h=1/128
7
Fig. 1 Le f t Computer simulations of paths x(t) with a = 0.5, b = 0.2, c = 0.3, d = 1.5, e = 0.1, f = 0.2. right with a = 0.6, b = 0.8, c = 0.3, d = 6.5, e = 1.1, f = 1.2
0
4
0
2
0.2
0.2
0
0.4
0.4
0.8
1
1.2
1.4
0.6
h=1/1024 h=1/512 h=1/256 h=1/128
x(t)
0.6
0.8
1
1.2
1.4
S. Zhou, C. Hu
Numerical approximation of stochastic differential...
step sizes. In the numerical experiments, the graphs are drawn with the same initial x(0) = 1, and the level axis represents time variable t. We identify the numerical solutions of EM method with a sufficiently small step size (h = 1/1024) as the “true solution”. We draw the numerical solutions obtained from the EM method to Eq. (4.2) with step size h = 1/128, h = 1/256 and h = 1/512 together with the “true solution”. There one can see that the numerical solutions on the left and right plots are close to the their “true” solutions, which exhibits convergence of numerical solutions. Acknowledgments The author expresses her sincere gratitude to two anonymous referees for their detailed comments and helpful suggestions. The financial support from the National Natural Science Foundation of China (Grant No.11301198;11422110).
References 1. Baker, C.T.H., Buckwar, E.: Exponential stability in p-th mean of solutions and of convergent Eulertype solutions of stochastic delay differential equations. J. Comput. Appl. Math. 184, 404–427 (2005) 2. Buckwar, E.: Introduction to the numerical analysis of stochastic delay differential equations. J. Comput. Appl. Math. 125, 297–307 (2000) 3. Mao, X.: Exponential stability of equidistant Euler-Maruyama approximations of stochastic differential delay equations. J. Comput. Appl. Math. 200, 297–316 (2007) 4. Wang, Z., Zhang, C.: An analysis of stability of Milstein for stochastic differential equation with delay. Comput. Math. Appl. 51, 1445–1452 (2006) 5. Zhou, S., Fang, Z.: Numerical approximation of nonlinear neutral stochastic functional differential equation. J. Appl. Math. Comput. 41, 427–445 (2013) 6. Wu, F., Mao, X.: Almost sure stability of Euler-Maruyama approximations for stochastic delay differential equations. Numer. Math. 115, 681–697 (2010) 7. Li, R., Hou, Y.: Convergence and stability of numerical solutions to SDDEs with Markovian switching. Appl. Math. Comput. 175, 1080–1091 (2006) 8. Liu, M., Cao, W., Fan, Z.: Convergence and stability of the semi-implicit Euler method for a linear stochastic differential delay equation. J. Comput. Appl. Math. 170, 255–268 (2004) 9. Mao, X., Yuan, C.: Stochastic Differential Equations with Markovian Switching. Imperial College Press, London (2006) 10. Rathinasamy, A., Balachandran, K.: Mean square stability of semi-implicit Euler method for linear stochastic differential equations with multiple delays and Markovian switching. Appl. Math. Comput. 206, 968–979 (2006) 11. Zhang, H., Gan, S., Hu, L.: The split-step backward Euler method for linear stochastic delay differential equations. J. Comput. Appl. Math. 225, 558–568 (2009) 12. Zhou, S., Wu, F.: Convergence of numerical solutions to neutral stochastic delay differential equation with Markovian switching. J. Comput. Appl. Math. 229, 85–96 (2009) 13. Wang, X., Gan, S.: The improved split-step backward Euler method for stochastic differential delay equations. Int. J. Comput. Math. 88, 2359–2378 (2011) 14. Tretyakov, M.V., Zhang, Z.: A fundamental mean-square convergence theorem for SDEs with locally Lipschitz coefficients and its applications. SIAM J. Numer. Anal. 51, 3135–3162 (2013) 15. Mao, X., Szpruch, L.: Strong convergence and stability of implicit numerical methods for stochastic differential equations with non-globally Lipschitz continuous coefficients. J. Comput. Appl. Math. 238, 14–28 (2013) 16. Szpruch, L., Mao, X., Higham, D.J., Pan, J.: Numerical simulation of a strongly nonlinear ait-sahaliatype interest rate model. BIT Numer. Math. 51, 405–425 (2011) 17. Mao, X., Szpruch, L.: Strong convergence rates for backward Euler-Maruyama method for non-linear dissipative-type stochastic differential equations with super-linear diffusion coefficients. Stoch. Int. J. Prob. Stoch. Proc. 85, 144–171 (2013) 18. Zhou, S.: Strong convergence and stability of backward Euler Maruyama scheme for highly nonlinear hybrid stochastic differential delay equation, Calcolo, doi:10.1007/s10092-014-0124-x
123
S. Zhou, C. Hu 19. Higham, D., Mao, X., Stuart, A.: Strong convergence of Euler-type methods for nonlinear stochastic differential equations. SIAM J. Numer. Anal. 40, 1041–1063 (2002) 20. Hutzenthaler, M., Jentzen, A., Kloeden, P.E.: Strong and weak divergence in finite time of Euler’s method for stochastic differential equations with nonglobally Lipschitz continuous coefficients. Proc. R. Soc. Lond. Ser. A Math. Phys. Eng. Sci. 467, 1563 (2011) 21. Hutzenthaler, M., Jentzen, A., Kloeden, P.: Strong convergence of an explicit numerical method for SDEs with non-globally Lipschitz continuous coefficients. Ann. Appl. Probab. 22, 1611–1641 (2012) 22. Kumar, C., Sabanis, S.: Strong convergence of Euler approximations of stochastic differential equations with delay under local Lipschitz condition. Stoch. Anal. Appl. 32, 207–228 (2014) 23. Mao, X.: Numerical solutions of stochastic differential delay equations under the generalized Khasminskii-type conditions. Appl. Math. Comput. 217, 5512–5524 (2011) 24. Milo˘sevi´c, M.: Highly nonlinear neutral stochastic differential equations with time-dependent delay and the Euler-Maruyama method. Math. Comput. Model 54, 2235–2251 (2011) 25. Hu, Y., Wu, F.: Stochastic Kolmogorov-type population dynamics with infinite distributed delays. Acta Appl. Math. 110, 1407–1428 (2010) 26. Wu, F., Xu, Y.: Stochastic Lotka-Volterra population dynamics with infinite delay. SIAM J. Appl. Math. 70, 641–657 (2009) 27. Mao, X.: Exponential Stability of Stochastic Differential Equation. Dekker, NewYork (1994) 28. Mao, X.: Stochastic Differential Equations and Its Application. Horwood, Chichester (1997)
123