利用 Cayley-Hamilton 定理計算矩陣函數

本文的閱讀等級：中級

達文西說：「簡單是最終極的細緻。」(Simplicity is the ultimate sophistication.) 就數學而言，簡單不意味平庸，反而是優雅的體現。本文從一個簡單的問題開始，並致力於發展簡單的解法。令 $A=\begin{bmatrix} 1&1\\ 0&2 \end{bmatrix}$ ，求 $A^{100}$ 。典型的問題通常有典型的解法，對角化 (diagonalization) 是目前最常用的冪矩陣算法。矩陣 $A$ 有相異特徵值 $1$ 和 $2$ ，對應特徵向量 $\begin{bmatrix} 1\\ 0 \end{bmatrix}$ 和 $\begin{bmatrix} 1\\ 1 \end{bmatrix}$ ，故可對角化為

$A=\begin{bmatrix} 1&1\\ 0&1 \end{bmatrix}\begin{bmatrix} 1&0\\ 0&2 \end{bmatrix}\begin{bmatrix} 1&1\\ 0&1 \end{bmatrix}^{-1}=S\Lambda S^{-1}$ 。

利用上式立得

$A^{100}=(S\Lambda S^{-1})^{100}=S\Lambda^{100}S^{-1}=\begin{bmatrix} 1&1\\ 0&1 \end{bmatrix}\begin{bmatrix} 1&0\\ 0&2^{100} \end{bmatrix}\begin{bmatrix} 1&1\\ 0&1 \end{bmatrix}^{-1}=\begin{bmatrix} 1&2^{100}-1\\ 0&2^{100} \end{bmatrix}$ 。

不過，遺憾的是對角化並不適用於所有的矩陣。若 $A=\begin{bmatrix} 1&2\\ 0&1 \end{bmatrix}$ ，則 $A$ 有兩個特徵值 $1$ ，但其特徵空間僅含一線性獨立向量 $\begin{bmatrix} 1\\ 0 \end{bmatrix}$ ，即知 $A$ 不可對角化 (見“可對角化矩陣與缺陷矩陣的判定”)。針對不可對角化矩陣，典型的方法是分解出 $A$ 的 Jordan 形式 $A=MJM^{-1}$ ，接著計算 $A^{100}=MJ^{100}M^{-1}$ (見“利用 Jordan form 解差分方程與微分方程”)。表面上，計算冪矩陣並不很困難，使用 Jordan 典型形式似乎有些「小題大作」。倘若不採用 Jordan 分解，那麼還有其他方法嗎？繼續閱讀前，建議讀者先花個幾分鐘想一想。

最簡單的辦法是修改 $A$ 矩陣使之可對角化。我們可以在 $A$ 的主對角元加入擾動量，譬如 $A_{\epsilon}=\begin{bmatrix} 1&2\\ 0&1+\epsilon \end{bmatrix}$ ，其中 $\epsilon$ 是一極小的數。如此一來， $A_{\epsilon}$ 有相異特徵值 $1$ 和 $1+\epsilon$ ，因此可對角化為

$\displaystyle A_{\epsilon}=\begin{bmatrix} 1&2\\ 0&1+\epsilon \end{bmatrix}=\begin{bmatrix} 1&2\\ 0&\epsilon \end{bmatrix}\begin{bmatrix} 1&0\\ 0&1+\epsilon \end{bmatrix}\frac{1}{\epsilon}\begin{bmatrix} \epsilon &-2\\ 0&1 \end{bmatrix}$ 。

按照先前方式計算冪矩陣，如下：

$\displaystyle\begin{aligned} A_{\epsilon}^{100}&=\begin{bmatrix} 1&2\\ 0&\epsilon \end{bmatrix}\begin{bmatrix} 1&0\\ 0&(1+\epsilon)^{100} \end{bmatrix}\frac{1}{\epsilon}\begin{bmatrix} \epsilon &-2\\ 0&1 \end{bmatrix}\\ &=\frac{1}{\epsilon}\begin{bmatrix} 1&2\\ 0&\epsilon \end{bmatrix}\begin{bmatrix} 1&0\\ 0&1+100\epsilon+o(\epsilon^2) \end{bmatrix}\begin{bmatrix} \epsilon &-2\\ 0&1 \end{bmatrix}\\ &=\frac{1}{\epsilon}\begin{bmatrix} 1&2+200\epsilon+o(\epsilon^2)\\ 0&\epsilon+100\epsilon^2+o(\epsilon^3) \end{bmatrix}\begin{bmatrix} \epsilon &-2\\ 0&1 \end{bmatrix}\\ &=\frac{1}{\epsilon}\begin{bmatrix} \epsilon&200\epsilon+o(\epsilon^2)\\ 0&\epsilon+100\epsilon^2+o(\epsilon^3) \end{bmatrix}\\ &=\begin{bmatrix} 1&200+o(\epsilon)\\ 0&1+100\epsilon+o(\epsilon^2) \end{bmatrix},\end{aligned}$

