Optimizing Talbot's contours for the inversion of the Laplace transform

Talbot's method for the numerical inversion of the Laplace transform consists of numerically integrating the Bromwich integral on a special contour by means of the trapezoidal or midpoint rules. In this paper we address the issue of parameter selection in the method, for the particular situation when parabolic PDEs are solved. In the process the well-known subgeometric convergence rate O(exp(-c root N)) of this method is improved to the geometric rate O(exp(-cN)), with N the number of nodes in the integration rule. The value of the maximum decay rate c is explicitly determined. Numerical results for two versions of the heat equation are presented. With the choice of parameters derived here, the rule of thumb is that to achieve an accuracy of 10(-l) at any given time, the associated elliptic problem has to be solved no more than l times.