[论文解读] Stability of fast algorithms for structured linear systems

Richard P. Brent|May 5, 2010

Matrix Theory and Algorithms参考文献 60被引用 13

一句话总结

本文研究了用于求解结构化线性系统（如托普利茨矩阵和汉克尔矩阵）的快速 $O(n^2)$ 算法的数值稳定性，重点关注在选主元过程中保持结构的方法。研究发现，尽管广义舒尔法和半正则方程法等算法具有弱稳定性，但快速正交分解方法通常因正交矩阵 $Q$ 的不稳定计算而失败，并提出了通过乔列斯基降阶和嵌入技术等方法来改进稳定性，但这些方法的计算成本更高。

ABSTRACT

We survey the numerical stability of some fast algorithms for solving systems of linear equations and linear least squares problems with a low displacement-rank structure. For example, the matrices involved may be Toeplitz or Hankel. We consider algorithms which incorporate pivoting without destroying the structure, and describe some recent results on the stability of these algorithms. We also compare these results with the corresponding stability results for the well known algorithms of Schur/Bareiss and Levinson, and for algorithms based on the semi-normal equations.

研究动机与目标

• 分析快速 $O(n^2)$ 算法在结构化线性系统中的数值稳定性，特别是针对低位移秩结构（如托普利茨矩阵和汉克尔矩阵）的系统。
• 比较莱文森法、舒尔/巴雷伊斯法和广义舒尔法与基于半正则方程和快速正交分解的算法在稳定性方面的差异。
• 阐明为何避免将 $Q$ 表示为正交矩阵乘积的快速算法尽管速度快，却常常不稳定。
• 探索在结构化分解中计算 $Q$ 或 $U$ 的稳定替代方法，包括基于嵌入的方法和改进的乔列斯基降阶技术。
• 澄清在结构化线性系统背景下弱稳定性与完全稳定性的区别，并强调证明快速算法稳定性所面临的挑战。

提出的方法

• 使用位移秩框架统一分析托普利茨、汉克尔和柯西矩阵等结构化矩阵。
• 应用布奇定义的稳定性与弱稳定性概念，评估快速算法的数值行为。
• 分析半正则方程方法：通过 $T^T T$ 的稳定乔列斯基分解求解 $\tilde{U}^T\tilde{U}x = T^Tb$，从而确保 $\tilde{U}$ 的弱稳定性。
• 考察将广义舒尔算法应用于 $T^T T$ 作为直接计算 $Q$ 的替代方法，指出该方法同样会使条件数的平方。
• 评估基于嵌入的方法，如 Chandrasekharan 和 Sayed 提出的方法，该方法将 $n \times n$ 矩阵嵌入到 $2n \times 2n$ 系统中，以实现稳定 $Q$ 的计算，但代价是常数因子更大。
• 将经典的 $O(n^3)$ 正交分解（如吉文斯、豪斯霍尔德）与 $O(n^2)$ 方法进行比较，表明后者因乔列斯基降阶对舍入误差的敏感性，无法稳定地形成 $Q$。

实验结果

研究问题

• RQ1 为何用于结构化线性系统的快速 $O(n^2)$ 算法尽管具有速度优势，却常常无法保持数值稳定性？
• RQ2 在求解结构化最小二乘问题时，广义舒尔算法的稳定性与半正则方程法相比如何？
• RQ3 能否构建用于计算结构化分解中正交因子 $Q$ 的稳定 $O(n^2)$ 算法？其代价是什么？
• RQ4 在半正则方程方法中，形成 $T^T T$ 时条件数的平方化会产生何种影响？这如何影响实际可解性？
• RQ5 为何乔列斯基降阶是快速正交分解算法在结构化矩阵中不稳定的根源？

主要发现

• 当结合 $T^T T$ 的稳定乔列斯基分解时，半正则方程法可实现 $O(n^2)$ 操作下对一般托普利茨系统的弱稳定性。
• 半正则方程法计算出的上三角矩阵 $\tilde{U}$ 满足 $\|T^T T - \tilde{U}^T \tilde{U}\| / \|T^T T\| = O_m(\varepsilon)$，表明其具有接近最优的向后误差。
• 通过乔列斯基降阶计算 $Q$ 的快速正交分解算法在数值上不稳定，因为降阶过程本身对舍入误差极为敏感。
• Chandrasekharan 和 Sayed 提出的基于嵌入的方法可实现 $T = LQU$ 分解的稳定计算，但其计算量高达 $59n^2 + O(n\log n)$，远高于 BBH 和半正则方程法的 $8n^2 + O(n\log n)$。
• 当条件数 $\kappa$ 落在 $1/\sqrt{\varepsilon} \leq \kappa \leq 1/\varepsilon$ 范围内时，迭代修正可能失效，需依赖双精度算术才能恢复精度。
• 本文结论认为，证明快速算法的稳定性比设计它们更具挑战性，且在未进行细致数值分析的情况下，若仅利用结构特性，稳定性往往无法保持。