GPU ACCELERATED ALGORITHMS FOR COMPUTING MATRIX FUNCTION VECTOR PRODUCTS WITH APPLICATIONS TO EXPONENTIAL INTEGRATORS AND FRACTIONAL DIFFUSION

The efficient computation of matrix function vector products has recently become an important area of research, driven in particular by two important applications: the numerical solution of fractional partial differential equations and the integration of large systems of ordinary differential equations. In this work we consider a problem that combines these two applications in the form of a numerical solution algorithm for fractional reaction-diffusion equations that, after spatial discretization, is advanced in time using the exponential Euler method. We focus on the efficient implementation of the algorithm on graphics processing units (GPUs), as we wish to make use of the increased computational power available with this hardware. We compute the matrix function vector products using the contour integration method in [N. Hale, N. J. Higham, and L. N. Trefethen, SIAM J. Numer. Anal., 46 (2008), pp. 2505 2523]. Multiple levels of preconditioning are applied to reduce the GPU memory footprint and to further accelerate convergence. We also derive an error bound for the convergence of the contour integral method that allows us to predetermine the appropriate number of quadrature points. Results are presented that demonstrate the effectiveness of the method for large two-dimensional problems, showing a speedup of more than an order of magnitude compared to a CPU-only implementation.