特殊矩陣 (19)：Hessenberg 矩陣

Posted on 04/02/2013 by ccjou

本文的閱讀等級：中級

令 A=[a_{ij}] 為一 n\times n 階矩陣。若 a_{ij}=0， i>j+1，則 A 稱為上 Hessenberg 矩陣，也就是說，A 的主對角下標元 (subdiagonal，即 a_{i+1,i}) 之下的所有元為零。若 A 的主對角上標元 (superdiagonal，即 a_{i,i+1}) 之上的所有元為零，則稱為下 Hessenberg 矩陣。此特殊矩陣因德國工程師黑森貝格 (Karl Adolf Hessenberg) 而得名。見下例，A 是上 Hessenberg 矩陣，B 是下 Hessenberg 矩陣，C 同時是上、下 Hessenberg 矩陣，稱為三對角 (tridiagonal) 矩陣 (見“特殊矩陣 (11)：三對角矩陣”)：

A=\begin{bmatrix} \ast&\ast&\ast&\ast&\ast\\ \ast&\ast&\ast&\ast&\ast\\ 0&\ast&\ast&\ast&\ast\\ 0&0&\ast&\ast&\ast\\ 0&0&0&\ast&\ast \end{bmatrix},~~B=\begin{bmatrix} \ast&\ast&0&0&0\\ \ast&\ast&\ast&0&0\\ \ast&\ast&\ast&\ast&0\\ \ast&\ast&\ast&\ast&\ast\\ \ast&\ast&\ast&\ast&\ast \end{bmatrix},~~C=\begin{bmatrix} \ast&\ast&0&0&0\\ \ast&\ast&\ast&0&0\\ 0&\ast&\ast&\ast&0\\ 0&0&\ast&\ast&\ast\\ 0&0&0&\ast&\ast \end{bmatrix}。

明顯地，對稱 Hessenberg 矩陣必定是三對角矩陣。下 Hessenberg 矩陣的轉置即為上 Hessenberg 矩陣，在不失一般性的原則下，以下討論限定於上 Hessenberg 矩陣。若上 Hessenberg 矩陣的所有主對角下標元皆不等於零，即 a_{i+1,i}\neq 0，i=1,\ldots,n-1，我們稱之為未約化 (unreduced) Hessenberg 矩陣，否則稱為約化 Hessenberg 矩陣。值得注意的是，約化 Hessenberg 矩陣的主對角分塊由未約化 Hessenberg 矩陣構成。上例中，若 a_{43}=0 且其他主對角下標元皆不為零，則

A=\begin{bmatrix} \ast&\ast&\ast&\vline&\ast&\ast\\ \ast&\ast&\ast&\vline&\ast&\ast\\ 0&\ast&\ast&\vline&\ast&\ast\\ \hline 0&0&0&\vline&\ast&\ast\\ 0&0&0&\vline&\ast&\ast \end{bmatrix}=\begin{bmatrix} A_{11}&A_{12}\\ 0&A_{22} \end{bmatrix}，

其中 A_{11} 和 A_{22} 都是未約化 Hessenberg 矩陣。下面我們介紹 Hessenberg 矩陣的一些性質，LU 分解和 QR 分解，以及一個相當有用的算法，稱為 Hessenberg 約化 (reduction)：利用 Householder 變換，任一實方陣皆正交相似於一 Hessenberg 矩陣。

 
矩陣乘積

兩個 Hessenberg 矩陣的乘積未必是 Hessenberg 矩陣，但 Hessenberg 矩陣與三角矩陣的乘積仍為 Hessenberg 矩陣。準確地說，若 A 為一上 Hessenberg 矩陣，T 為一上三角矩陣，則 AT 和 TA 都是上 Hessenberg 矩陣。

 
對角化

