## Abstract

The stochastic Euler scheme is known to converge to the exact solution of a stochastic differential equation (SDE) with globally Lipschitz continuous drift and diffusion coefficients. Recent results extend this convergence to coefficients that grow, at most, linearly. For superlinearly growing coefficients, finite-time convergence in the strong mean-square sense remains. In this article, we answer this question to the negative and prove, for a large class of SDEs with non-globally Lipschitz continuous coefficients, that Euler’s approximation converges neither in the strong mean-square sense nor in the numerically weak sense to the exact solution at a finite time point. Even worse, the difference of the exact solution and of the numerical approximation at a finite time point diverges to infinity in the strong mean-square sense and in the numerically weak sense.

## 1. Introduction

An important numerical scheme for simulating stochastic differential equations (SDEs) is Euler’s method (e.g. Kloeden & Platen 1992; Milstein 1995; Higham 2001). If the coefficients of an SDE are globally Lipschitz continuous, then standard results (e.g. ch. 10 and ch. 14 in Kloeden & Platen 1992) show con-vergence of the Euler approximation in the strong and numerically weak sense to the exact solution of the SDE. It remained an open question whether the Euler approximation also converges in the strong or numerically weak sense at a finite time point if the coefficients of the SDE are not globally Lipschitz continuous; see §1 in Higham *et al.* (2002) for a detailed description of this open problem. In this paper, we answer this question to the negative. More precisely, we prove, for a large class of SDEs with superlinearly growing coefficient functions, that both the distance in the strong *L*^{p}-sense and the distance between the *p*th absolute moments of the Euler approximation and of the exact solution of the SDE diverge to infinity for all . Thus, the Euler scheme does not produce an approximation in the strong or numerically weak sense of the exact solution of such an SDE.

For clarity of exposition, we concentrate in this section on the following prominent example. Let be the unique solution process of the one-dimensional SDE
1.1
for *t*≥0 where is a one-dimensional standard Brownian motion and where is a constant. The Euler approximation of , fixed, is defined recursively through and
1.2
for all *n*∈{0,1,2,…} and all .

Two results motivated us to try to prove convergence of the stochastic Euler approximation to the exact solution *X*_{T} of the SDE (1.1) in the strong mean-square sense as the number of time steps *N* goes to infinity. Gyöngy (1998) established pathwise convergence for SDEs with locally Lipschitz continuous coefficients. More precisely, theorem 1 in Gyöngy (1998) implies
1.3
almost surely pathwise convergence of the stochastic Euler approximation (1.2) to the exact solution of the SDE (1.1). This implies convergence of expectations of continuous and bounded functionals of the difference between the Euler approximation and the exact solution. Of course, the squared difference needed for mean-square convergence is not a bounded, but is an unbounded, functional. Another motivation was an instructive conditional result of Higham *et al.* (2002). They assume, in their theorem 2.2, local Lipschitz continuity of the coefficients of the SDE and boundedness of the *p*th moment of the Euler approximation and of the exact solution in the sense that
1.4
for one arbitrary . Under these assumptions, they establish in theorem 2.2 of Higham *et al.* (2002) strong mean-square convergence of Euler’s method to the exact solution of the SDE. Under assumption (1.4), they, in particular, establish strong mean-square convergence
1.5
of the Euler approximation (1.2) to the exact solution of the SDE (1.1). Boundedness of moments (1.4) of the Euler approximation, however, remained an open question. Higham *et al.* (2002, p. 1060) state that in ‘general, it is not clear when such moment bounds can be expected to hold for explicit methods with *f*,*g*∈*C*^{1}’(*f* and *g* denote in Higham *et al.* (2002) the drift function and the diffusion function, respectively).

In this article, we answer Higham *et al*.’s question in the case of the explicit Euler method and superlinearly growing coefficients of the SDE to the negative (see theorem 2.1 below for details) and establish strong *L*^{p}-divergence
1.6
of the stochastic Euler approximation (1.2) in finite time , for all (see equation (2.6) for details). In addition, numerically weak convergence (e.g. §9.4 in Kloeden & Platen 1992; not to be confused with stochastic weak convergence) fails to hold,
1.7
for all , see equation (2.6). In particular, theorem 2.1 implies that the absolute moments of the Euler approximation (1.2) at a finite time point diverge to infinity, i.e.
1.8
for all . Thus, the moment-bound assumption in theorem 2.2 of Higham *et al.* (2002) (see also inequality (1.4) here) is not satisfied for the SDE (1.1).

