Adaptive Euler methods for stochastic systems with non-globally Lipschitz coefficients

We present strongly convergent explicit and semi-implicit adaptive numerical schemes for systems of semi-linear stochastic differential equations (SDEs) where both the drift and diffusion are not globally Lipschitz continuous. Numerical instability may arise either from the stiffness of the linear operator or from the perturbation of the nonlinear drift under discretization, or both. Typical applications arise from the space discretization of an SPDE, stochastic volatility models in finance, or certain ecological models. Under conditions that include montonicity, we prove that a timestepping strategy which adapts the stepsize based on the drift alone is sufficient to control growth and to obtain strong convergence with polynomial order. The order of strong convergence of our scheme is (1 − ε)/2, for ε ∈ (0,1), where ε becomes arbitrarily small as the number of finite moments available for solutions of the SDE increases. Numerically, we compare the adaptive semi-implicit method to a fully drift-implicit method and to three other explicit methods. Our numerical results show that overall the adaptive semi-implicit method is robust, efficient, and well suited as a general purpose solver.


Introduction
Consider the d-dimensional semi-linear stochastic differential equation (SDE) of Itô type (1) dX(t) = [AX(t) + f (X(t)]dt + g(X(t))dW (t), t ∈ [0, T ]; X(0) ∈ R d , where T > 0, A ∈ R d×d , f : R d → R d , g : R d → R d×m , and W is an m-dimensional Wiener process.We suppose that the drift coefficient f and the diffusion coefficient g together satisfy polynomial bounds and a monotone condition permitting g to grow superlinearly as long as that growth is countered sufficiently strongly by f .Global Lipschitz bounds are not required.For example, consider f (x) = −x 2 with g(x) = x 3/2 or f (x) = −x 5 with g(x) = x 2 .Such applications arise in finance: for example the Lewis stochastic volatility model [17] which has a polynomial diffusion coefficient of order 3/2.It was shown in [11] that the explicit Euler-Maruyama method with constant stepsize fails to converge in the strong sense to solutions of (1) if either the drift or the diffusion coefficients grow superlinearly.Also, as noted in [4], fixed stepsize schemes may need to use very small stepsizes when the SDE being solved is stiff.We address these issues here by a semi-implicit scheme with adaptive timestepping.
In [14] a class of timestepping strategies, referred to as admissible, was motivated for the numerical discretisation of SDEs where the drift satisfies a one-sided Lipschitz coefficient and the diffusion satisfies a global Lipschitz bound.An admissible strategy uses the present value of the numerical trajectory to select the next timestep to avoid spuriously large drift responses.This is distinct from the error control approach in (for example) [4,5,13].
Timesteps selected by an admissible strategy are subject to upper and lower limits h max and h min in a fixed ratio ρ, with h max serving as a convergence parameter and h min serving to ensure that the simulation completes in a reasonable time.If the strategy attempts to select a timestep smaller than h min , then a backstop method is applied instead over a single step of length h min .It was proved in [14] that the explicit Euler-Maruyama method over a random mesh generated by an admissible timestepping strategy is strongly convergent in h max with order 1/2.The proof relied upon p th -moment bounds on the supremum of solutions of the underlying SDE.Note also the adaptive approach in [6] which is consistent with the admissibility condition of [14].
Here, we examine more general SDEs and consider simultaneously both explicit and semiimplicit Euler-Maruyama schemes.Due to the monotone condition on the drift and diffusion terms, our analysis must contend with only a finite number of available bounded SDE moments (see for example the estimates provided by parts (i) and (ii) of Lemma 15).Unlike in [14], we characterize precisely the backstop scheme and integrate it into the analysis in a way that is compatible with taking a random number of timesteps.In this way we show that a class of admissible timestepping strategies depending only on the drift response can be used to ensure that both the explicit and semi-implicit adaptive Euler-Maruyama schemes are strongly convergent to solutions of (1) with order (1 − ε)/2, in the sense that for any ε ∈ (0, 1), there exists C ε > 0, independent of h max such that where Y N is value of the numerical scheme at time T , and • is the l 2 norm.The reduction in the order of strong convergence in our main result (when compared to that in [14]) is a direct consequence of the loss of global Lipschitz continuity in the diffusion coefficient.If we reimpose global Lipschitz continuity on the diffusion, we recover a strong convergence order of 1/2, and if we decompose the drift of (1) so that A = 0, we recover the main result of [14]: see Remark 22 for more details of this.The nature of the monotone condition is such that a timestepping scheme which is admissible, and which can therefore successfully control the drift response, will also be sufficient to control the diffusion response.It is well documented that the structure of the drift function (both linear and nonlinear) under discretisation may have local dynamics that render the stability of equilibria vulnerable to the effects of perturbation, either stochastic or numerical [1,3,7,8,11].
Our method handles stiffness leading to potential instability in the discretisation in two distinct ways.Where there is a classic (deterministic) stiff linear system we are able to treat this term implicitly without sacrificing numerical efficiency.Adaptive timestepping then treats nonlinear stability under stochastic perturbation.Thus, we deal with each source of potential instability separately, as would a stochastic IMEX-type method.The use of an implicit approach to deal with the linear part of the drift avoids any consideration of potential interactions between it and the diffusion or between it and the nonlinear part of the drift.Note that the decomposition of the drift into the form AX(t) + f (X(t)) is determined by the modeller, and when A = 0 the convergence analysis in this article applies equally to a fully explicit method if desired.The literature already contains numerical schemes with fixed stepsizes that converge strongly to solutions of SDEs with coefficients that satisfy local Lipschitz and monotone conditions.Several of these extend the idea of taming as introduced in [12], which rescales the functional response of the drift coefficient in the scheme; they do so by allowing the entire stochastic Euler map to be rescaled by some combination of drift and diffusion responses.For example, see the balanced method introduced in [30] and the variant presented in [25], which are both strongly convergent in this setting.The projected Euler method of [2] handles runaway trajectories by projecting them back onto a ball of radius inversely proportional to the step size; hence the authors control moments of the numerical solution.It was shown in [24] that a drift-implicit discretisation could also ensure strong convergence in our setting.Finally we highlight [10], which treats SDEs and SPDEs with non-globally monotone coefficients.
In Section 5, we compare the numerical performance of a selection of these methods to that of the adaptive scheme presented in this article.We note this selection cannot be exhaustive and there are a growing number of variations; see for example [9,22,23,26,29,31].However our examples illustrate some of the drawbacks of fixed-step explicit schemes (when linear stability is an issue) and where for fixed, relatively large h, the taming perturbation which imposes convergence may change the dynamics of the solution.Compared to the fixed-step explicit methods, our numerical results show that the semi-implicit adaptive method gives consistently reliable numerical convergence results.It is also more efficient than the drift-implicit scheme for SODEs, though this comparison is less clear for the SPDE example.
The structure of the article is as follows.In Section 2 we describe the monotone condition and polynomial bounds that must be satisfied by f and g, and provide the p th -moment bounds satisfied by the solutions of (1) within that framework.In Section 3, we introduce the semiimplicit Euler-Maruyama method that, applied stepwise over a random mesh and combined with an appropriate backstop method, is the focus of the article.A mathematical definition for meshes produced by admissible timestepping strategies is provided, and conditional moment bounds for the SDE solution associated with these meshes are derived.In Section 4, we present our main convergence result and state several technical lemmas, with proofs provided in Section 6.In Section 5, we carry out a comparative numerical investigation of strongly convergent schemes from the selection discussed above.

Setting
Throughout the paper, N denotes the set {0, 1, 2, . ..}, • denotes the l 2 norm of a ddimensional vector, • F the Frobenius norm of a d × m-dimensional matrix, and for any x ∈ R d and i = 1, . . ., m, g i (x) denotes the i th column of the diffusion coefficient matrix g(x).For a, b ∈ R let a ∨ b denote max{a, b}.For any A ∈ R d×d , we let A 1/2 ∈ C d×d denote the matrix such that (A 1/2 ) 2 = A. Let (F t ) t≥0 be the natural filtration of W .To ensure the existence of a unique strong solution for (1) (in the sense of [21, Definition 2.2.1]) over the interval [0, T ], it suffices to place local Lipschitz and monotone conditions on f and g: for x, y ∈ R d with x ∨ y ≤ R, and there exists c ≥ 0 such that for some p ≥ 0 We also require a set of polynomial bounds on the derivatives of f and g, and hence on f and g themselves.The minimum value of p in (2) required to prove our main strong convergence result depends on the order of these bounds.Assumption 2. Suppose f : R d → R d and g : R d → R d×m are continuously differentiable with derivatives bounded as follows: for some c j , γ 0 , γ 1 ≥ 0; j = 1, . . ., 4, we have where Df (x) ∈ R d×d is the matrix of partial derivatives of f , and Dg i (x) ∈ R d×d is the matrix of partial derivatives of the i th column of g, and We require that some of the moments of the solutions of (1) are bounded over the interval [0, T ].Further, (2) in Assumption 1 implies (see, for example, Tretyakov & Zhang [30]) that there exists c ′ ≥ 0 such that (5) x, f (x This is a special case of Khasminskii's condition using the Lyapunov-type function V (x) = 1 + x 2 , and it guarantees the existence of a unique strong solution of (1) over [0, T ] for any T < ∞ (see [21,Theorem 2.3.5]),while also ensuring p th -moment bounds as follows: Lemma 3. Let (X(t)) t∈[0,T ] be the unique solution of (1).Suppose that (5) holds for some p ≥ 2 and (4) in Assumption 2 holds, then there exists M p,T < ∞ such that Proof.The proof of (6) follows from [23,Lemma 4.2], since the bound on g provided by (4) implies Eq. (4.2) in that article, which we reproduce here as To ensure sufficiently many bounded moments of the form (6) for our analysis to work, we now impose a lower bound on the value of p in (2) that depends on the order of the polynomial bounds on f and g.This bound is associated with a tolerance parameter ε which then appears in the the rate of strong convergence in Theorem 20.
Finally, note that the analysis in this article is also valid if the initial vector is random, F 0 -measurable, and E X(0) p < ∞.

An adaptive semi-implicit Euler scheme with backstop
The adaptive timestepping scheme under investigation in this article is based upon the semiimplicit Euler-Maruyama scheme over a random mesh {t n } n∈N on the interval [0, T ] given by ( 8) where △W n+1 := W (t n+1 ) − W (t n ), and {h n } n∈N is a sequence of positive random timesteps and {t n := n i=1 h i−1 } n∈N\{0} with t 0 = 0.For the setting described in Section 2, we show that, in order to ensure strong convergence with order (1 − ε)/2 of the method (8) for any ε ∈ (0, 1), it is sufficient to construct the stepsize sequence {h n } n∈N in the same way as in [14], demonstrating the applicability of this strategy to a significantly broader class of SDEs.We review the construction now.Definition 5. Suppose that each member of {t n := n i=1 h i−1 } n∈N\{0} , with t 0 = 0, is an F tstopping time: i.e. {t n ≤ t} ∈ F t for all t ≥ 0, where (F t ) t≥0 is the natural filtration of W .The filtration (F t ) t≥0 can be extended (see [21]) to any F t -stopping time τ by In particular this allows us to condition on F tn at any point on the random time-set {t n } n∈N .Remark 6.Throughout the article, the index of a random sequence reflects its F tn -measurability.For example, consider the timestep sequence {h n } n∈N : each h n is F tn -measurable.The only exception to this is {t n } n∈N , since each t n is F tn−1 -measurable.
Assumption 7. Suppose that each h n is F tn -measurable, and that there are constant minimum and maximum stepsizes h min and h max imposed in a fixed ratio ρ so that Definition 8.For each t ∈ [0, T ], define the random integer N (t) such that Set N := N (T ) and t N := T , so that T is always the last point on the mesh.
We note that N (t) is a.e. the index of the right endpoint of the step that contains t, i.e. t ∈ [t N (t) −1 , t N (t) ] a.s.Both t N (t) and t N (t) −1 are F t -stopping times, and N (t) is supported on the finite set {N Remark 9.In (8), note that each Therefore △W n+1 is a function of values of the Wiener process up to time t n , is not independent of F tn , and there is no reason to expect that △W n+1 ∼ N (0, h n I d×d ), where I d×d is the identity matrix.However since t n+1 is a bounded F tn -stopping time then, by Doob's optional sampling theorem (see, for example, [27]), In our main analysis, we use the following lemma on the boundedness of the moments of solutions of (1) conditioned at points on our adaptive mesh.The proof is a modification of that of [21,Theorem 2.4.1].
Lemma 10.Let (X(t)) t∈[0,T ] be a unique solution of (1), and suppose that (5) holds for some p ≥ 2. Then there exist constants ν 1 and ν 2 such that We are now in a position to define the scheme which is the subject of this article, and which combines a semi-implicit Euler scheme over an adaptive mesh, generated according to an admissible timestepping strategy, with a backstop method.
Then we define a semi-implicit Euler scheme with backstop as the sequence { Y n } n∈N by where {h n } n∈N satisfies the conditions of Assumption 7. The map ϕ : characterises the form a user-selected backstop method.We require that for some positive constants C 1 and C 2 , and ε ∈ (0, 1) the fixed parameter from Assumption 4.
Remark 12.Note that ϕ will satisfy (13) if the backstop method is subject to a mean-square consistency requirement that bounds the propagation of discretisation error over a single step.In practice, rather than checking (13) directly, we use as our backstop a method that is known to be strongly convergent of order 1/2 in this setting: for the numerical experiments in Section 5 we use the balanced method introduced by Tretyakov & Zhang [30], which satisfies a similar local accuracy bound (see [30, Eq. (2.9)]) and corresponds to the map Examples of h min , h max , and ρ for the scheme (12) with backstop characterised by (14) are given in Section 5.
Definition 13.Let { Y n } n∈N be a solution of (12) where f and g satisfy the conditions of Assumptions 1, 2, and 4. We say that {h n } n∈N is an admissible timestepping strategy for (12) if Assumption 7 is satisfied and there exist real non-negative constants For example (see [14]) the timestepping rule given by is admissible for (12).Choosing the larger of 1/ f ( Y n ) and Y n / f ( Y n ) helps maximize the stepsize while maintaining its admissibility.The backstop is needed since it may not always be possible to control Y n via timestep so that (15) holds.See Section 7 for a more detailed comment.

4.
Strong convergence of the adaptive scheme 4.1.Preliminary lemmas.These lemmas provide a regularity bound in time and an estimate on remainder terms from Taylor expansions of f and g.Proofs are given in Section 6.
Lemma 14.Let (X(t)) t∈[0,T ] be a solution of (1) with coefficients f and g satisfying the conditions of Assumptions 1, 2, and 4, and suppose that the sequence of random times {t n } n∈N is as in Definition 5 and satisfies the conditions of Assumption 7. Then for all n ∈ N and p ≥ 2, there exists an a.s.finite and F tn -measurable random variable Lp,n so that

Now consider the Taylor expansions of f and g
where the remainders R f and R gi are given in integral form by and z can be taken to read either f or g i .Furthermore let We now give F tn -conditional mean-square bounds on the integrals of these remainder terms.
Lemma 15.Let (X(t)) t∈[0,T ] be a solution of (1) with coefficients f and g satisfying the conditions of Assumptions 1, 2, and 4. Let {t n } n∈N be as in Definition 5, satisfying the conditions of Assumption 7.
Then for any ε ∈ (0, 1) there is an a.s.finite and F tn -measurable random variable Λε,n > 0, and a constant Λ ε < ∞, the latter independent of {h n } n∈N and h max , such that Remark 16.The notational convention used in Part (iv) of Lemma 15 is extended throughout the paper to F tn -adapted sequences for which there exists a finite uniform upper bound on the expectation of each term.

Main results.
In this section, we provide a lemma giving a bound on the one-step error for the combined scheme, followed by our main theorem which uses discrete Gronwall inequalities to produce an order of strong convergence.
Lemma 17.Let (X(t)) t∈[0,T ] be a solution of (1) with drift and diffusion coefficients f and g satisfying the conditions of Assumptions respectively, such that where ε ∈ (0, 1) is the fixed parameter from Assumption 4 and constants Q, Γ 1 are given by where c, R 2 , C 2 satisfy (2) in Assumption 1, (13) in Definition 11,and (15) in Definition 13 respectively.
Proof.For h n selected at time t n , for some n ∈ N, by an admissible timestepping strategy, there are two possible cases (denoted (I) and (II)), first, h min < h n ≤ h max and second, h n = h min .We consider each in turn.
(I) In this case we rely on the bound (15) provided by the admissibility of the timestepping scheme.When h min < h n ≤ h max , Y n+1 is derived from Y n using ( 8), and we have Expand f and g as Taylor series around X(t n ) over the interval of integration, and write Therefore which is Let Q be as defined in (18).Then, using that h max ≤ 1 and the inequality 2 x, y ≤ x 2 + y 2 , we find We develop bounds on The terms after C n+1 on the RHS of the inequality have zero conditional expectation, by (11) in Remark 9, and the fact that E n and each Rgi are conditionally independent with respect to F tn , with the latter an Itô integral with zero conditional expectation.
By (2) in Assumption 1, By (15) in Definition 13,and (4) in Assumption 2 we have )), a.s.Next, by Jensen's inequality and part (i) of Lemma 15, We also have from part (iii) of Lemma 15 Applying parts (i)-(iii) of Lemma 15 then gives Here Y n+1 is generated from Y n via an application of the backstop method over a single step of length h min .This corresponds to a single application of the map ϕ and therefore the relation ( 13) is satisfied a.s.
To combine the two cases (I) and (II), define the a.s.finite and F tn -measurable random variables Γ(m) where C 1 and C 2 are as given in (13).Since Q, Λε,n ≥ 0 (the latter in the a.s.sense), we have for any h n selected by an admissible adaptive timestepping strategy, (19) we may apply (6) to (20) to show that, under Assumption 4, where M 2+2γ1,T and M 2(γ0+γ1)+2,T are finite constants satisfying (6) for p = 2, 2γ 0 + 2 respectively.It can be similarly shown under Assumption 4 that there exist finite constants Γ and E Λε,n ≤ Λ ε .
The bound (22) characterises the propagation of error in mean-square over a single step of the combined semi-implicit Euler scheme with backstop (12), and holds regardless of whether or not the timestepping strategy requires an application of the semi-implicit scheme or the backstop scheme.
Assumption 18.Finally, h max is chosen so that there exists a constant δ ∈ [0, 1] that does not depend on h max , such that where Q is defined in (18) and Γ 1 is defined in (19).It follows that there exists γ < ∞ such that Although these conditions are required in the proof of Theorem 20 we have observed no practical implications in our numerical experiments.
The accumulation of error in mean-square for (12), and hence the order of strong convergence, can now be estimated.Theorem 20.Let (X(t)) t∈[0,T ] be a solution of (1) with drift and diffusion coefficients f and g satisfying the conditions of Assumptions 1, 2, and 4. Let { Y n } n∈N be a solution of (12) in Definition 11, with initial value Y 0 = X 0 and admissible timestepping strategy {h n } n∈N satisfying the conditions of Definition 13 and Assumption 18.Then, if ε ∈ (0, 1) is the fixed parameter from Assumption 4, there exists C ε,m,ρ,T > 0, independent of h max such that where E 2 c (t) is as defined in Definition 19, and in particular where N is as given in Definition 8.

Remark 21.
By setting A = 0 in (1) we obtain strong convergence of identical order of a backstopped fully explicit Euler-Maruyama adaptive method.
Since t N (t) is a F t -stopping time, the event {N (t) ≤ n} ∈ F tn and therefore the indicator variable I {N (t) ≥n+1} is F tn -measurable.Thus we have + 126 Λε,n h n 3−ε I {N (t) ≥n+1} , a.s.

A comparative numerical review of some available schemes
Given the semi-linear SDE (36 with solution u : [0, T ] × Ω → R d , we compare our semi-implicit adaptive numerical method to a number of different fixed-step schemes (with time step h) that we outline below.Some numerical examples for an explicit adaptive scheme are given in [14].The majority of recent developments concentrate on a perturbation of the flow (or solution) of order h 1/2 or higher, however the first method we present is the classic implicit approach.We do not consider an exhaustive list of taming-type schemes and there are other variants available, see for example [9,22,23,26,29,31].
Our examples illustrate some of the drawbacks of explicit schemes, for example where linear stability is an issue.1. Drift implicit scheme [24] This is given for (36) by Although strong convergence has been proved (see [24]), at each step a nonlinear system of the form needs to be solved for Y n+1 for some vector b.Even for the deterministic case there is no guarantee the nonlinear solver will converge to the correct root (see [28,Chapter 4]).We observe in our numerical experiments that both a standard Newton method and the matlab nonlinear solver fsolve (or fzero in one-dimension) may fail to converge.In the event of a step where this occurs we use as a backstop an alternative explicit method, in this article taken to be the balanced method (see below).The drift implicit scheme with this backstop method is denoted by Drift Implicit in the figures of this Section.Finally, note that for several examples in this section the implicit solver may be made more efficient by exploiting a known closed-form solution for the nonlinear system (37).Such solutions are not in general available and so we do not make use of them in our comparative analysis here.

Tamed [25]
A tamed version which may be used when the solutions of (36) have a limited number of finite moments [31] Strong convergence of order 1/2 is achieved by setting β = 1/2.We denote this method Tamed.3. Balanced Method [30] is given for (36) by .
This was proved to be strongly convergent with order 1/2 (including for additive noise) and is denoted in the figures of this Section as BM.

4.
Projected EM [2] uses the standard EM method when Y n is inside a ball of radius inversely proportional to the step step size.However, outside of the ball, the numerical solution is projected onto the ball.We have for

This is denoted below as Projected EM.
We provide a comparative illustration of the combined effect of semi-implicitness and adaptivity using five examples ranging from geometric Brownian motion to a system of SDEs arising from the spatial discretisation of an SPDE.Recall that our use of a semi-implicit method controls instabilities from a linear operator and the adaptive timestepping controls the discretisation of the nonlinear structure.Stiffness is manifested in the structure of each of these equations in different ways: ranging from the linearity only (in geometric Brownian motion) to both in the linear operator and nonlinearities for a discretisation of an SPDE.
To examine strong convergence in h max for the SDE examples below we solve with M = 1000 samples to estimate the root mean square error (RMSE) at a final time ] and we estimate the standard deviation from 20 groups of 50 samples included on the error plots.Reference solutions are computed with 10 6 uniform steps on [0, T ].For efficiency we compare the RMSE against the average computing time over the 1000 samples (denoted cputime).Unless otherwise stated we take ρ = 10 throughout.5.1.Geometric Brownian Motion.The classic example to illustrate linear mean square stability is geometric Brownian motion (38) du(t) = ru(t)dt + σu(t)dW (t), u(0) = u 0 , t ≥ 0.
If r + σ 2 /2 < 0 it is straightforward to see that E (u(t) 2 → 0 as t → ∞ and that the (fixed step) explicit Euler method is only stable if 0 < ∆t < −2(r + σ 2 /2)/r 2 .The drift and diffusion are both linear functions, so there is no need for either taming or adaptivity to control growth from a nonlinear term; indeed in this example the semi-implicit adaptive and fully drift implicit schemes co-coincide if A = r and f (u) = 0. However it is instructive to compare the explicit schemes to the implicit schemes (Adaptive IEM and Drift Implicit).We take r = −8 and σ = 3 so that the explicit Euler method is unstable for ∆t = 0.25 and ∆t = 0.5.In Fig. 1 we plot the error squared, |u(t n ) − Y n | 2 , of two sample paths one with h max = 0.25 (a) and h max = 0.5 (b).Although the tamed and projected schemes control growth from the linear instability, we observe that this control can come at a price of bounded oscillations.5.2.1D stochastic Ginzburg-Landau.The 1D stochastic Ginzburg-Landau SDE is a classic example with a cubic nonlinearity in the drift and linear diffusion (39) We take here parameter values as in [22,Example 4.7], a = 0.1, b = 1 and c = 0.2, x(0) = 2, and solve to T = 1.We see in Fig. 2 that all the methods demonstrate convergence and that Adaptive IEM and Projected EM are similar in convergence and efficiency.Neither the adaptive nor drift-implicit schemes used the backstop method.5.3.Stochastic volatility system.We consider an extension of the 3/2-volatility model to two dimensions as in [26]  We see in Fig. 3 that all the methods demonstrate convergence but that BM and Tamed have a larger error constant and the adaptive method Adaptive IEM is the most efficient.The initial data taken was X(0) = [2, 2] T and the backstop method was not used for either drift implicit or adaptive methods (as for (39)).In Fig. 4 we examine the time steps h n for a single noise path with the same value of h max = 0.01 but with ρ = 10 and ρ = 100 corresponding to h min = 10 −3 and h min = 10 −4 .In (a) we have initial data X(0) = [20,20] T and in (b) X(0) = [200, 200] T .In (a) we observe that with ρ = 100 the backstop method is not used at all (but it is required for ρ = 10) whereas in (b), to deal with the larger initial data, the backstop is used in both cases.As time progresses the time-step taken increases until it reaches h max .
These observations suggest that practitioners who apply a standard explicit or semi-implicit Euler-Maruyama scheme over a uniform mesh with a step size sufficiently small (e.g.close to h min with large ρ) may rarely encounter the spurious coefficient responses that underlie the lack of strong convergence for the scheme.
for some m ∈ N \ {0} and where W j (t) are mutually independent standard Brownian motions.We introduce a grid of d + 2 uniformly spaced points x k = k∆x on [0, 1] for k = 1, . . ., d + 2. Then the finite difference approximation in space leads to a system of d SDEs: where we denote u := (u 1 , u 2 , . . ., u d ) T , u k (t) ≈ u(x k , t) and the noise process is  The tri-diagonal matrix A is the standard finite difference approximation to the Laplacian with zero Dirichlet boundary conditions given by For further details on the finite difference approximation of SPDEs see [20].To be able to compare different system sizes d we scaled the Euclidean norm in the standard way (see for example [20]) to approximate the function space norm L 2 ([0, 1]).This system of SDEs displays linear stiffness (similar to the geometric Brownian motion) and nonlinear stiffness arising from the drift and diffusion coefficients.The parameter ǫ then determines the degree of linear stiffness.
To examine convergence to the SDE system (42) we take ǫ = 0.1, T = 1 for d = m = 10 (Fig. 5) and d = m = 100 (Fig. 6).For the smaller system size (d = 10) in Fig. 5 (a) we see all methods converging for h max sufficiently small, although Projected EM initially diverges.We see in (b) that Adaptive IEM is more efficient than Drift Implicit and is similar to the other explicit schemes.When m = 100 Projected EM is no longer seen to converge on this range of h and so is not plotted in Fig. 6.In (a) we now see that the Drift Implicit is more accurate for a given h max and from (b) that it is comparable in efficiency with Adaptive IEM.

Proofs of Technical Results
In this section we frequently use the inequality a + b p ≤ 2 p ( a p + b p ), where a, b ∈ R d , and p ∈ R + , which follows from a  By the triangle inequality and [21, Theorem 1.7.1](with conditioning on F tn ), Next, we apply (4), Lemma 10, and the fact that A p F < ∞, to get a.s.
Therefore, since |s − t n | ≤ 1, we can define an a.s.finite and F tn -measurable random variable so that (16) holds.
Part (iii) holds as a special case of Part (i).Part (iv) follows by an application of the Cauchy-Schwarz inequality, followed by Jensen's inequality for the functions (•) 1/(2 q−1 ) and (•) (both of which are concave over R + , by the second derivative test), to get which is finite under the conditions of Assumption 4: p satisfies (7), and therefore by ( 6) the finiteness of E Ῡq,n is ensured by (44) and that of E L2 n is ensured by (43).
Remark 22.By making q successive applications of the Cauchy-Schwarz inequality in the proof of Lemma 15 we separate the expectation of dependent random factors in R f and R g in such a way that the highest possible order of h n is achieved in the estimate, given the available finite moment bounds.This is necessary to ensure a polynomial order of strong convergence in the statement of Theorem 20.If the diffusion coefficient g is globally Lipschitz continuous then the resulting uniform bound on each Dg i (x) F , along with stronger moment bounds of the form (6), sidesteps that requirement.In this case the statement of Lemma 15, and hence the statement of Theorem 20, would hold with ε = 0 (and order constant independent of q, and therefore ε), giving an order of strong convergence of 1/2 for the semi-implicit method with backstop (12), using an admissible timestepping strategy.If we then set A = 0 in (1), our method becomes explicit and we recover the main result of [14].

Conclusion
The discretisation of SDEs with non-Lipschitz drift and diffusion coefficients is a challenging numerical problem.We have proved strong convergence for both adaptive semi-implicit and explicit Euler schemes, and presented numerical results that indicate the semi-implicit variant is well suited as a general purpose solver, being more robust than several competing explicit fixed-step methods and more efficient than the drift implicit method.
Both the drift implicit and the adaptive scheme make use of a backstop method which is triggered when the adaptive timestepping strategy attempts to select a timestep at the minimum stepsize h min .Our numerical experiments indicate that, for an appropriate choice of ρ, h min may be achieved only rarely (if at all).It may be possible to characterise the probability of this occurrence and, if it can be bounded appropriately, a strong convergence result may be possible for a numerical method of the form (8) that does not rely on a backstop method (provided T is reached in a finite number of steps).A step in this direction may be found in [15].
SDEs where the drift coefficient is both positive and non-globally Lipschitz continuous are not covered by the analysis in this article, though adaptive meshes have been used to reproduce positivity of solutions with high probability and a.s.stability and instability of equilibria in [16] (informed by the approach of Liu & Mao [18]).We are unaware of any strong convergence results for such equations.
Finally, since our analysis relies upon the boundedness of A F , and since the error constant in the strong convergence estimate increases without bound with the number of independent noise terms m, the results of the article do not automatically extend to SPDEs.This setting is now considered in [19].

Figure 2 .Figure 3 .
Figure 2. Convergence and efficiency of methods applied to the stochastic Ginzburg Landau equation (39).We compare RMSE at T = 1 against h max in (a) and efficiency (RMSE vs cputime) in (b).

Figure 5 .
Figure 5. Convergence and efficiency for the methods applied to the finite difference approximation of the SPDE given by (42) with d = 10 (a) RMSE vs h max (with reference line of slope 0.5) and (b) the efficiency (RMSE vs cputime).

Figure 6 .
Figure 6.Convergence and efficiency for the methods applied to the finite difference approximation of the SPDE given by (42) with d = 100 (a) RMSE vs h max (with reference line of slope 0.5) and (b) the efficiency (RMSE vs cputime).