上式中 $o(\epsilon^k)$ 代表最小次為 $k$ 的多項式。令 $\epsilon\to 0$ ，則 $A_{\epsilon}^{100}\to A^{100}=\begin{bmatrix} 1&200\\ 0&1 \end{bmatrix}$ 。這個方法的主要缺點在於引入代數符號運算。當 $A$ 的尺寸增大時，計算過程變得尤其複雜。例如， $A=\begin{bmatrix} 1&2&0\\ 0&1&2\\ 0&0&1 \end{bmatrix}$ ，令 $\epsilon_1$ 和 $\epsilon_2$ 代表兩個相異極小數，不難想像 $A_{\epsilon_1,\epsilon_2}=\begin{bmatrix} 1&2&0\\ 0&1+\epsilon_1&2\\ 0&0&1+\epsilon_2 \end{bmatrix}$ 有相當繁瑣的對角化展開式。

除了對角化，還有別的方法嗎？欲跳脫思路的陷阱，我們應暫時遺忘可對角化與不可對角化之間的對立，並將注意力集中在問題本身──冪矩陣。仔細觀察後發現 $A=\begin{bmatrix} 1&2\\ 0&1 \end{bmatrix}$ 滿足 $(A-I)^2=0$ 。這個事實稱為 Cayley-Hamilton 定理，它指出：對於任一 $n\times n$ 階矩陣 $A$ ，令 $p(t)=\det(A-tI)$ 為 $A$ 的特徵多項式 (或 $p(t)=\det(tI-A)$ 亦可)，則 $p(A)=0$ ，也就是說，矩陣 $A$ 被其特徵多項式消滅 (見“Cayley-Hamilton 定理”)。此例 $A$ 的特徵多項式是

$p(t)=\begin{vmatrix} 1-t&2\\ 0&1-t \end{vmatrix}=(1-t)^2=t^2-2t+1$ ，

就有 $p(A)=(A-I)^2=A^2-2A+I=0$ 。利用恆等式 $A^2=2A-I$ ，我們可以計算更高次的冪矩陣：

$\begin{aligned} A^3&=A(A^2)=A(2A-I)=2A^2-A=2(2A-I)-A=3A-2I,\\ A^4&=A(A^3)=A(3A-2I)=3A^2-2A=3(2A-I)-2A=4A-3I,\\ \vdots \end{aligned}$

重複此過程，任一冪矩陣 $A^k$ 皆可以表示成 $A$ 和 $I$ 的線性組合。從另一個角度來看，設想 $A^{100}$ 除以 $(A-I)^2$ 的表達式：

$A^{100}=q(A)(A-I)^2+c_1A+c_0I$ ，

其中 $q(A)$ 代表商， $c_1$ 和 $c_0$ 是未定係數。因為 $(A-I)^2=0$ ，可知 $A^{100}=c_1A+c_0I$ 。接下來要解出 $c_1$ 和 $c_0$ 。考慮多項式除法：

$t^{100}=q(t)(t-1)^2+c_1t+c_0$ 。

為了從 $t^{100}$ 的表達式得到兩個條件式 (因為有兩個未知係數)，求上式等號兩邊的導數，即有

$100t^{99}=q^{\prime}(t)(t-1)^2+2q(t)(t-1)+c_1=\left(q^{\prime}(t)(t-1)+2q(t)\right)(t-1)+c_1$ 。

將 $t=1$ 代入上面兩式可得

$c_1+c_0=1,~~c_1=100$ ，

解出 $c_1=100$ ， $c_0=-99$ ，所以 $A^{100}=100A-99I=\begin{bmatrix} 1&200\\ 0&1 \end{bmatrix}$ 。

