數值線性代數

數值線性代數（numerical linear algebra），又稱應用線性代數（applied linear algebra）是一門研究在電腦上進行線性代數計算，特別是矩陣運算演算法的學科，是數值分析的一個分支。電腦用浮點數運算，無法精確表示無理數數據，因此電腦演算法應用於數據矩陣時，有時會產生較大的誤差。數值線性代數利用向量和矩陣的性質開發演算法，以最大限度地減少這誤差，同時還要確保演算法儘可能高效。

數值線性代數試圖利用精度有限的電腦解決連續的數學問題，因此在自然科學和社會科學中的應用同連續數學一樣廣泛。這些問題包括圖像處理、信號處理、金融工程學、材料科學模擬、結構生物學、數據挖掘、生物資訊學、流體動力學和其他很多領域。矩陣方法尤其用於有限差分法、有限元法和微分方程建模。尼克·特雷費森和大衛·鮑三世（David Bau III）注意到數值線性代數的廣泛應用，認為它「同微積分和微分方程一樣是數學科學的基礎」，儘管是個相對較小的領域。由於矩陣和向量的很多性質也適用於函數和算子，因此數值線性代數也可視作泛函分析中注重實用演算法的一支。

數值線性代數中的常見問題如LU分解、QR分解、奇異值分解、特徵分解等，然後可用於解答常見的線性代數問題，如求解線性方程組、定位特徵值、最小二乘優化等。數值線性代數的核心問題是開發在有限精度電腦上應用真實數據時不會引入誤差的演算法，這通常通過迭代法來實現，而非直接方法。

## 歷史

數值線性代數是由約翰·馮·諾伊曼、艾倫·圖靈、詹姆斯·哈迪·威爾金森、阿爾斯通·斯科特·豪斯霍爾德、喬治·福賽思、海因茨·魯蒂紹爾等電腦先驅開發的，目的是將最早的電腦應用於連續數學問題，如彈道問題與偏微分方程求解。約翰·馮·諾依曼和赫爾曼·戈德斯坦在1947年的工作中，首次認真嘗試將演算法應用於真實數據計算時的誤差降至最低。隨著技術發展，研究者越來越有能力在超大型高精度矩陣解決複雜問題，這一領域也在不斷壯大，而一些數值演算法也因並行計算等技術，使其成為解決科學問題的使用方法而日益突出。

## 矩陣分解

### 分割矩陣

對應用線性代數中的很多問題，將矩陣視作列向量的組合是很有用的。例如，求解線性系統$x = A^{-1}b$時，與其將x理解為$A^{-1}$與b的積，不如將x視作b在A的列形成的基上線性展開的係數向量。將矩陣視作列之組合，也是矩陣演算法的一種實用方法。這是因為矩陣演算法常常包含兩個嵌套的循環，一個是A的列循環，一個是A的行循環。例如，對於矩陣$A^{m \times n}$和向量$x^{n \times 1}$、$y^{m \times 1}$，可以用列劃分的方法計算$y:=Ax+y$

```matlab
for q = 1:n
  for p = 1:m
    y(p) = A(p,q)*x(q) + y(p)
  end
end
```

### 奇異值分解

矩陣$A^{m \times n}$的奇異值分解是$A = U \Sigma V^\ast$，其中U、V是酉矩陣，$\Sigma$是對角矩陣。$\Sigma$的對角項稱作A的奇異值。由於奇異值是$AA^\ast$特徵值的平方根，所以奇異值分解同特徵值分解有緊密聯繫。這意味著，計算奇異值分解的大多數方法都與特徵值方法類似，最常用的方法可能涉及豪斯霍德變換。

### QR分解

矩陣$A^{m \times n}$的QR分解是$Q^{m \times m}$與$R^{m \times n}$，使得$A=QR$，其中Q是正交矩陣，R是上三角矩陣。計算QR分解的兩種主要演算法是格拉姆-施密特正交化與豪斯霍德變換。QR分解常用於解決線性最小二乘問題與特徵值問題（通過迭代QR演算法）。

### LU分解

矩陣A的LU分解由下三角矩陣L和上三角矩陣U構成，使得$A=LU$。矩陣U通過上三角化過程生成，涉及將A左乘一系列矩陣$M_1,\ldots,M_{n-1}$，以形成積$M_{n-1} \cdots M_1 A = U$，因此等價地$L = M_1^{-1} \cdots M_{n-1}^{-1}$。

### 特徵值分解

矩陣$A^{m \times m}$的特徵值分解是$A = X \Lambda X^{-1}$，其中X的列是A的特徵向量，$\Lambda$是對角矩陣，對角項是A的相應特徵值。目前還沒有直接得到任意矩陣特徵值分解的直接方法，因為程序不可能在有限時間內精確求得任意多項式的根，所以任何通用的特徵值求解器必是迭代的。

## 演算法

### 高斯消元法

從數值線性代數的視角看，高斯消元法是將A進行LU分解的一種方法。高斯消元法通過將A與一系列矩陣$L_{m-1} \cdots L_2 L_1 A = U$左乘來實現，直到U是上三角矩陣，L是下三角矩陣，且$L \equiv L_1^{-1}L_2^{-1} \cdots L_{m-1}^{-1}$。眾所周知，高斯消元法極不穩定，用於有效數字較多的矩陣時誤差巨大。最簡單的解決方法是引入主元，從而得到更穩定的改進高斯消元法。

### 線性方程組的解