Note that the strong divergence (1.6) and the weak divergence (1.7) of the Euler approximation is not a special property of equation (1.1). We establish this divergence for a large class of SDEs with superlinearly growing coefficients in §2. Moreover, our estimates are easily adapted to prove divergence of other numerical schemes such as the Milstein scheme. For simplicity, we restrict ourselves to the Euler scheme. Presence of noise, however, is essential. In the deterministic case, the Euler scheme does converge, and both equations (1.6) and (1.7) fail to hold.

Next, we relate the divergence result (1.8) regarding finite time intervals to a divergence result of Mattingly *et al.* (2002) regarding infinite time intervals (see also Roberts & Tweedie 1996, theorem 3.2; Higham *et al.* 2003, lemma 4.1). Their lemma 6.3 shows, for the SDE (1.1), that the second moment of the Euler approximation diverges in infinite time for any fixed discretization step size, i.e.
1.9
for every fixed . Note that this divergence result is rather different to our divergence result (1.8). First, the time discretization step size does not converge to 0 in equation (1.9). Secondly, the Euler approximation diverges even pathwise with positive probability on an infinite time interval. More precisely, lemma 6.3 in Mattingly *et al.* (2002) proves
1.10

for every fixed . The divergence (1.9) is then an immediate consequence of this pathwise divergence result. On a finite time interval, the Euler approximation does converge pathwise to the exact solution according to Gyöngy (1998) (see equation (1.3) here). So here, the divergence (1.8) of the moments is not a consequence of the pathwise behaviour.

In the next step, we compare the divergence results (1.6)–(1.8) in this article with a few further related results in the literature. It is a classical result (e.g. Kloeden & Platen 1992; Milstein 1995) that the Euler approximation converges in the cases of globally Lipschitz continuous drift and diffusion coefficient functions. In contrast to our results on superlinearly growing coefficients, Yan (2002) proves numerically weak convergence of the Euler scheme if both drift and diffusion functions have, at most, linear growth. Zhang (2006) and Berkaoui *et al.* (2008) prove strong convergence for a class of drift and diffusion functions with a singularity. Yuan & Mao (2008) obtain the rate of strong *L*^{2}-convergence of the stochastic Euler scheme for locally Lipschitz continuous coefficients if these grow, at most, linearly. An explicit bound for the strong *L*^{1}-error is established in Higham & Mao (2005) for the mean-reverting square-root process (also known as the Cox–Ingersoll–Ross process), which has linearly bounded coefficients. Bernard & Fleury (2001) establish convergence in probability under weak hypotheses such as local Lipschitz continuity. Pathwise convergence results can be found in Gyöngy (1998), Fleury (2005), Jentzen & Kloeden (2009) or in Jentzen *et al.* (2009). A number of authors obtain convergence for modified Euler schemes. Milstein & Tretyakov (2005) consider a modified Euler scheme for non-globally Lipschitz continuous coefficients. They obtain numerically weak convergence by discarding trajectories of the approximation, which leave a sufficiently large sphere. Lamba *et al.* (2007) prove strong convergence of an Euler scheme with an adaptive time-stepping algorithm for locally Lipschitz continuous coefficients. A completely different adaptive time-stepping algorithm with a focus on long-time approximation is proposed by Lemaire (2007). Finally, note that in contrast to the explicit Euler method (1.2), the implicit Euler method converges in the root mean-square sense according to Higham *et al.* (2002) (see also Hu 1996; Talay 2002; Szpruch & Mao 2010 and references therein for more convergence results on implicit numerical methods for SDEs).

The rest of this article is organized as follows. Our main result (theorem 2.1) and several examples for the divergence of the Euler scheme are provided in §2. Simulations in §3 illustrate this divergence. The proof of theorem 2.1 is given in §4.

## 2. Main result and examples

Throughout this section, assume that the following setting is fulfilled. Fix and let be a probability space with normal filtration . Additionally, let be a one-dimensional standard -Brownian motion and let be a /-measurable mapping. Moreover, let be two -measurable functions, such that the SDE
2.1
for *t*∈[0,*T*] has a solution. More precisely, we assume the existence of a predictable stochastic process satisfying
2.2
and
2.3
for all *t*∈[0,*T*]. The drift function *μ*(⋅) is the infinitesimal mean of the process *X* and the diffusion function *σ*(⋅) is the infinitesimal standard deviation of the process *X*. The Euler approximation for the SDE (2.1)—denoted by -measurable mappings , *n*∈{0,1,…,*N*}, —is given by and
2.4
for all *n*∈{0,1,…,*N*−1}, and all *ω*∈*Ω*. Now we formalize the main result of this article that asserts the divergence of Euler’s method for the SDE (2.1) if at least one coefficient grows superlinearly.

