This Licentiate Thesis contains contributions in two challenging topics in matrix computations: the parallel multishift QR algorithm with aggressive early deflation (AED) and the matrix exponential of essentially nonnegative matrices. They are both structure-preserving algorithms in numerical linear algebra. We focus on performance in the former problem and on accuracy in the latter one.
The solution of matrix eigenvalue problems is a fundamental topic in numerical linear algebra. The QR algorithm which computes the Schur decomposition of a matrix is by far the most important approach for solving dense nonsymmetric eigenvalue problems. Recently a novel parallel QR algorithm has been developed by incorporating some modern techniques such as small-bulge multishift and AED. The novel parallel approach significantly outperforms the pipelined QR algorithm in ScaLAPACK v1.8.0 and earlier versions. But AED becomes a computational bottleneck in the new parallel QR algorithm. We develop multilevel AED algorithms which indeed decrease the total amount of communications and further improve the performance of the parallel QR algorithm. We also identify and remedy some anomalies in the original ScaLAPACK implementation. The improved version of the new parallel QR algorithm is now available as a part of ScaLAPACK version 2.0. Both performance models and numerical experiments demonstrate the efficiency of the new approach.
The computation of the matrix exponential is also a classical topic in numerical linear algebra. The exponential of essentially nonnegative matrices (real square matrices with nonnegative off-diagonal entries) is an important special case since it has applications in continuous-time Markov processes and positive linear dynamical systems. Because of the nonnegativity, this problem is usually well-conditioned in the sense of componentwise perturbations. We derive lower and upper bounds of the matrix exponential, and combining the scaling and squaring method with aggressively truncated Taylor expansion. More precisely, we establish new à priori error estimates and use them to propose an efficient strategy to balance the scale factor and the order of expansion. This leads to new algorithms with componentwise high relative accuracy. Rounding error analyses and numerical experiments confirm the efficiency and accuracy of our algorithms.