# Session C6 - Stochastic Computation

Organizers: Mike Giles (University of Oxford, UK) - Arnulf Jentzen (Swiss Federal Institute of Technology Zurich, Switzerland) - Klaus Ritter (University of Kaiserslautern, Germany)

Workshop website

View session abstracts PDF

## Talks

December 18, 14:30 ~ 15:00 - Room A12

## Multi-Index Monte Carlo: When Sparsity Meets Sampling

### King Abdullah University of Science and Technology, Saudi Arabia   -   abdullateef.hajiali@kaust.edu.sa

We present a novel Multi-Index Monte Carlo (MIMC) method for weak approximation of stochastic models that are described in terms of differential equations either driven by random measures or with random coefficients. The MIMC method is both a stochastic version of the combination technique introduced by Zenger, Griebel and collaborators and an extension of the Multilevel Monte Carlo (MLMC) method first described by Heinrich and Giles. Inspired by Giles's seminal work, instead of using first-order differences as in MLMC, we use in MIMC high-order mixed differences to reduce the variance of the hierarchical differences dramatically. This in turn yields new and improved complexity results, which are natural generalizations of Giles's MLMC analysis, and which increase the domain of problem parameters for which we achieve the optimal complexity, $\mathcal{O}({\rm TOL}^{-2})$. We propose a systematic construction of optimal sets of indices for MIMC based on properly defined profits that in turn depend on the average cost per sample and the corresponding weak error and variance. Under standard assumptions on the convergence rates of the weak error, variance and work per sample, the optimal index set turns out to be of Total Degree (TD) type. In some cases, using optimal index sets, MIMC achieves a better rate for the computational complexity than does the corresponding rate when using Full Tensor sets. Moreover, we present a simple numerical example that illustrates the method and its computational benefits.

Joint work with Fabio Nobile (École Polytechnique Fédérale de Lausanne, Switzerland) and Raul Tempone (King Abdullah University of Science and Technology, Saudi Arabia).

View abstract PDF

December 18, 15:00 ~ 15:30 - Room A12

## Multilevel Monte Carlo for the simulation of dilute polymers

### University of Oxford, U.K.   -   mike.giles@maths.ox.ac.uk

Polymers immersed in a fluid can be modelled as chains connected by finitely-extensible bonds with nonlinear elastic potentials, subject to random forcing. In this work we discuss the simulation of the system of coupled SDEs which results from the modelling, and the use of multilevel Monte Carlo (MLMC) to efficiently estimate expectations arising from the associated equilibrium distribution. One important element is the use of adaptive timestepping which can be incorporated fairly easily into MLMC. Another is the use of a new multilevel coupling idea developed by Glynn and Rhee (2014) for the simulation of equilibrium expectations in the context of contracting Markov Chains.

Joint work with Endre Suli (Oxford), James Whittle (Oxford) and Shenghan Ye (Oxford).

View abstract PDF

December 18, 15:35 ~ 16:25 - Room A12

## Multi Level Monte Carlo for Coulomb Collisions in a Plasma

### University of California at Los Angeles, United States   -   caflisch@math.ucla.edu

This talk will describe the application of MLMC to solution of stochastic differential equations (SDEs) describing Coulomb collisions in a plasma. The collisions are between "test particles" and a Maxwellian distribution of "field Particles". The SDEs are three-dimensional and the Levy areas in the Milstein method are not tractable, but acceleration is achieved by reduction in the variance of the simulation with coarsest time step and by the antithetic method or a "Ito linearization method". Preliminary results are presented on MLMC for systems that depend on mean fields.

Joint work with Lee Ricketson (Courant Institute, NYU, USA), Mark Rosin (Pratt Institute, USA), Andris Dimits (Lawrence Livermore National Labs, USA) and Bruce Cohen (Lawrence Livermore National Labs, USA).

View abstract PDF

December 18, 17:00 ~ 17:30 - Room A12

## Weak approximation of stochastic differential equations by a multilevel Monte Carlo method using mean square adaptive numerical integration

### University of Oslo, Norway   -   haakonah1@gmail.com