若 \lambda 是未約化 Hessenberg 矩陣 A 的特徵值，則 \lambda 的幾何重數，即 \dim N(A-\lambda I)，等於 1。證明於下。因為 A-\lambda I 不可逆，\mathrm{rank}(A-\lambda I)<n。再有，A-\lambda I 的前 n-1 個行向量包含非零的主對角下標元 a_{i+1,i}，這些行向量構成一線性獨立集，因此 \mathrm{rank}(A-\lambda I)\ge n-1。合併以上結果，\mathrm{rank}(A-\lambda I)=n-1。使用秩─零度定理，可得 \dim N(A-\lambda I)=n-\mathrm{rank}(A-\lambda I)=1。

 
如果未約化 Hessenberg 矩陣 A 包含相重的特徵值 \lambda，則 \lambda 的幾何重數必小於代數重數，故知 A 是不可對角化矩陣 (見“可對角化矩陣與缺陷矩陣的判定”)。從 Jordan 形式來說，令 J 是未約化 Hessenberg 矩陣 A 的 Jordan 形式：

J=J(\lambda_1)\oplus\cdots J(\lambda_k)，

其中 \lambda_1,\ldots,\lambda_k 是 A 的相異特徵值，J(\lambda_i) 表示對應 \lambda_i 的「超級」Jordan 分塊 (見“Jordan 形式大解讀 (上)”)。因為特徵值 \lambda_i 的幾何重數等於 1，故僅有一個基本 Jordan 分塊 J(\lambda_i)，換句話說，\lambda_i 的指標 (index，即最大基本 Jordan 分塊的階數) 等於代數重數。譬如，未約化 Hessenberg 矩陣 A 有特徵值 2,3,3,4,4,4，立知 A 的 Jordan 形式為

J=\begin{bmatrix} 2&\vline&0&0&\vline&0&0&0\\ \cline{1-4} 0&\vline&3&1&\vline&0&0&0\\ 0&\vline&0&3&\vline&0&0&0\\ \cline{3-8} 0&&0&0&\vline&4&1&0\\ 0&&0&0&\vline&0&4&1\\ 0&&0&0&\vline&0&0&4 \end{bmatrix}。

利用上述性質可推論實對稱、未約化三對角矩陣，形如

\begin{bmatrix} a_{1}&b_1&0&\cdots&0&0\\ b_{1}&a_{2}&b_{2}&\cdots&0&0\\ 0&b_{2}&a_{3}&\cdots&0&0\\ \vdots&\vdots&\vdots&\ddots&\vdots&\vdots\\ 0&0&0&\cdots&a_{n-1}&b_{n-1}\\ 0&0&0&\cdots&b_{n-1}&a_n \end{bmatrix},~~b_j\neq 0~~j=1,\ldots,n-1，

必有相異的特徵值，否則它是不可對角化矩陣，這與實對稱矩陣可正交對角化相矛盾 (見“實對稱矩陣可正交對角化的證明”)。

 
行列式

為方便說明，考慮 4\times 4 階上 Hessenberg 矩陣 A。令 A_k 代表 A 的 k\times k 階領先主子陣 (它們也都是 Hessenberg 矩陣)。利用餘因子公式計算行列式，針對第四行展開，可得

\begin{aligned} \begin{vmatrix} a_{11}&a_{12}&a_{13}&a_{14}\\ a_{21}&a_{22}&a_{23}&a_{24}\\ 0&a_{32}&a_{33}&a_{34}\\ 0&0&a_{43}&a_{44} \end{vmatrix}&=a_{44}\begin{vmatrix} a_{11}&a_{12}&a_{13}\\ a_{21}&a_{22}&a_{23}\\ 0&a_{32}&a_{33} \end{vmatrix}-a_{34} \begin{vmatrix} a_{11}&a_{12}&a_{13}\\ a_{21}&a_{22}&a_{23}\\ 0&0&a_{43} \end{vmatrix}\\ &+a_{24} \begin{vmatrix} a_{11}&a_{12}&a_{13}\\ 0&a_{32}&a_{33}\\ 0&0&a_{43} \end{vmatrix}-a_{14} \begin{vmatrix} a_{21}&a_{22}&a_{23}\\ 0&a_{32}&a_{33}\\ 0&0&a_{43} \end{vmatrix}\\ &=a_{44}\det A_3-a_{34}a_{43}\det A_2+a_{24}a_{43}a_{32}\det A_1-a_{14}a_{43}a_{32}a_{21}\det A_0, \end{aligned}

