A Unified Framework for Numerically Inverting Laplace Transforms

We introduce and investigate a framework for constructing algorithms to invert Laplace transforms numerically. Given a Laplace transform \hat{f} of a complex-valued function of a nonnegative real-variable, f, the function f is approximated by a finite linear combination of the transform values; i.e., we use the inversion formula f(t) \approx f_n (t) \equiv \frac{1}{t} \sum_{k = 0}^{n}\omega_{k}\hat{f}\biggl(\frac{\alpha_{k}}{t}\biggr),\quad 0 < t < \infty, where the weights ω k and nodes α k are complex numbers, which depend on n, but do not depend on the transform \hat{f} or the time argument t. Many different algorithms can be put into this framework, because it remains to specify the weights and nodes. We examine three one-dimensional inversion routines in this framework: the Gaver-Stehfest algorithm, a version of the Fourier-series method with Euler summation, and a version of the Talbot algorithm, which is based on deforming the contour in the Bromwich inversion integral. We show that these three building blocks can be combined to produce different algorithms for numerically inverting two-dimensional Laplace transforms, again all depending on the single parameter n. We show that it can be advantageous to use different one-dimensional algorithms in the inner and outer loops.