In this talk, we present a multilevel Monte Carlo (MLMC) method for weak approximation of stochastic differential equations (SDE) that uses an a posteriori mean square error (MSE) adaptive Euler--Maruyama step-size control in the numerical integration of SDE realizations. MSE adaptivity is useful for weak approximation MLMC methods since it provides a reliable and efficient way to control the statistical error of the weak approximation MLMC estimator. For a large set of low-regularity weak approximation problems, the adaptive Euler--Maruyama method produces output whose weak error is bounded by $\mathcal{O}(\epsilon)$ at the cost $\mathcal{O}({\epsilon^{-2} |\log(\epsilon)|^4})$. This is a lower asymptotic cost than what can typically be obtained by the uniform time-step Euler-Maruyama MLMC method on the given set of problems. The cost reduction is illustrated in numerical studies.

Joint work with Juho Häppölä (Applied Mathematics and Computational Sciences, KAUST, Thuwal, Saudi Arabia) and Raúl Tempone (Applied Mathematics and Computational Sciences, KAUST, Thuwal, Saudi Arabia).

View abstract PDF

December 18, 17:30 ~ 18:00 - Room A12

## Higher order QMC Galerkin Discretization for parametric operator equations

### The University of New South Wales, Australia   -   josef.dick@unsw.edu.au

We discuss quasi-Monte Carlo methods to approximate the expected values of linear functionals of Petrov-Galerkin discretizations of parametric operator equations which depend on a possibly infinite sequence of parameters. Such problems arise in the numerical solution of differential and integral equations with random field inputs. We analyze the regularity of the solutions with respect to the parameters in terms of the rate of decay of the fluctuations of the input field. If $p\in (0,1]$ denotes the summability exponent'' corresponding to the fluctuations in affine-parametric families of operators, then we prove that deterministic interlaced polynomial lattice rules'' of order $\alpha = \lfloor 1/p \rfloor+1$ in $s$ dimensions with $N$ points can be constructed using a fast component-by-component algorithm, in $\mathcal{O}(\alpha\,s\, N\log N + \alpha^2\,s^2 N)$ operations, to achieve a convergence rate of $\mathcal{O}(N^{-1/p})$, with the implied constant independent of $s$. This dimension-independent convergence rate is superior to the rate $\mathcal{O}(N^{-1/p+1/2})$ for $2/3\leq p\leq 1$, which was recently established for randomly shifted lattice rules under comparable assumptions. In our analysis we use a non-standard Banach space setting and introduce smoothness-driven product and order dependent (SPOD)'' weights for which we develop a new fast CBC construction.