其中我們定義 \det A_0\equiv 1。推廣至一般情況，n\times n 階上 Hessenberg 矩陣 A_n 的行列式可用下列遞歸表達式計算[1]：

\det A_n=a_{nn}\det A_{n-1}+\displaystyle \sum_{i=1}^{n-1}\left((-1)^{n-i}a_{in}\prod_{j=i}^{n-1}a_{j+1,j}\det A_{i-1}\right)。

運用這個行列式公式 (或其他公式)，不難驗證 n\times n 階 Hessenberg 矩陣

F_n=\begin{bmatrix} 2&1&1&\cdots&1\\ 1&2&1&\cdots&1\\ 0&1&2&\cdots&1\\ \vdots&\vdots&\vdots&\ddots&\vdots\\ 0&0&0&\cdots&2 \end{bmatrix}

的行列式即是費布納西數列：\det F_n=f_{n+2}，其中 f_0=0，f_1=1，f_n=f_{n-1}+f_{n-2} (見“費布納西數列的表達式”)。

 
LU 分解

如果 Hessenberg 矩陣 A 可逆且所有的 k\times k 階領先主子陣也都可逆，則 A 有唯一的 LU 分解 (見“LU 分解”)：

\begin{bmatrix} a_{11}&a_{12}&a_{13}&\cdots&a_{1n}\\ a_{21}&a_{22}&a_{23}&\cdots&a_{2n}\\ 0&a_{32}&a_{33}&\cdots&a_{3n}\\ \vdots&\vdots&\vdots&\ddots&\vdots\\ 0&0&0&\cdots&a_{nn} \end{bmatrix}=\begin{bmatrix} 1&0&\cdots&0&0\\ l_1&1&\cdots&0&0\\ 0&l_2&\cdots&0&0\\ \vdots&\vdots&\ddots&\vdots&\vdots\\ 0&0&\cdots&l_{n-1}&1 \end{bmatrix}\begin{bmatrix} u_{11}&u_{12}&u_{13}&\cdots&u_{1n}\\ 0&u_{22}&u_{23}&\cdots&u_{2n}\\ 0&0&u_{33}&\cdots&u_{3n}\\ \vdots&\vdots&\vdots&\ddots&\vdots\\ 0&0&0&\cdots&u_{nn} \end{bmatrix}，

其中下三角矩陣 L 是上 Hessenberg 矩陣，U 是上三角矩陣 (主對角元 u_{ii} 皆不為零)。乘開上式，可解出 (1) u_{1j}=a_{1j},~~j=1,\ldots,n，(2) l_i=a_{i+1,i}/u_{ii},~~i=1,\ldots,n-1，(3) u_{i+1,j}=a_{i+1,j}-l_iu_{ij},~~i=1,\ldots,n-1,~~j=i+1,\ldots,n。相較一般 n\times n 階矩陣的 LU 分解的計算量 O(n^3)，Hessenberg 矩陣的計算量為 O(n^2)。

 
QR 分解

在一般情況下，Householder 變換應用於 QR 分解所耗費的運算量少於 Givens 旋轉 (見“QR 分解的數值計算方法比較”)，不過也有例外發生，Hessenberg 矩陣即為其一。原因在於 Hessenberg 矩陣已經具有許多零元，我們只要消去主對角下標元即可獲得 QR 分解的上三角矩陣 R。見下例，使用 Givens 旋轉 (見“Givens 旋轉於 QR 分解的應用”)，可得

Q_{21}A=\left[\begin{array}{crcc} c&-s&0&0\\ s&c&0&0\\ 0 & 0& 1&0 \\ 0 & 0& 0& 1 \end{array}\right]\begin{bmatrix} a_{11}&\ast&\ast&\ast\\ a_{21}&\ast&\ast&\ast\\ 0&\ast&\ast&\ast\\ 0&0&\ast&\ast \end{bmatrix}=\begin{bmatrix} \ast&\ast&\ast&\ast\\ 0&\ast&\ast&\ast\\ 0&\ast&\ast&\ast\\ 0&0&\ast&\ast \end{bmatrix}，

