Exponential time differencing for stiff systems with nondiagonal linear part

Abstract

Exponential time differencing methods provide instability-free explicit schemes for systems with fast decaying or oscillating modes (stiff systems), without limitation on the time step size. Moreover, with these methods, one can drastically diminish the error accumulation rate for numerical simulation of conservative systems. The methods yield an especially large performance gain for PDEs with high order of spatial derivatives. Simultaneously, the problem of analytical calculation of coefficients of exponential time differencing schemes becomes laborious or unsolvable in the case of a nondiagonal form of the principal linear part of equations. We introduce an approach, where the scheme coefficients are obtained from the direct numerical integration of certain auxiliary problems over a short time interval - one scheme step size. The approach is universal and its implementation is illustrated with four examples: analytically solvable system of two first-order ODEs, one-dimensional reaction-diffusion system under time-dependent conditions, two-dimensional reaction-diffusion system under time-independent and time-dependent conditions, and one-dimensional Cahn-Hilliard equation with constant coefficients. The employment of an exponential time differencing method of the two-step Runge-Kutta type yields a simulation performance gain for the diffusion-type equation, with program optimization made. Without program optimization, the performance gain increases by one order with respect to the spatial step size for each order of the highest spatial derivative, and appears starting from the third order of the derivative. With the tested method, one can extend the study of an analogue of the Anderson localization to two- and three-dimensional active media and achieve an acceptable performance for the direct numerical simulations of dynamics of the probability density function for active Brownian particles.