Joint work with F. Y. Kuo (UNSW), Q. T. Le Gia (UNSW), D. Nuyens (KU Leuven) and Ch. Schwab (ETH Z\"urich).

View abstract PDF

December 18, 18:00 ~ 18:30 - Room A12

## Optimal mesh hierarchies in Multilevel Monte Carlo methods

### Kungliga Tekniska Högskolan, KTH, Sweden   -   schwerin@csc.kth.se

We will discuss how to choose optimal mesh hierarchies in Multilevel Monte Carlo (MLMC) simulations based on uniform discretization methods with general approximation orders and computational costs. We will compare optimized geometric and non-geometric hierarchies and discuss how enforcing some domain constraints on parameters of MLMC hierarchies affects the optimality of these hierarchies. We also discuss the optimal tolerance splitting between the bias and the statistical error contributions and its asymptotic behavior.

Joint work with N. Collier, A. L. Haji-Ali, F. Nobile, and R. Tempone.

View abstract PDF

December 19, 14:30 ~ 15:00 - Room A12

## Quadrature for self-affine distributions on $\mathbf R^d$

### University of Passau, Germany   -   thomas.mueller-gronbach@uni-passau.de

We study numerical integration of $q$-times differentiable functions on $\mathbf R^d$ for a probability measure that is self-similar with respect to $m$ affine contraction mappings $S_1,\dots,S_m$ on $\mathbf R^d$ and corresponding probability weights $\rho_1,\dots,\rho_m$. Under mild conditions on the contractions we provide lower bounds for the worst case errors of deterministic as well as randomized algorithms in terms of the worst case (average) number of function evaluations that are used. The matching upper bounds are obtained by composite quadrature rules, which are easy to implement and are based on divide and conquer strategies that are adapted to the structure of the self-similarity. The optimal order of convergence is characterized in terms of the similarity dimension of the contractions.

Joint work with Steffen Dereich (University of Muenster).

View abstract PDF

December 19, 15:00 ~ 15:30 - Room A12

## A multilevel stochastic collocation method for PDEs with random inputs

### University of Warwick, United Kingdom   -   arethat24@gmail.com

By employing a hierarchy of both spatial approximations and interpolation operators in stochastic parameter space, we develop a multilevel version of stochastic collocation methods for random partial differential equations, leading to a significant reduction in computational cost. We provide a convergence and cost analysis of the new algorithm, and demonstrate the gains possible on a typical random diffusion model problem.

Joint work with Max Gunzburger (Florida State University, USA), Peter Jantsch (University of Tennessee, USA) and Clayton Webster (Oak Ridge National Laboratory, USA).

View abstract PDF

December 19, 15:35 ~ 16:25 - Room A12

## Multilevel Monte-Carlo Methods for hyperbolic PDEs with random input data

### SAM, ETH Zurich, Switzerland   -   schwab@math.ethz.ch

We consider random scalar hyperbolic conservation laws (RSCLs) in spatial dimension d with bounded random flux functions which are P-a.s. Lipschitz continuous with respect to the state variable.

There exists a unique random entropy solution (i.e., a strongly measurable mapping from a probability space into $C([0,T];L^1({\mathbb R}^d))$) with finite second moments.

We present a convergence analysis of a Multi-Level Monte-Carlo Front-Tracking (MLMCFT) algorithm.

It is based on pathwise'' application of the Front-Tracking Method for deterministic SCLs.

We compare the MLMCFT algorithms to the Multi-Level Monte-Carlo Finite-Volume methods. Due to the absence of a CFL time step restriction in the pathwise front tracking scheme, we can prove favourable complexity estimates: in spatial dimension $d\geq 2$, the mean field of the random entropy solution can be approximated numerically with (up to logarithmic terms) the same complexity as the solution of one instance of the deterministic problem, on the same mesh.

We then present results on large scale simulations of MLMC for wave propagation in heterogeneous media with log-Gaussian random coefficients. Here, conventional explicit timestepping schemes encounter the CFL constraint which, due to the lognormal Gaussian constitutive parameter, is random. A novel probabilistic complexity analysis and and adaptive load balancing algorithm achieve near linear strong scaling on up to 40K processors.

Joint work with Siddhartha Mishra (ETH), Nils Henrik Risebro (Oslo), Jonas Sukys (ETH) and Franziska Weber (Oslo).

View abstract PDF

December 19, 17:00 ~ 17:30 - Room A12

## Analytical approximations of BSDEs with non-smooth driver

### Ecole Polytechnique, France   -   emmanuel.gobet@polytechnique.edu

We provide and analyse analytical approximations of BSDEs in the limit of small non-linearity and short time to maturity, in the case of non-smooth drivers. We identify the first and the second order approximations within this asymptotics and consider two topical financial applications: the two interest rates problem and the Funding Value Adjustment. In high dimensional diffusion setting, we show how to compute explicitly the first order formula by taking advantage of recent proxy techniques. Numerical tests up to dimension 10 illustrate the efficiency of the numerical schemes. We also show that third order expansions may fail.

Joint work with Stefano Pagliarani (Ecole Polytechnique, France).

View abstract PDF

December 19, 17:30 ~ 18:00 - Room A12

## Simulation of forward-reverse stochastic representations for conditional diffusions

### Weierstrass Institute for Applied Analysis and Stochastics, Germany   -   christian.bayer@wias-berlin.de

We derive stochastic representations for the finite dimensional distributions of a multidimensional diffusion on a fixed time interval, conditioned on the terminal state. The conditioning can be with respect to a fixed point or more generally with respect to some subset. The representations rely on a reverse process connected with the given (forward) diffusion as introduced in Milstein et al. [Bernoulli 10(2):281--312, 2004] in the context of a forward-reverse transition density estimator. The corresponding Monte Carlo estimators have essentially root-N accuracy, hence they do not suffer from the curse of dimensionality. We provide a detailed convergence analysis and give a numerical example involving the realized variance in a stochastic volatility asset model conditioned on a fixed terminal value of the asset. We show that the algorithm can be applied for inference under incomplete data for the calculation of the expectation step in the EM algorithm. Finally, we give some applications.

Joint work with John Schoenmakers (Weierstrass Institute, Germany).

View abstract PDF

December 19, 18:00 ~ 18:30 - Room A12

## Customized fully implementable numerical schemes for FBSDEs

### University of Edinburgh, UK   -   l.szpruch@ed.ac.uk

In this talk we introduce a family of explicit numerical approximations for the forward backward stochastic differential equations (FBSDEs). We show that newly developed methodology allows to analyse BSDEs with drivers having polynomial growth and that are also monotone in the state variable. This offers a probabilistic scheme for wide class of reaction-diffusion PDEs. Proposed schemes preserve qualitative properties of the solutions to the FBSDEs for all ranges of time-steps. We conclude the talk by presenting a new efficient algorithm that allows to approximate conditional expectations in BSDEs setting. This leads to fully implementable numerical scheme.

Joint work with Goncalo dos Reis (Univeristy of Edinburgh), Arnaud Lionnet (University of Oxford) and Plamen Turkedjiev (Ecole Polytechnique).

View abstract PDF

December 20, 14:30 ~ 15:00 - Room C11

## On a SDE with no polynomial convergence rate for strong approximation at the final time

### University of Passau, Germany   -   larisa.yaroslavtseva@uni-passau.de

We consider the problem of strong approximation of the solution of a stochastic differential equation (SDE) at the final time based on point evaluations of the driving Brownian motion at a uniform grid. We present an example of a SDE with smooth and bounded coefficients for which no sequence of such approximations can achieve a polynomial rate of convergence. This generalizes a result from [1], which only covers the Euler scheme.

[1] Hairer,M., Hutzenthaler,M., Jentzen,A., Loss of regularity for Kolmogorov equations, Annals of Probability (to appear).

Joint work with Thomas Mueller-Gronbach (University of Passau, Germany).

View abstract PDF

December 20, 15:00 ~ 15:30 - Room C11

## Explicit numerical schemes for SDEs driven by Levy noise and for Stochastic Evolution Equations

### University of Edinburgh, UK   -   s.sabanis@ed.ac.uk

The idea of 'tamed' Euler schemes, which was pioneered by Hutzenthaler, Jentzen and Kloeden [1] and Sabanis [2], led to the development of a new generation of explicit numerical schemes

- for SDEs driven by Levy noise with superlinear coefficients and,

- for stochastic evolutions equations with super-linearly growing operators appearing in the drift.

Moreover, high order schemes (such as Milstein) are established (with optimal rates of convergence) by the natural extension of the aforementioned ideas. Theoretical results on this topic along with relevant simulation outputs will be presented during this talk.

[1] M. Hutzenthaler, A. Jentzen, P.E. Kloeden, Strong convergence of an explicit numerical method for SDEs with non-globally Lipschitz continuous coefficients. Ann. Appl. Probab. 22 (2012) 1611--1641.

[2] S. Sabanis, A note on tamed Euler approximations, Electron. Commun. Probab. 18 (2013), no. 47, 1–-10.

Joint work with Istvan Gyongy (University of Edinburgh), David Siska (University of Edinburgh), Chaman Kumar (University of Edinburgh) and Konstantinos Dareiotis (University of Edinburgh).

View abstract PDF

December 20, 15:30 ~ 16:00 - Room C11

## A perturbation formula as universal tool for strong approximations of stochastic differential equations

### University of Duisburg-Essen, Germany   -   martin.hutzenthaler@uni-due.de

The main object of this talk is a pathwise perturbation formula for stochastic differential equations (SDEs) which expresses the distance between the solution and any It\^o process in terms of the distances of the local characteristics. Together with suitable integrability properties, this is a convenient tool for proving strong convergence rates of various approximations of SDEs. For example, this yields a sufficient condition for local Lipschitz continuity in the strong sense in the initial value, a sufficient condition for explicit numerical approximations of finite-dimensional SDEs and a sufficient condition for spatial discretizations of nonlinear SPDEs. We illustrate these conditions with example SDEs from finance, physics and biology.

Joint work with Sonja Cox (University of Amsterdam, Netherlands), Arnulf Jentzen (ETH Zurich, Switzerland) and Xiaojie Wang (Changsha, China).

View abstract PDF

December 20, 16:00 ~ 16:30 - Room C11

## How to simulate stochastic differential equations without discretizing time

### Rutgers, USA   -   nawaf.bourabee@rutgers.edu

Numerical methods for SDEs are primarily based on pathwise approximation, basically if the Brownian motion is given in advance then the solution to the SDE (in a pathwise sense) looks like the solution to an ODE. Since the ODE solution is unique and smooth, it makes sense to discretize in time. In fact, solving the ODE is quite often thought of as an interpolation problem, where the largest time step is used to fit the approximation to the solution. The situation in the SDE context is a bit different: the solution to the SDE is continuous, but not differentiable in general, and an ensemble of trajectories originates from any initial point. Passing to the continuous limit is also nontrivial in the SDE context because existence of a pathwise unique solution imposes stringent conditions on the coefficients and domain of the SDE, which can be difficult to satisfy, and in fact, are not needed in practice. In actual computations one mainly considers the approximation of statistics of the law of the solution (mean first passage times, multi-time expectations, invariant measures, etc.), which is the main aim of SDE solvers, and the existence of a pathwise unique solution is not necessary for convergence in law. In this talk I will numerically illustrate why standard integrators are sometimes inconvenient to use, and present a different approach motivated by the issues encountered when simulating SDEs coming from quantitative finance, population dynamics, chemical kinetics, epidemiology, and polymeric fluids. This talk is based on joint work with Eric Vanden-Eijnden at NYU.

Joint work with Eric Vanden-Eijnden (NYU).

View abstract PDF

December 20, 17:00 ~ 17:30 - Room C11

## Weak approximation of the Heston model: non-smooth payoffs

### University of Mannheim, Germany   -   neuenkirch@kiwi.math.uni-mannheim.de

In this talk, we will study the weak approximation of the Heston price process for payoff functions, which are only measurable and bounded. The main tool for our analysis will be the explicit knowledge of the characteristic function of the Heston price, since we can not rely on the seminal work of Bally and Talay (1995). The latter requires smooth coefficients and Gaussian tails for the underlying SDE, which is not fulfilled for the Heston model.

Joint work with Martin Altmayer (University of Mannheim, Germany).

View abstract PDF

December 20, 17:30 ~ 18:00 - Room C11

## On a mild Ito formula for stochastic partial differential equations (SPDEs) and on weak convergence rates for SPDEs with nonlinear diffusion coefficients

### ETH Zurich, Switzerland   -   arnulf.jentzen@sam.math.ethz.ch

Strong convergence rates for (temporal, spatial, and noise) numerical approximations of semilinear stochastic evolution equations (SEEs) with smooth and regular nonlinearities are well understood in the scientific literature. Weak convergence rates for numerical approximations of such SEEs have been investigated since about 11 years and are far away from being well understood: roughly speaking, no essentially sharp weak convergence rates are known for parabolic SEEs with nonlinear diffusion coefficient functions; see Remark 2.3 in [A. Debussche, Weak approximation of stochastic partial differential equations: the nonlinear case, Math. Comp. 80 (2011), no. 273, 89--117] for details. In this talk we solve the weak convergence problem emerged from Debussche's article and establish esssentially sharp weak convergence rates for different type of temporal and spatial numerical approximations of semilinear SEEs with nonlinear diffusion coefficient functions. Our solution to the weak convergence problem does not use Malliavin calculus. Rather, a key ingredient in our solution to the weak convergence problem emerged from Debussche's article is a new -- somehow mild -- Ito type formula for solutions and numerical approximations of semilinear SEEs.

Joint work with Daniel Conus (Lehigh University, USA), Giuseppe Da Prato (Scuola Normale Superiore di Pisa, Italy), Ryan Kurniawan (University of Zurich and ETH Zurich, Switzerland), and Michael Röckner (Bielefeld University, Germany).

View abstract PDF