### Theorem 2.1

*Assume that the setting above is fulfilled with* *and let C≥1, β>α>1 be constants such that
*
2.5
*for all |x|≥C. Then, there exists a constant* *and a sequence of non-empty events* *,* *, with* *and* *for all ω∈Ω*_{N} *and all* *. Moreover, if the exact solution* *of the SDE (*2.1*) satisfies* *for one* *, then
*
2.6

Roughly speaking, the theorem asserts that in the presence of noise, there exists a sequence of events of at least exponentially small probability on which the Euler approximations grow at least double-exponentially fast. Consequently, as being *double-exponentially large* over-compensates that the events have at least *exponentially small* probability, the *L*^{1}-norm of the Euler approximations are unbounded in . The proof of the double-exponential growth of the Euler approximations is based on elementary calculations. For details, the reader is referred to §4 where the proof of theorem 2.1 can be found.

Condition (2.5) should be read as follows. Either the drift function grows in a higher polynomial order than linearly and the diffusion function grows slower than that, or the diffusion function grows in a higher polynomial order than linearly and the drift function grows slower than that. More formally, it suffices to show that either
2.7
for all |*x*|≥*C* or
2.8
for all |*x*|≥*C* and some constants *β*>1, *β*>*α*≥0, *C*>0. Note that our estimates need *β*>1. A drift function of the form for all is too small for our estimates. The assumption that the diffusion function does not vanish on the starting point ensures the presence of noise in the first time step.

In the remainder of this section, we apply theorem 2.1 to a selection of examples. Note that the coefficients in the following examples satisfy an appropriate one-sided linear-growth condition. Therefore, the *p*th absolute moment of the exact solution is finite in each example for every according to theorem 2.4.1 in Mao (1997). For all examples, we check that assumption (2.7) is satisfied. Thus, both the distance in the strong *L*^{p}-sense and the distance between the *p*th absolute moments of the Euler approximation and of the exact solution diverges to infinity for every in each of the following examples.

### (a) The introductory example

The example in §1 is
2.9
for *t*∈[0,*T*]. The dominating coefficient is the drift function with dominating exponent *β*=3. The diffusion function has exponent 0, and we may choose *α*=0. The constant *C* in condition (2.7) can be chosen as *C*=1.

### (b) The stochastic Ginzburg–Landau equation

The Ginzburg–Landau equation is from the theory of superconductivity. It has been introduced by Ginzburg & Landau (1950) to describe a phase transition. Its stochastic version with multiplicative noise can be written as
2.10
for *t*∈[0,*T*] where *η*≥0, . Its solution is known explicitly (e.g. §4.4 in Kloeden & Platen 1992),
2.11
for *t*∈[0,*T*]. Since
2.12
for all , the constants in condition (2.7) can be chosen as *β*=3, *α*=1, .

### (c) The stochastic Verhulst equation

The Verhulst equation is an ordinary differential equation and is a simple model for a population with competition between individuals. Its stochastic version with multiplicative noise can be written as
2.13
for *t*∈[0,*T*], where . Its solution is known explicitly (e.g. §4.4 in Kloeden & Platen 1992) and is given by
2.14
for *t*∈[0,*T*]. The dominant exponent of the drift function is 2 and of the diffusion function is 1. Thus, we may choose *β*=2 and *α*=1. The constant *C* in condition (2.7) can be chosen as .

### (d) Feller diffusion with logistic growth

The branching process with logistic growth (e.g. Lambert 2005) is a stochastic Verhulst equation with Feller noise. It solves
2.15
for *t*∈[0,*T*], where . There is no explicit solution for this equation. However, it features the following self-duality:
2.16
for all *t*∈[0,*T*] and all , where refers to the starting point *X*_{0}=*x* (Hutzenthaler & Wakolbinger 2007). As constants for condition (2.7) serve *β*=2, *α*=1, .

### (e) Protein kinetics

The proportion *x* of one form of a certain protein can be modelled by an ordinary differential equation whose appropriate stochastic version is given by
for *t*∈[0,*T*], where *X*_{0}=*x*_{0}∈(0,1), *η*∈[0,1] and , see §7.1 in Kloeden & Platen (1992). Here, theorem 2.1 applies with *β*=3, *α*=2 and *C* sufficiently large.

## 3. Simulations