除了計算冪矩陣，Cayley-Hamilton 定理是否也可以用來計算一般的矩陣函數 $f(A)$ ？(關於矩陣函數的定義，請見“矩陣函數 (上)”和“矩陣函數 (下)”。) 答案是肯定的。利用泰勒展開式和 Cayley-Hamilton 定理 $p(A)=0$ ，任何 $n\times n$ 階矩陣 (不論是否可對角化) $A$ 的良好定義函數 $f(A)$ 都可表示成矩陣多項式：

$f(A)=q(A)p(A)+c_{n-1}A^{n-1}+\cdots+c_1A+c_0I=c_{n-1}A^{n-1}+\cdots+c_1A+c_0I$ 。

設 $A$ 有相異特徵值 $\lambda_1,\ldots,\lambda_m$ ，且 $\beta_1,\ldots,\beta_m$ 為對應的相重數 (亦稱代數重數)，則 $A$ 的特徵多項式為

$p(t)=(t-\lambda_1)^{\beta_1}\cdots(t-\lambda_m)^{\beta_m}$ 。

考慮

$\begin{aligned} f(t)&=q(t)p(t)+c_{n-1}t^{n-1}+\cdots+c_1t+c_0\\ &=q(t)(t-\lambda_1)^{\beta_1}\cdots(t-\lambda_m)^{\beta_m}+c_{n-1}t^{n-1}+\cdots+c_1t+c_0, \end{aligned}$

對於特徵值 $\lambda_i$ ，我們可以使用 $\beta_i$ 個限制條件： $f(\lambda_i),~f^{\prime}(\lambda_i),\ldots,~f^{(\beta_i-1)}(\lambda_i)$ 。因為 $\sum_{i=1}^m\beta_i=n$ ，所以共有 $n$ 個條件式可以用來決定係數 $c_0,c_1,\ldots,c_{n-1}$ ，以矩陣形式表示如下：

$\begin{bmatrix} 1&\lambda_1&\lambda_1^2&\lambda_1^3&\cdots&\lambda_1^{n-1}\\ \vdots&\vdots&\vdots&\vdots&&\vdots\\ 1&\lambda_m&\lambda_m^2&\lambda_m^3&\cdots&\lambda_m^{n-1}\\ -&-&-&-&-&-\\ \vdots&\vdots&\vdots&\vdots&&\vdots\\ 0&1&2\lambda_i&3\lambda_i^2&\cdots&(n-1)\lambda_i^{n-2}\\ \vdots&\vdots&\vdots&\vdots&&\vdots\\ -&-&-&-&-&-\\ \vdots&\vdots&\vdots&\vdots&&\vdots\\ 0&0&2&6\lambda_i&\cdots&(n-1)(n-2)\lambda_i^{n-3}\\ \vdots&\vdots&\vdots&\vdots&&\vdots\\ -&-&-&-&-&-\\ \vdots&\vdots&\vdots&\vdots&&\vdots \end{bmatrix}\begin{bmatrix} c_0\\ c_1\\ c_2\\ c_3\\ \vdots\\ \vdots\\ \vdots\\ c_{n-1} \end{bmatrix}=\begin{bmatrix} f(\lambda_1)\\ \vdots\\ f(\lambda_m)\\ -\\ \vdots\\ f^{\prime}(\lambda_i)\\ \vdots\\ -\\ \vdots\\ f^{\prime\prime}(\lambda_i)\\ \vdots\\ -\\ \vdots \end{bmatrix}$ ，

其中由特徵值構成的 $n\times n$ 階係數矩陣是可逆的，因此保證 $c_0,c_1,\ldots,c_{n-1}$ 唯一存在，也就有唯一的 $f(A)$ 。係數矩陣為可逆矩陣的理由如下：(一) 最上面分塊包含相異特徵值的 $m$ 個列屬於 Vandermonde 矩陣 $V$ 的一部分，因此是線性獨立的 (見“特殊矩陣 (8)：Vandermonde 矩陣”)；(二) 隨後的各分塊皆可表示為 $VD$ ，其中 $D$ 是可逆對角矩陣，因此也是線性獨立的；(三) 各分塊的領先行位置相異，如第一分塊起始於第一行，第二分塊起始於第二行，故各分塊彼此獨立。

