[Paper Review] Residual, restarting and Richardson iteration for the matrix exponential, revised

Mike A. Botchev|arXiv (Cornell University)|Dec 23, 2011

Matrix Theory and Algorithms38 references3 citations

TL;DR

This paper introduces a novel residual concept for the matrix exponential based on interpreting the solution as the outcome of an initial-value problem (IVP) for a system of ODEs. By defining the residual as the violation of the ODE system by an approximate solution, the method enables reliable, computable stopping criteria for Krylov subspace and Chebyshev methods, and leads to new residual-based iterative solvers, including a restarted Krylov-Richardson method that improves convergence efficiency.

ABSTRACT

A well-known problem in computing some matrix functions iteratively is the lack of a clear, commonly accepted residual notion. An important matrix function for which this is the case is the matrix exponential. Suppose the matrix exponential of a given matrix times a given vector has to be computed. We develop the approach of Druskin, Greenbaum and Knizhnerman (1998) and interpret the sought-after vector as the value of a vector function satisfying the linear system of ordinary differential equations (ODE) whose coefficients form the given matrix. The residual is then defined with respect to the initial-value problem for this ODE system. The residual introduced in this way can be seen as a backward error. We show how the residual can be computed efficiently within several iterative methods for the matrix exponential. This completely resolves the question of reliable stopping criteria for these methods. Further, we show that the residual concept can be used to construct new residual-based iterative methods. In particular, a variant of the Richardson method for the new residual appears to provide an efficient way to restart Krylov subspace methods for evaluating the matrix exponential.

Motivation & Objective

• To address the longstanding lack of a reliable residual notion for iterative matrix function computation, particularly for the matrix exponential.
• To provide a computable and meaningful backward error measure for approximate solutions of the matrix exponential via an ODE interpretation.
• To enable robust and efficient stopping criteria in Krylov subspace and polynomial iterative methods for the matrix exponential.
• To develop new residual-based iterative algorithms, including a restarted variant of the Richardson method for Krylov methods.
• To extend the residual framework to other matrix functions such as the matrix cosine and sine.

Proposed method

• Define the matrix exponential as the solution of an initial-value problem (IVP): y'(t) = -Ay(t), y(0) = v, with exact solution y(t) = exp(-tA)v.
• Introduce the residual r_k(t) = -A y_k(t) - y_k'(t) for an approximate solution y_k(t), which measures the violation of the ODE system.
• Use the residual norm as a stopping criterion for iterative methods, ensuring backward stability and reliable convergence monitoring.
• Develop a residual-based Richardson iteration for the matrix exponential, which can be used to efficiently restart Krylov subspace methods.
• Apply the residual concept to Krylov subspace methods (e.g., Arnoldi/Lanczos) and Chebyshev polynomial methods, showing efficient residual computation.
• Extend the framework to matrix cosine and sine functions by identifying corresponding second-order ODEs and defining residuals accordingly.

Experimental results

Research questions

• RQ1 How can a reliable and computable residual be defined for iterative methods solving the matrix exponential?
• RQ2 Can the residual concept be used to construct new iterative solvers for the matrix exponential?
• RQ3 How does the residual-based stopping criterion improve the robustness and efficiency of Krylov subspace methods?
• RQ4 Can the residual-based Richardson method serve as an effective preconditioner or restart strategy for Krylov methods?
• RQ5 To what extent can the residual framework be generalized to other matrix functions like cosine and sine?

Key findings

• The proposed residual, defined as the residual of the ODE system, serves as a backward error and enables reliable stopping criteria for iterative matrix exponential solvers.
• The residual can be efficiently computed within Krylov subspace and Chebyshev polynomial methods, allowing real-time convergence monitoring.
• A new residual-based Richardson iteration is proposed, which, when combined with Krylov methods, enables efficient and effective restarting, reducing the number of outer iterations.
• Numerical experiments show that the Krylov-Richardson method with the SaI strategy converges faster than standard methods, particularly when the initial vector v satisfies boundary conditions.
• The convergence rate improves significantly when the initial vector v is compatible with the boundary conditions of the underlying PDE, as demonstrated in a 3D diffusion problem.
• The residual framework is extendable to other matrix functions such as cosine and sine by identifying their governing second-order ODEs and defining residuals accordingly.