Efficient weak second-order stochastic Runge–Kutta methods for Itô stochastic differential equations

Adaptive time-stepping with high-order embedded Runge-Kutta pairs and rejection sampling provides efficient approaches for solving differential equations. While many such methods exist for solving deterministic systems, little progress has been made for stochastic variants. One challenge in developing adaptive methods for stochastic differential equations (SDEs) is the construction of embedded schemes with direct error estimates. We present a new class of embedded stochastic Runge-Kutta (SRK) methods with strong order 1.5 which have a natural embedding of strong order 1.0 methods. This allows for the derivation of an error estimate which requires no additional function evaluations. Next we derive a general method to reject the time steps without losing information about the future Brownian path termed Rejection Sampling with Memory (RSwM). This method utilizes a stack data structure to do rejection sampling, costing only a few floating point calculations. We show numerically that the methods generate statistically-correct and tolerance-controlled solutions. Lastly, we show that this form of adaptivity can be applied to systems of equations, and demonstrate that it solves a stiff biological model 12.28x faster than common fixed timestep algorithms. Our approach only requires the solution to a bridging problem and thus lends itself to natural generalizations beyond SDEs.

Tamed Runge-Kutta methods for SDEs with super-linearly growing drift and diffusion coefficients

Some efficient stochastic Runge–Kutta (SRK) methods for the strong as well as for the weak approximation of solutions of stochastic differential equations (SDEs) with improved computational complexity are considered. Their convergence is analyzed by a concise colored rooted tree approach for both, Ito as well as Stratonovich SDEs. Further, order conditions for the coefficients of order 1.0 and 1.5 strong SRK methods as well as for order 2.0 weak SRK methods are given. As the main novelty, the computational complexity of the presented order 1.0 strong SRK method and the order 2.0 weak SRK method depends only linearly on the dimension of the driving Wiener process. This is a significant improvement compared to well known methods where the computational complexity depends quadratically on the dimension of the Wiener process.

It is possible to construct fully implicit Runge-Kutta methods like Gaus-Legendre, Radau-IA, Radau-IIA, Lobatto-IIIA, -IIIB, and -IIIC methods of arbitrary high order of convergence. The aim of this paper is to find a new adaptive time stepping for these classes which is based on the embedding technique. Adaptive time step control with embedding is well-known for Runge-Kutta methods, and therefore new embedded methods of order s-1 for the above classes of fully implicit Runge-Kutta methods are constructed. Since these fully implicit methods need the solution of a huge non-linear system of equations different approaches for non-linear equations are discussed and compared. It can be observed that non-linear solvers like the usually used simplified Newton method have a step size restriction if they are applied on higher order methods. We apply our new methods on some lower dimensional ODEs to show that our approach leads to an efficient method.

On strong convergence of explicit numerical methods for stochastic delay differential equations under non-global Lipschitz conditions

In this paper we present a family of balanced stochastic Runge-Kutta (BSRK) methods for stiff systems of stochastic differential equations (SDEs). This class of methods will have general framework for balancing the stochastic terms in arbitrary SRK methods of stochastic strong order p to get some balanced derivative-free methods of arbitrary stochastic order p. We also investigate the mean-square stability (MS-stability) properties for SDEs with stiffness in their drift and diffusion parts for some selected BSRK methods. Using such new methods, contrary to most existing derivative based balanced methods, in addition to having superior stability properties, one can get some significant reduction in computational cost. Numerical examples are presented to demonstrate the effectiveness of these methods for the pathwise approximation of solution of Ito SDEs.

This paper presents the cyber-physical system (CPS) of a numerical method (the widely-used Euler- Maruyama method) and establishes a foundational theory of the CPSs of numerical methods for stochastic differential equations (SDEs), which transforms the way we understand the relationship between the numerical method and the underlying dynamical system. Unlike in the literature where they are treated as separate systems linked by inequalities, the CPS is a seamless integration of the SDE and the numerical method and we construct a new and general class of stochastic impulsive differential equations (SiDEs) that can serve as a canonic form of the CPSs of numerical methods. By the CPS approach, we show the equivalence and intrinsic relationship between the stability of the SDE and the stability of the numerical method using the Lyapunov stability theory we develop for our class of SiDEs. Applying our established theory, we present the CPS Lyapunov inequality that is the necessary and sufficient condition for meansquare stability of the CPS of the Euler-Maruyama method for linear SDEs. The proposed CPS and theory initiate the study of systems numerics and provoke many open and interesting problems for future work.

Second-order stabilized explicit Runge–Kutta methods for stiff problems

The standard approach to applying IRK methods in the solution of two-point boundary value problems involves the solution of a non-linear system ofn×s equations in order to calculate the stages of the method, wheren is the number of differential equations ands is the number of stages of the implicit Runge-Kutta method. For two-point boundary value problems, we can select a subset of the implicit Runge-Kutta methods that do not require us to solve a non-linear system; the calculation of the stages can be done explicitly, as is the case for explicit Runge-Kutta methods. However, these methods have better stability properties than the explicit Runge-Kutta methods. We have called these new formulas two-point explicit Runge-Kutta (TPERK) methods. Their most important property is that, because their stages can be computed explicity, the solution of a two-point boundary value problem can be computed more efficiently than is possible using an implicit Runge-Kutta method. We have also developed a symmetric subclass of the TPERK methods, called ATPERK methods, which exhibit a number of useful properties.