數值線性代數通常將矩陣視作列向量的組合。對於求解線性方程組$x = A^{-1}b$，傳統代數方法是將x理解為$A^{-1}$與b的積。數值線性代數則將x理解為b在A的列形成的基上現行展開的係數向量。

格努矩陣A及向量x、b的特徵，可以用很多不同的分解解決線性問題。若$A=QR$，則等價於$Rx = Q^\ast b$。這和矩陣分解一樣容易計算。將A進行特徵分解：$A = X \Lambda X^{-1}$，並試圖找到使$b=Ax$的b，且$b' = X^{-1}b,\ x' = X^{-1}x$，則有$b' = \Lambda x'$。這同使用奇異值分解求解線性系統密切相關，因為矩陣的奇異值是特徵值的絕對值，也等價于格拉姆矩陣$X^{*} X$特徵值絕對值的平方根。將A進行LU分解得$A=LU$，則$Ax=b$可用三角矩陣$Ly=b,\ Ux=y$求解。

### 最小二乘優化

矩陣分解提出了許多解線性方程組$r=b-Ax$的方法，其中我們尋求最小的r，如回歸問題。QR演算法計算A的既約QR分解並重排，得到$\widehat{R}x = \widehat{Q}^\ast b$。此上三角方程組可用於解x。SVD還提供了一種獲得線性最小二乘的演算法：計算既約SVD分解$A = \widehat{U}\widehat{\Sigma}V^\ast$，並計算向量$\widehat{U}^\ast b$，就可以將最小二乘問題簡化為簡單的對角方程組。QR分解和SVD分解可以產生最小二乘解，這意味著除了用經典正規方程法解最小二乘問題外，還可用格拉姆-施密特法、豪斯霍德法等方法。

## 條件和穩定性

假設問題是函數$f: X \to Y$，其中X是數據的賦范向量空間，Y是解的賦范向量空間。對於數據點$x \in X$，若x的微小擾動會導致$f(x)$的值發生很大變化，則稱之為病態的（ill-conditioned）。可以定義條件數量化之，代表問題的完備程度：

$$\widehat{\kappa} = \lim_{\delta \to 0} \sup_{\| \delta x \| \leq \delta} \frac{\| \delta f \|}{\| \delta x \|}.$$

不穩定性是指依賴浮點運算的電腦演算法的結果與問題的精確解產生差距的趨勢。矩陣包含多位有效數字的真實數據時，很多解線性方程組或最小二乘優化的算法會產生不準確的解。為病態問題創建穩定演算法是數值線性代數的核心問題，例如豪斯霍德三角化是線性方程組的很穩健的解法，而解最小二乘法的正規方程法不穩定，是奇異值分解更常用的一個原因。有些矩陣分解方法可能不穩定，但可以直接修改而變得穩定，例如格拉姆-施密特法可以很容易地改成改進格拉姆-施密特法；高斯消元法不穩定，引入主元後變得穩定。

## 迭代法

迭代法成為數值線性代數重要組成有兩個原因。其一，很多重要的數值問題沒有直接解，如求任意矩陣的特徵值和特徵向量。其二，對任意$m \times m$矩陣的非迭代演算法需要$O(m^3)$時間，而矩陣只含$m^2$個數字，這個下限非常高。迭代法可利用矩陣的特點縮短時間。例如對稀疏矩陣，迭代演算法可跳過很多直接方法必須的步驟，即便它們在矩陣高度結構化的情形下是多餘的。

數值線性代數中，很多迭代法的核心是將矩陣投影到低維克雷洛夫子空間，這樣就可從低維空間開始迭代計算類似矩陣的等效特徵，依次向高維移動以逼近高維矩陣的特徵。對對稱的A，解線性方程$Ax=b$時，經典迭代法是共軛梯度法。A若不對稱，則迭代法有廣義最小殘量方法和CGN（正規方程共軛梯度法）。若A對稱，則解特徵值、特徵向量問題時，可以用蘭佐斯算法；A不對稱時，可以用阿諾德迭代法。

## 軟體

有幾種程式語言用數值線性代數優化技術，旨在實現數值線性代數演算法。有MATLAB、Analytica、Maple和Mathematica。其他語言也有提供數值線性代數例程與優化的庫，如C語言及Fortran有BLAS、LAPACK等庫，Python有NumPy，Perl有Perl Data Language。R語言中的很多數值線性代數命令依賴於LAPACK這些更基本的庫。

## 參見

- 數值線性代數是數值分析的子領域。
- 高斯消去法，數值線性代數中一種重要的演算法。
- BLAS和LAPACK，高度優化的電腦程式庫，以Fortran程式語言寫成，可以實現數值線性代數中最基本的演算法。
- 數值分析軟體列表
- 數值程序庫列表
- 高斯-賽德爾疊代
- 雅可比疊代

## 閱讀更多

- J. H. Wilkinson and C. Reinsch, "Linear Algebra, volume II of Handbook for Automatic Computation" SIAM Review 14, 658 (1972).
- Golub, Gene H.; van Loan, Charles F. (1996), Matrix Computations, 3rd edition, Johns Hopkins University Press, ISBN 978-0-8018-5414-9

## 外部連結

- Freely available software for numerical algebra on the web , composed by Jack Dongarra and Hatem Ltaief, University of Tennessee
- NAG Library of numerical linear algebra routines
- nla_grp (Research group in the United Kingdom)
- siagla (Activity group about numerical linear algebra in the Society for Industrial and Applied Mathematics)
- The GAMM Activity Group on Applied and Numerical Linear Algebra