下面以一個例子展示計算步驟。考慮 $A=\begin{bmatrix} 1&2&0\\ 0&1&0\\ 0&0&3 \end{bmatrix}$ ，求 $e^{A}$ 。

(1) 求出特徵多項式 $p(t)=(t-1)^2(t-3)$ ，可知 $A$ 有特徵值 $1, 1, 3$ 。

(2) 令 $e^A=c_2A^2+c_1A+c_0I$ 。寫出函數 $f(t)=e^t$ 的多項式除法表達式：

$f(t)=e^{t}=q(t)(t-1)^2(t-3)+c_2t^2+c_1t+c_0$ 。

求 $f(t)$ 的導數，

$f^{\prime}(t)=e^t=g(t)(t-1)+2c_2t+c_1$ ，

其中 $g(t)$ 是整理合併後的多項式。

(3) 將 $t=1$ 和 $t=3$ 代入 $f(t)$ ，並將 $t=1$ 代入 $f^{\prime}(t)$ ，可得限制條件式：

$\begin{aligned} e&=c_2+c_1+c_0\\ e^3&=9c_2+3c_1+c_0\\ e&=2c_2+c_1, \end{aligned}$

或表示為矩陣形式：

$\begin{bmatrix} 1&1&1\\ 1&3&9\\ 0&1&2 \end{bmatrix}\begin{bmatrix} c_0\\ c_1\\ c_2 \end{bmatrix}=\begin{bmatrix} e\\ e^3\\ e \end{bmatrix}$ 。

由上式解出

$\displaystyle c_0=\frac{1}{4}e^3-\frac{3}{4}e,~~c_1=-\frac{1}{2}e^3+\frac{5}{2}e,~~c_2=\frac{1}{4}e^3-\frac{3}{4}e$ ，

也就得到

$\displaystyle\begin{aligned} e^A&=\left(\frac{1}{4}e^3-\frac{3}{4}e\right)A^2+\left(-\frac{1}{2}e^3+\frac{5}{2}e\right)A+\left(\frac{1}{4}e^3-\frac{3}{4}e\right)I\\ &=\left(\frac{1}{4}e^3-\frac{3}{4}e\right)\begin{bmatrix} 1&4&0\\ 0&1&0\\ 0&0&9 \end{bmatrix}+\left(-\frac{1}{2}e^3+\frac{5}{2}e\right)\begin{bmatrix} 1&2&0\\ 0&1&0\\ 0&0&3 \end{bmatrix}+\left(\frac{1}{4}e^3-\frac{3}{4}e\right)\begin{bmatrix} 1&0&0\\ 0&1&0\\ 0&0&1 \end{bmatrix}\\ &=\begin{bmatrix} e&2e&0\\ 0&e&0\\ 0&0&e^3 \end{bmatrix}.\end{aligned}$

相較於對角化 $f(A)=Sf(\Lambda)S^{-1}$ 或 Jordan 典型形式 $f(A)=Mf(J)M^{-1}$ ，利用 Cayley-Hamilton 定理計算矩陣函數 $f(A)$ 有三個明顯的優點：(一) 不須計算特徵向量或廣義特徵向量；(二) 不須計算特徵向量矩陣 $S$ 或 $M$ 的逆矩陣；(三) 不須計算 Jordan 矩陣 $J$ ，自然也就免除 $f(J)$ 的繁瑣算式 (見“矩陣函數 (下)”)。然而，此法必須使用一般矩陣乘法算出所有的冪矩陣 $A^k$ ， $k<n$ 。提醒讀者，二個 $n\times n$ 階矩陣相乘需要 $O(n^3)$ 個乘法運算，因此乘法運算總量是 $O(n^4)$ 。換句話說，本文介紹的這個方法僅適用於小尺寸矩陣。

1 Response to 利用 Cayley-Hamilton 定理計算矩陣函數

胡霖 says:

03/21/2018 at 3:34 pm

请问周老师良好定义函数是怎样的函数呢？

	陳倍恩 on 線性代數的第一堂課──矩陣乘法的定義
	輕鬆談如何教學二項式定理？… on 牛頓的二項式定理 (上)
	madhouse on 高斯消去法
	WishMobile on 翻轉 LU 分解
	周子傑 on Karush-Kuhn-Tucker (KKT) 條件
	Cloud Huang on 線性泛函與伴隨