The weak approximation of the solution of a system of Stratonovich stochastic differential equations with a m–dimensional Wiener process is studied. Therefore, a new class of stochastic Runge–Kutta methods is introduced. As the main novelty, the number of stages does not depend on the dimension m of the driving Wiener process which reduces the computational effort significantly. The colored rooted tree analysis due to the author is applied to determine order conditions for the new stochastic Runge–Kutta methods assuring convergence with order two in the weak sense. Further, some coefficients for second order stochastic Runge–Kutta schemes are calculated explicitly.

In the multidimensional case, second-order weak Runge-Kutta methods for stochastic differential equation (SDE) need simulation of correlated random variables, unless the diffusion matrix of SDE satisfies the commutativity condition. In this paper, we show that this can be avoided for some types of diffusion matrices and test functions important for applications.

A numerical method for some stochastic differential equations with multiplicative noise

In this paper we study a class of explicit pseudo two-step Runge-Kutta (EPTRK) methods for first-order ODEs for parallel computers. We investigate linear stability and derive methods with enlarged stability regions. In numerical experiments on a shared memory computer we compare a parallel variable step size EPTRK implementation with the effcient sequential Runge-Kutta method dopri5.

Reacting flow calculations are stiff for time-stepping. Recently, Zhong derived three sets of semi-implicit Runge-Kutta(SIRK) methods for split ordinary differential equations in the form of u' = /(«) 4g(u), where the nonstiff term / is treated explicitly while the stiff term g is simultaneously treated implicitly. The purpose of this work is twofold: first, to simulate multidimensional detonation waves using the high-order Rosenbrock semi-implicit Runge-Kutta time-stepping scheme(SIRK-3C); secondly to quantitatively compare the performance of SIRK with the second-order timesplitting schemes in dealing with the hyperbolic conservation laws with stiff relaxation terms. In the simulation of two-dimensional detonations, a full set of elementary hydrogen-oxygen kinetics with nine species are linear-implicitly treated while the explicit discretization of the basic advection equations are achieved by the third-order ENO schemes. Several model problems of hyperbolic systems with stiff relaxation are considered in the comparison of the Time-Splitting and the SIRK-3C scheme. A set of low-storage SIRK(LSSIRK) schemes is derived and tested to be high-order accurate and strongly A-stable. The results show that high-order SIRK time-stepping methods are suitable for stiff reactive flow simulations like the multi-dimensional detonations. INTRODUCTION In the reactive flow field solutions and the combustion related problems, the existence of several nonequilibrium states impose additional difficulties in the solutions of the reactive Euler equations with stiff relaxation terms. Since the smallest time scale is introduced from the chemical reaction kinetics, the second-order time-splitting method of Strang^ allows an independent calculation of the stiff ODE for the source term via implicit methods in one step, and an explicit highorder convective calculation in the following time steps. By the main source of stiffness from the basic hyperbolic conservation laws, it can accomplish a •Graduate Student, Member AIAA 'Assistant Professor, Mechanical and Aerospace Engineering Department, Member AIAA Copyright ©1997 by American Institute of Aeronautics and Astronautics, Inc. AH rights reserved. robust and easy-to-implement computation of reactive flow of current interests. Despite the second-order in time accuracy and the robustness of these methods, there are two main difficulties with the stiffness which have been minimally handled. The spurious solutions are one of the often observed phenomena, which are stable and free of oscillations, and yet may be completely incorrect.'-' Some type of resolution technique is needed to treat this known difficulties.^'!] The degrading of the second-order time accuracy to the first-order is another and quite recent observation when there exist solutions with many small-scaled structures of high-degree of complication. Jin® has modified the Strang splitting into the secondorder time-splitting methods which preserve its highorder accuracy when the considered solutions are complicated. New third-order semi-implicit Rosenbrock type Runge-Kutta scheme, based on the original version of Zhong'---', is computationally efficient with the linearized source term integration, thus requiring no iterations during the implicit treatment of the stiff source term. Third-order ENO schemes discretize the convective fluxes such that both implicit and explicit terms are treated at each of the three Runge-Kutta stages. This high-order Runge-Kutta scheme(SIRK-3C) still possesses the typical under-resolved characteristics of unphysical spurious numerical results of the TimeSplitting method. However the reduction to lower order if the small relaxation time is not temporally wellresolved, is successfully removed. Multi-dimensional Detonations Detonation waves are multi-dimensional and unstable phenomenon in nature as demonstrated by the early experiments of Urtiew and Oppenheinv-'. Existence of the triple points'^, consisting of an incident shock, a reflecting shock and a Mach stem, is the main characteristic of the reacting region behind the propagating shock front, and these detaching triple points from the leading front further contributes to the rolling up of vortices of opposite strength. One distinctive feature of the instability process observed in experiments is the formation of regular cell structures as triple points collide as they come together in the incident shock and move away from each other in the newly formed Mach stem. Figure[13] depicts this process of triple point col-

Convergence and stability of the backward Euler method for jump–diffusion SDEs with super-linearly growing diffusion and jump coefficients