其中 Q_{21} 代表消去 (2,1) 元的 Givens 旋轉矩陣，且 c=a_{11}/\sqrt{a_{11}^2+a_{21}^2}，s=-a_{21}/\sqrt{a_{11}^2+a_{21}^2}。繼續對右下 3\times 3 階 Hessenberg 子陣執行消去運算可將主對角下標元 a_{i+1,i} 悉數消滅，過程如下：

\begin{bmatrix} \ast&\ast&\ast&\ast\\ \ast&\ast&\ast&\ast\\ 0&\ast&\ast&\ast\\ 0&0&\ast&\ast \end{bmatrix}\xrightarrow[]{~Q_{21}~}\begin{bmatrix} \ast&\ast&\ast&\ast\\ 0&\ast&\ast&\ast\\ 0&\ast&\ast&\ast\\ 0&0&\ast&\ast \end{bmatrix}\xrightarrow[]{~Q_{32}~}\begin{bmatrix} \ast&\ast&\ast&\ast\\ 0&\ast&\ast&\ast\\ 0&0&\ast&\ast\\ 0&0&\ast&\ast \end{bmatrix}\xrightarrow[]{~Q_{43}~}\begin{bmatrix} \ast&\ast&\ast&\ast\\ 0&\ast&\ast&\ast\\ 0&0&\ast&\ast\\ 0&0&0&\ast \end{bmatrix}

所以，Q_{43}Q_{32}Q_{21}A=R。旋轉矩陣是正交矩陣 (orthogonal matrix)，Q_{ij}^T=Q_{ij}^{-1}，正交矩陣乘積仍是正交矩陣，因此 A 有 QR 分解式 A=QR=(Q_{21}^TQ_{32}^TQ_{43}^T)R。如同 LU 分解，一般 n\times n 階矩陣的 QR 分解的計算量是 O(n^3)，Hessenberg 矩陣的計算量則為 O(n^2)。

 
Hessenberg 約化

如果方陣 A 相似於一三角矩陣 T，即存在一可逆矩陣 P 使得 P^{-1}AP=T，則 T 的主對角元就是 A 的特徵值。現實上，尋找使矩陣三角化的變換矩陣 P 是一個相當困難的工作。退而求其次，下面我們介紹如何運用 Householder 變換 (見“特殊矩陣 (4)：Householder 矩陣”) 求一正交矩陣 P 使得 P^TAP 為一 Hessenberg 矩陣，這裡 A 限定為實矩陣。

 
我們用 n=5 為例說明計算步驟。將 A 以分塊矩陣表示：

A=\begin{bmatrix} \ast&\vline&\ast&\ast&\ast&\ast\\ \hline \ast&\vline&\ast&\ast&\ast&\ast\\ \ast&\vline&\ast&\ast&\ast&\ast\\ \ast&\vline&\ast&\ast&\ast&\ast\\ \ast&\vline&\ast&\ast&\ast&\ast \end{bmatrix}=\begin{bmatrix} a_{11}&\mathbf{b}^T\\ \mathbf{c}&A_1 \end{bmatrix}。

步驟一：消去向量 \mathbf{c} 的末三元。設計主對角分塊正交矩陣 P_1=\begin{bmatrix} 1&0\\ 0&H_1 \end{bmatrix}，其中 4\times 4 階 Householder 矩陣 H_1 使得 H_1\mathbf{c}=\begin{bmatrix} \ast\\ 0\\ 0\\ 0 \end{bmatrix}。注意，H_1 滿足 H_1=H_1^T=H_1^{-1}，同樣地，P_1=P_1^T=P_1^{-1}，因此它們是對稱正交矩陣。所以，