In this section, we present two numerical simulations that illustrate the divergence of Euler’s method as formulated in theorem 2.1. For this purpose, we choose an equation with an explicit solution to compare with. Consider the SDE
3.1
for *t*∈[0,3], where is a constant. This equation agrees with the stochastic Ginzburg–Landau equation (2.10) in the case *η*=0,*λ*=1,*x*_{0}=1. Its explicit solution is given by equation (2.11) with *η*=0,*λ*=1,*x*_{0}=1.

Table 1 shows Monte Carlo simulations of the second moment both of the exact solution and of the Euler approximation for five values of . For each value of , the table contains 10 simulation runs for the Euler approximation. The number of time steps is *N*=10^{3}.

In the case , the parameters of the simulation are such that the Euler approximation produces the value ‘NaN’ (NaN is the Institute of Electrical and Electronic Engineers (IEEE) arithmetic representation for ‘not-a-number’, here this is because of an operation ‘Inf-Inf’ where Inf is the IEEE arithmetic representation for positive infinity). In contrast to this, our Monte Carlo simulations in the cases are all finite. In these cases, the probability of the event on which the Euler approximation diverges seems to be rather small. The last but one column in table 1 exemplifies a parameter setting in which the event of large growth has a probability such that in 10^{5} runs, in some Monte Carlo simulations, no explosion occurs and in some Monte Carlo simulations, at least one explosion occurs.

The values in table 1 are either within distance 2 of the true value or NaN. This is owing to the double-exponential growth of the deterministic system for some initial values. If the simulation starts to grow, then it reaches NaN very quickly. We encountered similar behaviour for other exponents greater than 2. In order to see double-exponential growth of the Euler approximation in a plot, we consider an exponent close to 1. More precisely, we plot the Monte Carlo simulation of the first absolute moment of the Euler approximation of *X*_{10}, where
3.2
for *t* in [0,10], see figure 1 (see also the simulation values in table 2). Note that the graph resembles an exponential function and that the *y*-axis is logarithmic. Thus, the growth of the graph in figure 1 is indeed close to a double-exponential function. In addition, we should mention that some fine-tuning was necessary to obtain suitable parameters for which the simulated absolute moment grows but is not NaN.

## 4. Proofs

### Lemma 4.1

*Let* *be a probability space and let* *be an* -*measurable mapping that is standard normal distributed. Then*,
4.1
*for all* .

### Proof of lemma 4.1.

We have 4.2 for all . Moreover, a similar estimate yields 4.3 for all . ■

### Proof of theorem 2.1.

By assumption, we have . Therefore, there exists a real number , such that
4.4
Then, we define
4.5
and consider the sets given by
4.6
for all . Note that
4.7
for all due to the definition of *r*_{N}, .

Let now and *ω*∈*Ω*_{N} be arbitrary. Then, we claim
4.8
for all *n*∈{1,2,…,*N*}. We prove equation (4.8) by induction on *n*∈{1,2,…,*N*}. In the base case *n*=1, we have
4.9
due to the definition (2.4) of and due to the definition (4.6) of *Ω*_{N}. For the induction step *n*→*n*+1, we assume that equation (4.8) holds for one *n*∈{1,2,…,*N*−1}. In particular, we then obtain
4.10
since *r*_{N}≥2. Additionally, we have
and
4.11
due to definitions (2.4) and (4.6). Therefore, the growth condition (2.5) and inequality (4.10) imply
4.12
where the last step follows from inequality (4.7). The induction hypothesis hence yields
4.13
This proves inequality (4.8) for *n*+1. Therefore, equation (4.8) indeed holds for all *n*∈{1,2,…,*N*}. In particular—as and *ω*∈*Ω*_{N} were arbitrary—we obtain
4.14
for all *ω*∈*Ω*_{N} and all .

Furthermore, definition (4.4) gives

for all . Therefore, lemma 4.1 yields 4.15 for all . This shows the existence of a constant such that 4.16 for all .

Combining equations (4.14) and (4.16) then gives 4.17 and Jensen’s inequality hence shows 4.18 Since by assumption, we finally conclude 4.19 Equations (4.14) and (4.16) finally complete the proof of theorem 2.1. ■

## Acknowledgements

This work has been partially supported by the Collaborative Research Centre 701 ‘Spectral Structures and Topological Methods in Mathematics’ funded by the German Research Foundation. Moreover, the authors thank three anonymous referees for their helpful comments.

- Received July 6, 2010.
- Accepted November 16, 2010.

- This journal is © 2010 The Royal Society