×

Computing fundamental matrix decompositions accurately via the matrix sign function in two iterations: the power of Zolotarev’s functions. (English) Zbl 1383.15012

Summary: The symmetric eigenvalue decomposition and the singular value decomposition (SVD) are fundamental matrix decompositions with many applications. Conventional algorithms for computing these decompositions are suboptimal in view of recent trends in computer architectures which require minimizing communication together with arithmetic costs. Spectral divide-and-conquer algorithms, which recursively decouple the problem into two smaller subproblems, can achieve both requirements. Such algorithms can be constructed with the polar decomposition playing two key roles: it forms a bridge between the symmetric eigendecomposition and the SVD, and its connection to the matrix sign function naturally leads to spectral-decoupling. For computing the polar decomposition, the scaled Newton and QDWH iterations are two of the most popular algorithms, as they are backward stable and converge in at most nine and six iterations, respectively. Following this framework, we develop a higher-order variant of the QDWH iteration for the polar decomposition. The key idea of this algorithm comes from approximation theory: we use the best rational approximant for the scalar sign function due to G. Zolotarev [Abh. St. Petersb. 30 (Russisch) (1877; JFM 09.0343.02), Reprint in Collected Works, Vol. II, Izdat. Akad. Nauk SSSR, Moscow, 1932, pp. 1–59 (in Russian)]. The algorithm exploits the extraordinary property enjoyed by the sign function that a high-degree Zolotarev function (best rational approximant) can be obtained by appropriately composing low-degree Zolotarev functions. This lets the algorithm converge in just two iterations in double-precision arithmetic, with the whopping rate of convergence seventeen. The resulting algorithms for the symmetric eigendecompositions and the SVD have higher arithmetic costs than the QDWH-based algorithms, but are better suited for parallel computing and exhibit excellent numerical backward stability.

MSC:

15A23 Factorization of matrices
65F15 Numerical computation of eigenvalues and eigenvectors of matrices
65F30 Other matrix algorithms (MSC2010)
65G50 Roundoff error

Citations:

JFM 09.0343.02

Software:

mftoolbox; JDQZ; FEAST; JDQR

References:

[1] M. Abramowitz and I. A. Stegun, {\it Handbook of Mathematical Functions}, Dover, New York, 1965. · Zbl 0171.38503
[2] N. I. Akhiezer, {\it Elements of the Theory of Elliptic Functions}, AMS, Providence, RI, 1990. · Zbl 0694.33001
[3] N. I. Akhiezer, {\it Theory of Approximation}, Dover, New York, 1992. · Zbl 0802.34054
[4] Z. Bai and J. Demmel, {\it Using the matrix sign function to compute invariant subspaces}, SIAM J. Matrix Anal. Appl., 19 (1998), pp. 205-225. · Zbl 0914.65035
[5] Z. Bai, J. Demmel, J. Dongarra, A. Ruhe, and H. van der Vorst, {\it Templates for the Solution of Algebraic Eigenvalue Problems: A Practical Guide}, SIAM, Philadelphia, 2000. · Zbl 0965.65058
[6] Z. Bai, J. Demmel, and M. Gu, {\it An inverse free parallel spectral divide and conquer algorithm for nonsymmetric eigenproblems}, Numer. Math., 76 (1997), pp. 279-308. · Zbl 0876.65021
[7] G. Ballard, J. Demmel, and I. Dumitriu, {\it Minimizing Communication for Eigenproblems and the Singular Value Decomposition}, Technical Report 237, LAPACK Working Note, 2010.
[8] G. Ballard, J. Demmel, O. Holtz, and O. Schwartz, {\it Minimizing communication in numerical linear algebra}, SIAM J. Matrix Anal. Appl., 32 (2011), pp. 866-901. · Zbl 1246.68128
[9] G. Ballard, J. Demmel, and N. Knight, {\it Avoiding Communication in Successive Band Reduction}, Technical Report UCB/EECS-2013-131, University of California, Berkeley, 2013.
[10] A. N. Beavers, Jr. and E. D. Denman, {\it A computational method for eigenvalues and eigenvectors of a matrix with real eigenvalues}, Numer. Math., 21 (1973), pp. 389-396. · Zbl 0296.65017
[11] B. Beckermann, {\it Optimally Scaled Newton Iterations for the Matrix Square Root}, talk at the Advances in Matrix Functions and Matrix Equations workshop, Manchester, UK, 2013.
[12] B. Beckermann and L. Reichel, {\it Error estimates and evaluation of matrix functions via the Faber transform}, SIAM J. Numer. Anal., 45 (2009), pp. 3849-3883. · Zbl 1204.65041
[13] D. Braess, {\it On rational approximation of the exponential and the square root function}, in Rational Approximation and Interpolation, P. R. Graves-Morris, E. B. Saff, and R. S. Varga, eds., Lecture Notes in Math. 1105, Springer, Berlin, Heidelberg, New York, 1984, pp. 89-99. · Zbl 0556.41013
[14] D. Braess, {\it Nonlinear Approximation Theory}, Springer, Berlin, Heidelberg, New York, 1986. · Zbl 0656.41001
[15] R. Byers and H. Xu, {\it A new scaling for Newton’s iteration for the polar decomposition and its backward stability}, SIAM J. Matrix Anal. Appl., 30 (2008), pp. 822-843. · Zbl 1170.65019
[16] B. C. Carlson and J. L. Gustafson, {\it Asymptotic expansion of the first elliptic integral}, SIAM J. Math. Anal., 16 (1985), pp. 1072-1092. · Zbl 0593.33002
[17] E. W. Cheney, {\it Introduction to Approximation Theory}, Chelsea, New York, 1982. · Zbl 0535.41001
[18] T. F. Cox and M. A. A. Cox, {\it Multidimensional Scaling}, 2nd ed., CRC Press, Boca Raton, FL, 2000. · Zbl 1147.68460
[19] I. S. Dhillon and B. N. Parlett, {\it Orthogonal eigenvectors and relative gaps}, SIAM J. Matrix Anal. Appl., 25 (2004), pp. 858-899. · Zbl 1068.65046
[20] V. Druskin, S. Güttel, and L. Knizhnerman, {\it Near-optimal perfectly matched layers for indefinite Helmholtz problems}, SIAM Rev., 58 (2016), pp. 90-116. · Zbl 1344.35009
[21] S. W. Gaaf and M. E. Hochstenbach, {\it Probabilistic bounds for the matrix condition number with extended Lanczos bidiagonalization}, SIAM J. Sci. Comput., 37 (2015), pp. S581-S601. · Zbl 1325.65051
[22] A. George and K. Ikramov, {\it Is the polar decomposition finitely computable?}, SIAM J. Matrix Anal. Appl., 17 (1996), pp. 348-354. · Zbl 0853.65046
[23] A. George and K. Ikramov, {\it Addendum: Is the polar decomposition finitely computable?}, SIAM J. Matrix Anal. Appl., 18 (1997), p. 264. · Zbl 0870.65037
[24] G. H. Golub and C. F. Van Loan, {\it Matrix Computations}, 4th ed., The Johns Hopkins University Press, Baltimore, MD, 2013. · Zbl 1268.65037
[25] A. A. Gončar, {\it Zolotarev problems connected with rational functions}, Math. USSR Sb., 7 (1969), pp. 623-635. · Zbl 0203.07302
[26] J. C. Gower and G. B. Dijksterhuis, {\it Procrustes Problems}, Oxford University Press, Oxford, UK, 2004. · Zbl 1057.62044
[27] M. Gu and S. C. Eisenstat, {\it A divide-and-conquer algorithm for the bidiagonal SVD}, SIAM J. Matrix Anal. Appl., 16 (1995), pp. 79-92. · Zbl 0821.65019
[28] M. Gu and S. C. Eisenstat, {\it A divide-and-conquer algorithm for the symmetrical tridiagonal eigenproblem}, SIAM J. Matrix Anal. Appl., 16 (1995), pp. 172-191. · Zbl 0815.65050
[29] S. Güttel, E. Polizzi, P. Tang, and G. Viaud, {\it Zolotarev Quadrature Rules and Load Balancing for the FEAST Eigensolver}, MIMS EPrint 2014.39, The University of Manchester, UK, 2014. · Zbl 1321.65055
[30] N. Hale, N. J. Higham, and L. N. Trefethen, {\it Computing \(A^α, \log(A)\), and related matrix functions by contour integrals}, SIAM J. Numer. Anal., 46 (2008), pp. 2505-2523. · Zbl 1176.65053
[31] N. Halko, P. G. Martinsson, and J. A. Tropp, {\it Finding structure with randomness: Probabilistic algorithms for constructing approximate matrix decompositions}, SIAM Rev., 53 (2011), pp. 217-288. · Zbl 1269.65043
[32] P. Henrici, {\it Applied and Computational Complex Analysis. Vol. 1. Power Series, Integration, Conformal Mapping, Location of Zeros}, Wiley, New York, 1974. · Zbl 0313.30001
[33] N. J. Higham, {\it Functions of Matrices: Theory and Computation}, SIAM, Philadelphia, 2008. · Zbl 1167.15001
[34] N. J. Higham and P. Papadimitriou, {\it A parallel algorithm for computing the polar decomposition}, Parallel Comput., 20 (1994), pp. 1161-1173. · Zbl 0815.65061
[35] N. J. Higham, {\it Computing the polar decomposition–with applications}, SIAM J. Sci. Statist. Comput., 7 (1986), pp. 1160-1174. · Zbl 0607.65014
[36] N. J. Higham, {\it The matrix sign decomposition and its relation to the polar decomposition}, Linear Algebra Appl., 212/213 (1994), pp. 3-20. · Zbl 0816.15012
[37] M. E. Hochstenbach, {\it Eigenvalue Tools}; available via http://www.win.tue.nl/ hochsten/eigenvaluetools/. · Zbl 1139.65028
[38] S. Huss-Lederman, E. S. Quintana-Ortí, X. Sun, and Y.-J. Y. Wu, {\it Parallel spectral division using the matrix sign function for the generalized eigenproblem}, Internat. J. High Speed Com., 11 (2000), pp. 1-14. · Zbl 0973.65028
[39] B. Iannazzo, {\it The geometric mean of two matrices from a computational viewpoint}, Numer. Linear Algebra Appl., 23 (2016), pp. 208-229. · Zbl 1374.65083
[40] T. Kaneko, S. Fiori, and T. Tanaka, {\it Empirical arithmetic averaging over the compact Stiefel manifold}, IEEE Trans. Signal Process., 61 (2013), pp. 883-894. · Zbl 1394.62063
[41] A. D. Kennedy, {\it Fast evaluation of Zolotarev coefficients}, in Proceedings of the Third International Workshop on Numerical Analysis and Lattice QCD, Edinburgh, Scotland, 2003, pp. 169-189. · Zbl 1082.33012
[42] A. D. Kennedy, {\it Approximation theory for matrices}, Nuclear Phys. B Proc. Suppl., 128 (2004), pp. 107-116.
[43] C. Kenney and A. Laub, {\it A hyperbolic tangent identity and the geometry of Padé sign function iterations}, Numer. Algorithms, 7 (1994), pp. 111-128, doi:10.1007/BF02140677. · Zbl 0812.65035
[44] C. Kenney and A. Laub, {\it The matrix sign function}, IEEE Trans. Automat. Control, 40 (1995), pp. 1330-1348. · Zbl 0830.93022
[45] B. Laszkiewicz and K. Zietak, {\it Numerical experiments with algorithms for the ADI and Zolotarev coefficients}, Appl. Math. Comput., 206 (2008), pp. 298-312. · Zbl 1154.65307
[46] J. Mao, {\it Optimal orthonormalization of the strapdown matrix by using singular value decomposition}, Comput. Math. Appl., 12 (1986), pp. 353-362. · Zbl 0606.65027
[47] Y. Nakatsukasa, {\it Algorithms and Perturbation Theory for Matrix Eigenvalue Problems and the Singular Value Decomposition}, Ph.D. thesis, University of California, Davis, 2011.
[48] Y. Nakatsukasa, Z. Bai, and F. Gygi, {\it Optimizing Halley’s iteration for computing the matrix polar decomposition}, SIAM J. Matrix Anal. Appl., 31 (2010), pp. 2700-2720. · Zbl 1215.65081
[49] Y. Nakatsukasa and N. J. Higham, {\it Backward stability of iterations for computing the polar decomposition}, SIAM J. Matrix Anal. Appl., 33 (2012), pp. 460-479. · Zbl 1252.65083
[50] Y. Nakatsukasa and N. J. Higham, {\it Stable and efficient spectral divide and conquer algorithms for the symmetric eigenvalue decomposition and the SVD}, SIAM J. Sci. Comput., 35 (2013), pp. A1325-A1349. · Zbl 1326.65049
[51] I. Ninomiya, {\it Best rational starting approximations and improved Newton iteration for the square root}, Math. Comp., 24 (1970), pp. 391-404. · Zbl 0213.16402
[52] B. N. Parlett, {\it The Symmetric Eigenvalue Problem}, SIAM, Philadelphia, 1998. · Zbl 0885.65039
[53] P. P. Petrushev and V. A. Popov, {\it Rational Approximation of Real Functions}, Cambridge University Press, Cambridge, UK, 2011. · Zbl 1209.41002
[54] M. J. D. Powell, {\it Approximation Theory and Methods}, Cambridge University Press, Cambridge, UK, 1981. · Zbl 0453.41001
[55] J. D. Roberts, {\it Linear model reduction and solution of the algebraic Riccati equation by use of the sign function}, Internat. J. Control, 32 (1980), pp. 677-687. · Zbl 0463.93050
[56] H. Rutishauser, {\it Betrachtungen zur Quadratwurzeliteration}, Monatsh. Math., 67 (1963), pp. 452-464. · Zbl 0143.39203
[57] D. E. Sukkari, H. Ltaief, and D. E. Keyes, {\it A High Performance QDWH-SVD Solver Using Hardware Accelerators}, Technical report, KAUST Repository, Thuwal, Saudi Arabia, 2015. · Zbl 1369.65058
[58] L. N. Trefethen, J. A. C. Weideman, and T. Schmelzer, {\it Talbot quadratures and rational approximations}, BIT, 46 (2006), pp. 653-670. · Zbl 1103.65030
[59] L. N. Trefethen, {\it Approximation Theory and Approximation Practice}, SIAM, Philadelphia, 2013. · Zbl 1264.41001
[60] L. N. Trefethen and J. A. C. Weideman, {\it The exponentially convergent trapezoidal rule}, SIAM Rev., 56 (2014), pp. 385-458. · Zbl 1307.65031
[61] J. van den Eshof, {\it Nested Iteration Methods for Nonlinear Matrix Problems}, Ph.D. thesis, Proefschrift Universiteit Utrecht, 2003.
[62] J. van den Eshof, T. Lippert, A. Frommer, K. Schilling, and H. van der Vorst, {\it Numerical methods for the QCD overlap operator: I. Sign-function and error bounds}, Comput. Phys. Comm., 146 (2002), pp. 203-224. · Zbl 1008.81508
[63] F. V. Zee, R. van de Geijn, and G. Quintana-Orti, {\it Restructuring the tridiagonal and bidiagonal QR algorithm for performance}, ACM Trans. Math. Software, 40 (2014), article 18. · Zbl 1322.65051
[64] Z. Zhang, H. Zha, and W. Ying, {\it Fast parallelizable methods for computing invariant subspaces of Hermitian matrices}, J. Comput. Math., 25 (2007), pp. 583-594. · Zbl 1140.65315
[65] E. I. Zolotarev, {\it Application of elliptic functions to questions of functions deviating least and most from zero}, Zap. Imp. Akad. Nauk. St. Petersburg, 30 (5), 1877. Reprinted in his Collected Works, Vol. II, Izdat. Akad. Nauk SSSR, Moscow, 1932, pp. 1-59 (in Russian).
This reference list is based on information provided by the publisher or from digital mathematics libraries. Its items are heuristically matched to zbMATH identifiers and may contain data conversion errors. In some cases that data have been complemented/enhanced by data from zbMATH Open. This attempts to reflect the references listed in the original paper as accurately as possible without claiming completeness or a perfect matching.