P_1^TAP_1=\begin{bmatrix} 1&0\\ 0&H_1 \end{bmatrix}\begin{bmatrix} a_{11}&\mathbf{b}^T\\ \mathbf{c}&A_1 \end{bmatrix}\begin{bmatrix} 1&0\\ 0&H_1 \end{bmatrix}=\begin{bmatrix} a_{11}&\mathbf{b}^TH_1\\ H_1\mathbf{c}&H_1A_1H_1 \end{bmatrix}=\begin{bmatrix} \ast&\vline&\ast&\ast&\ast&\ast\\ \hline \ast&\vline&\ast&\ast&\ast&\ast\\ 0&\vline&\ast&\ast&\ast&\ast\\ 0&\vline&\ast&\ast&\ast&\ast\\ 0&\vline&\ast&\ast&\ast&\ast \end{bmatrix}。

步驟二：對 4\times 4 階子陣 A_2=H_1A_1H_1 重複上述程序。設計 P_2=\begin{bmatrix} I_2&0\\ 0&H_2 \end{bmatrix}，其中 3\times 3 階 Householder 矩陣 H_2 使得 \begin{bmatrix} 1&0\\ 0&H_2 \end{bmatrix} A_2\begin{bmatrix} 1&0\\ 0&H_2 \end{bmatrix}=\begin{bmatrix} \ast&\vline&\ast&\ast&\ast\\ \hline \ast&\vline&\ast&\ast&\ast\\ 0&\vline&\ast&\ast&\ast\\ 0&\vline&\ast&\ast&\ast \end{bmatrix}，故

P_2^TP_1^TAP_1P_2=\begin{bmatrix} \ast&\vline&\ast&\vline&\ast&\ast&\ast\\ \hline \ast&\vline&\ast&\vline&\ast&\ast&\ast\\ \hline 0&\vline&\ast&\vline&\ast&\ast&\ast\\ 0&\vline&0&\vline&\ast&\ast&\ast\\ 0&\vline&0&\vline&\ast&\ast&\ast \end{bmatrix}。

步驟三：接著繼續對右下 3\times 3 階子陣重複消去程序，便告結束。給定一 n\times n 階矩陣，執行 n-2 個步驟即可得上 Hessenberg 矩陣 P^TAP，其中 P=P_1P_2\cdots P_{n-2} 是正交矩陣，因為 P^{-1}=P_{n-2}^{-1}\cdots P_2^{-1}P_1^{-1}=P_{n-2}^T\cdots P_2^TP_1^T=(P_1P_2\cdots P_{n-2})^T=P^T。

 
如果對下 Hessenberg 矩陣 (P^TAP)^T=P^TA^TP 重複上述約化程序，可得到一個三對角矩陣 Q^TP^TA^TPQ，其中 PQ 是正交矩陣。換句話說，對於任一實方陣 A，存在一正交矩陣 U 使得 U^TAU 為一三對角矩陣。

 
表面上，Hessenberg 矩陣似乎不具備豐富的 (理論) 性質，但它是數值線性代數中極為重要的特殊矩陣。本文解說了如何利用 Householder 變換將給定實方陣化約成一上 Hessenberg 矩陣，以及運用 Givens 旋轉快速求出 Hessenberg 矩陣的 QR 分解。因為這兩個特性，Hessenberg 矩陣遂成為今日普遍採用的一種特徵值數值計算方法──QR 迭代法──的骨幹。日後我將另文介紹這個優異的特徵值算法。

 
參考來源：
[1] Nathan D. Cahill, John R. E’Errico, Darren A. Narayan, and Jack Y. Narayan, Fibonacci Determinants, The College Mathematics Journal, Vol. 33, 2002, pp. 221-225.

Advertisement
This entry was posted in 特殊矩陣, 線性代數專欄 and tagged , , , , , , , . Bookmark the permalink.

Leave a Reply

Fill in your details below or click an icon to log in:

Gravatar
WordPress.com Logo

You are commenting using your WordPress.com account. ( Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. ( Log Out /  Change )

Cancel

Connecting to %s