特殊矩陣 (19):Hessenberg 矩陣

本文的閱讀等級:中級

A=[a_{ij}] 為一 n\times n 階矩陣。若 a_{ij}=0i>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 0i=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 為一上三角矩陣,則 ATTA 都是上 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_kA 的相異特徵值,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 代表 Ak\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=0f_1=1f_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.

Advertisements
本篇發表於 特殊矩陣, 線性代數專欄 並標籤為 , , , , , , , 。將永久鏈結加入書籤。

發表迴響

在下方填入你的資料或按右方圖示以社群網站登入:

WordPress.com Logo

您的留言將使用 WordPress.com 帳號。 登出 / 變更 )

Twitter picture

您的留言將使用 Twitter 帳號。 登出 / 變更 )

Facebook照片

您的留言將使用 Facebook 帳號。 登出 / 變更 )

Google+ photo

您的留言將使用 Google+ 帳號。 登出 / 變更 )

連結到 %s