Geometrically convergent simulation of the extrema of Lévy processes

We develop a novel approximate simulation algorithm for the joint law of the position, the running supremum, and the time of the supremum of a general Lévy process at an arbitrary finite time. We identify the law of the error in simple terms. We prove that the error decays geometrically in Lp (for any p≥1) as a function of the computational cost, in contrast with the polynomial decay for the approximations available in the literature. We establish a central limit theorem and construct nonasymptotic and asymptotic confidence intervals for the corresponding Monte Carlo estimator. We prove that the multilevel Monte Carlo estimator has optimal computational complexity (i.e., of order ϵ−2 if the mean squared error is at most ϵ2) for locally Lipschitz and barrier-type functions of the triplet and develop an unbiased version of the estimator. We illustrate the performance of the algorithm with numerical examples