1. 引言 : 亥維賽的妙招

眾所周知, 非齊次的常係數線性微分方程往往可以通過待定係數法求解, 這是通常教科書上介紹的方法, 如 。

記得從前我在給大一學生上微積分講到微分方程時, 常有愛動腦筋的學生問我 : 「為什麼要 用待定係數法求解, 這個方法從天而降, 讓我感覺自己很蠢。數學太深奧了！」 確實, 作為老師, 我講待定係數法也很難受。數學是講道理的, 我按照書上的講當然沒錯, 可是好像也說不出什麼道理。我都不記得當時我是怎麼回答他的！ 直到有一天我看到亥維賽(Oliver Heaviside, 1850$\sim$1925)的一個妙招, 我才覺得終於得救了, 下一次講微分方程時我一定要分享給學生！

我是從彭羅斯 (Rogers Penrose) 的一本書 [493-494] 中學到這個妙招的, 摘引如下:

亥維賽的洞見是, 微分運算元通常可以 如普通數一樣處理, 這個事實對求解某些類型的微分方程非常有用。 我們來看一個例子, 考慮微分方程 \begin{equation} y+\frac{\mathrm{d}^2y}{\mathrm{d}x^2}=x^5.\label{1} \end{equation} 我們想要求出一個特解。亥維賽的方法是, 像對待普通數那樣對待微分運算元 $\mathrm{d}/\mathrm{d}x$。 為讓它看起來更合「情理」, 我們用一個單獨的字母 $D$ 來表示這個運算元 \begin{equation} D=\frac{\mathrm{d}}{\mathrm{d}x}.\label{2} \end{equation} 兩次作用 $D$ 所得到的平方運算元 $D^2$ 表示 兩次微分, 即二階求導運算元 $\mathrm{d}^2/\mathrm{d}x^2$； $D^3$ 表示 三階求導運算元 $\mathrm{d}^3/\mathrm{d}x^3$；依次類推, 於是我們的方程成為 $y+D^2y=x^5$, 我們可以表示為 \begin{equation} (1+D^2)y=x^5.\label{3} \end{equation} 我們可以通過「除以係數 $1+D^2$」而 形式地「解」方程, 從而得到解 $y\!=\!(1\!+\!D^2)^{-1}x^5$。 將 $(1+D^2)^{-1}$ 展開成「$D$ 的冪級數」 : \begin{equation} (1+D^2)^{-1}=1-D^2+D^4-D^6+\cdots,\label{4} \end{equation} 從而我們求出一個 (正確的！) 特解 : \begin{align} y=&(1+D^2)^{-1}x^5\nonumber\\ =&(1-D^2+D^4-D^6+\cdots)x^5\nonumber\\ =&(1-D^2+D^4)x^5\qquad\qquad(\text{$D^6 x^5=0$, 等等})\label{5}\\ =&x^5-(x^5)"+(x^5)^{""}\nonumber\\ =&x^5-20x^3+120x.\nonumber \end{align}

如果注意恰當的規則, 那麼可以使得這個形式步驟無懈可擊 --- 不過, 亥維賽在首次運用它時卻遭到了強烈的反對。

2. 尋找關鍵點

據加州大學伯克利大學的數學教授伍鴻熙講 (), 他從數學大師陳省身那裡學到的最重要的一點, 就是

從整體而非局部地去把握問題, 不考慮 複雜的技術細節, 即抓住關鍵點。 $\cdots$ 他 (陳省身先生) 對待問題的人生哲學, 就是每件事情都有一個 關鍵點。如果你抓住了這個關鍵點, 那麼剩下的遲早會解決。

回到上述亥維賽解方程的巧妙方法, 成功的關鍵, 又在哪裡呢？

一個初步的分析指出, 成功的關鍵在於 : 方程 \eqref{1} 右邊的非齊次項 $x^5$ 是一個多項式, 而微分運算元 $D$ 的高次冪作用在上面等於0, 因此 \eqref{4} 右邊給出的 $D$ 的級數作用於 $x^5$ 後, 約化為 一個 有限截斷 (即 $1-D^2+D^4$) 的作用, 從而變無限為有限, 化分析為代數！

3. 隱藏的本質

可以設想, 亥維賽的方法當初之所以會遭到反對, 最主要的一點, 是無窮和 \eqref{4} 不好接受。 ^{1 1} 因此, 通常的微分方程教科書並不介紹這個妙招。 經典的北大教材 第一版 倒是介紹了, 但是作為加星號的一節 (「運算元法和拉普拉斯變換法簡介」)。也許是反響不好, 在第二版中, 這一節又被刪掉了！ 正如彭羅斯所指出的, 這是可以嚴格化的。 ^{2 2} 對此, 包括維納 (Norbert Wiener) 在內的許多數學家都作出了貢獻, 其關鍵字是~operational calculus, 有興趣的讀者可以進一步 瞭解。 但實際上, 有一個更簡單的方法, 完全不涉及無窮和！

我們只要觀察到 : $$(1+D^2)(1-D^2+D^4)=1+D^6.$$ 這意味著, 在次數小於 $6$ 的複係數多項式函數空間 ${\mathbb{C}}_6[x]$ 上 (有 $D^6=0$), 其實有 運算元乘積 $$(1+D^2)(1-D^2+D^4)=1+D^6=1,$$ 即, $(1-D^2+D^4)$ 是 $(1+D^2)$ 的右逆, 為簡單起見, 我們仍然記為 $(1+D^2)^{-1}$！ 從而, 對任意的 $f(x)\in {\mathbb C}_6[x]$, 可以直接寫出方程 $(1+D^2)y=f(x)$ 的解為 $$y=(1+D^2)^{-1}f(x)=(1-D^2+D^4)f(x).$$ 特別的, 對 $f(x)=x^5$, 我們即得到原來方程 $(1+D^2)y=x^5$ 的解 $$y=(1-D^2+D^4)x^5=x^5-20x^3+120x.$$

現在我們看出, 為了求解方程 $(1+D^2)y=x^5$, 我們真正需要的, 只是求出 $1+D^2$ 在包含 $x^5$ 的某個函數空間 (這裡就是 ${\mathbb{C}}_6[x]$) 上的右逆！ 亥維賽的方法, 相當於通過對 $(1+D^2)^{-1}$ 的形式冪級數 \eqref{4} 取有限截斷而得到這個右逆。 但通過進一步的分析, 我們發現有更簡單的方法。

4. 化微分方程為代數方程

為了說明這個更一般的方法, 我們不妨假設所考慮的, 是 一個更一般的方程 \begin{equation} P(D)y=x^5,\label{6} \end{equation} 其中 $P=P(x)$ 是一個複係數多項式。 根據前面的分析容易看出, 為了求出方程 \eqref{6} 的一個特解, 我們只要求出 $P(D)$ 在 ${\mathbb{C}}_6[x]$ 上的一個右逆。 由於 $D$ 在 ${\mathbb{C}}_6[x]$ 上的極小多項式是 $x^6$, 所以我們只要求出滿足同餘條件 ^{3 3} 這個條件等價於 $x^6\mid P(x)U(x)-1$, 即 $x^6$ 整除 $P(x)U(x)-1$。 \begin{equation} P(x)U(x)\equiv 1 \pmod{x^6}\label{7} \end{equation} 的一個複係數多項式 $U=U(x)$, 即可得到 $P(D)$ 在 $\mathbb{C}_6[x]$ 上的一個右逆 $U(D)$。 這是因為, 在 $P(x)U(x)-1=x^6Q(x)$ 中令 $x=D$ 即可得到 $$P(D)U(D)-I=D^6Q(D)=0\cdot Q(D)=0,$$ 即 $P(D)U(D)=I$。

因此, 我們要做的, 就是求解多項式同餘方程 \eqref{7}。 幸運的是, 這種類型的問題是古人早已解決的。 特別是, 當我們用整數代替多項式時, 相應的問題即整數的同餘方程 ^{4 4} 這個條件等價於 $b\mid ax-1$, 即 $b$ 整除 $ax-1$。 \begin{equation} ax\equiv 1 \pmod{b},\label{8} \end{equation} 其中 $a,b$ 是給定的正整數, $x$ 是要求的整數。

為幫助讀者更好地掌握多項式同餘方程 \eqref{7} 的解法, 我們下面先回顧一下整數同餘方程 \eqref{8} 的解法。

5. 解整數同餘方程的求一術

對整數同餘方程 \eqref{8}, 我們有經典的解法, 而且至少有三種等價的表述 : 分別是歐幾裡得演算法、連分數法與求一術。 這裡我們介紹的是第三種 : 求一術。

求一術是北宋數學家秦九韶 (1208$\sim$1268) 的發明, 他命名為大衍求一術 (至於「大衍」的意思, 在《數書九章》序中, 秦九韶 把這一方法與《周易》「大衍之數」附會), 是所謂「大衍總數術」的關鍵一步。 清代數學家黃宗憲後來進一步簡化了秦九韶的方法, 我們現在介紹的, 就是這個簡化的 版本。

秦九韶--黃宗憲的方法可用矩陣表述。首先寫出一個 2 行 2 列的陣 $$\begin{bmatrix} ~a~&~1~\\ b&0 \end{bmatrix},$$ 其中第一列 $a,b$ 都是源自方程 \eqref{8}, 而第二列的兩個元素 $1,0$ 則是添加進來的 (用於探測未知數 $x$)。

求一術演算法如下 (見 或 ) : 對第一列的數 $a,b$ 用帶餘除法 (較大的數除以較小的數), 設得到的商為 $q$, 則較大的數那一行減去較小數的那一行對應元素的 $q$ 倍； 於是新得到的矩陣的第一列兩個元素替換為第一次帶餘除法的除數與餘數, 重複之前的操作, 直到某一步帶餘除法得到的餘數為 1 (演算法結束, 用紅色標記), 此時 $1$ 的正 右方的數, 即為所求的 $x$ (用藍色示意)。

這個演算法的理論基礎可見最後一節的練習 1。

作為例子, 我們用秦九韶--黃宗憲的方法來求 \begin{equation} 5x\equiv1\pmod 7\label{9} \end{equation} 的一個解。

解 : 求一術步驟如下 : \begin{align*} &\begin{bmatrix} ~5~&\,\,~1~\,\,\\ 7&\,\,0\,\, \end{bmatrix} \xrightarrow[ 下行減去上行的{\color{blue}1}倍]{ 7={\color{blue}1}\cdot 5+2} \begin{bmatrix} ~5~&1\\ 2&-1 \end{bmatrix} \xrightarrow[ 上行減去下行的{\color{blue}2}倍]{5={\color{blue}2}\cdot 2+{\color{red}1}} \begin{bmatrix} {\color{red}1}&{\color{blue}3}\\ ~2~&-1 \end{bmatrix} \end{align*} 根據求一術, $x$ 在 $1$ 的右邊, 即 $x=3$。這是很容易驗證的 : $$5\cdot3=15\equiv1\pmod 7.$$

當然, 你或許以為我是把問題搞複雜了, 你甚至在一開始就試出來 $x=3$ 是一個解。 然而, 正如吳文俊先生多次強調的 (例如, 見 ), 中國古代數學講的是一種演算法。 簡單的例子你用技巧可以解決, 但如果換成一個稍微複雜的例子, 如解方程 $$250x\equiv1\pmod {2017},$$ 你可能就無計可施了！ ^{5 5} 這讓我們回想起著名數學家、 數學教育家波利亞 (George P\'{o}lya) 的一句名言 (「傳統的數學教授」一節) : 「方法與技巧的差別何在? 方法乃是放之四海而皆準的技巧。 (What is the difference between method and device? A method is a device which you used twice.)」 中國古代數學, 注重的是方法 (「法」「術」同義) 而非技巧。 讀者若想領教秦九韶--黃宗憲的求一術的威力, 不妨用上面的方程 $250x\equiv1\pmod {2017}$ 一試！

可以看出, 上述求一術的基礎是 :

整數的帶餘除法 : 設 $a$ 和 $b$ 是兩個整數, 其中 $b\gt 0$, 則存在唯一的整數 $q$ 和 $r$ 使得 \begin{equation*} a=qb+r, \end{equation*} 其中 $r$ 滿足 $0\leq\,r\lt b$。

6. 解多項式同餘方程的求一術

注意到, 類似的, 對多項式我們也有帶餘除法 :

多項式的帶餘除法 : 設 $a(x)$ 和 $b(x)$ 是兩個多項式, 其中 $b(x)\neq0$, 則存在唯一的多項式 $q(x)$ 和 $r(x)$ 使得 \begin{equation*} a(x)=q(x)b(x)+r(x), \end{equation*} 其中 $r(x)$ 滿足 $\deg r(x)\lt \deg b(x)$, 這裡 $\deg r(x),\deg b(x)$ 分別表示多項式 $r(x),b(x)$ 的次數(degree, 縮寫為 $\deg$).

因此, 同樣有求解一般多項式同餘方程 (假定其中 $a(x),b(x)$ 互質, 這是方程有解的充分必要條件) \begin{equation} a(x)u(x)\equiv1\pmod {b(x)}\label{10} \end{equation} 的求一術 :

首先寫出一個 2 行 2 列的陣 $$\begin{bmatrix} a(x)&1\\ b(x)&~0~ \end{bmatrix}.$$

對第一列的兩個多項式 $a(x),b(x)$ 用帶餘除法 (次數高的多項式除以次數低的), 設得到的商為 $q(x)$, 則次數高的那一行減去次數低的那一行對應元素的 $q(x)$ 倍； 於是新得到的矩陣的第一列兩個元素替換為第一次帶餘除法的除式與餘式, 重複之前的操作, 直到某一步帶餘除法得到的餘式為某個非零常數 $c$ (演算法結束), 此時 $c$ 的正右方的多項式除以 $c$ (「歸一化」), 即為所求的 $u(x)$。

評注 : 這裡用非零常數 $c$ 取代了整數同餘方程求一術中的目標"1", 可以這麼理解, 非零常數是多項式環中的單位 (可逆元), 而 $\pm1$ 也是整數環中的單位。 (在整數的情況, 由於我們限定了餘數大於等於0, 所以把 $-1$ 的情況就排除了。 如果我們限定餘數是關於 $0$ 最對稱的一組剩餘系, 求一術就成為了 「求正負一術」。) 所以求一術的目標「得一」(因而又名「得一術」), 實則是「得單位」。 因此, 最恰當的稱謂既不是「求一術」, 也不是「得一術」, 而是「求逆術」。 不論是整數還是多項式的情況, 我們所討論的方程 \eqref{8} \eqref{9} 實際上就是 求一個給定元素 (在某個商環中) 的逆。

作為例子, 我們來看一開頭彭羅斯所給出的微分方程 \eqref{3} 所確定的多項式同餘方程 : \begin{equation*} (x^2+1)u(x)\equiv1\pmod {x^6}. \end{equation*}

解 : 求一術步驟如下 : \begin{align*} \begin{bmatrix} x^2+1&\,\,~1~\,\,\\ x^6&\,\,0\,\, \end{bmatrix} \xrightarrow[ 下行減去上行的{\color{blue}{(x^4-x^2+1)}}倍]{x^6={\color{blue}{(x^4-x^2+1)}}\cdot (x^2+1)-1} \begin{bmatrix} x^2+1&~1~\\ {\color{red}{-1}}&{\color{blue}{-(x^4-x^2+1)}} \end{bmatrix}. \end{align*} 根據求一術, $u(x)$ 等於 $-1$ 的右邊的多項式除以 $-1$, 即 \begin{equation*} u(x)=\frac{-(x^4-x^2+1)}{-1}=x^4-x^2+1; \end{equation*} 這是很容易驗證的 : $$(x^2+1)\cdot(x^4-x^2+1)=x^6+1\equiv1\pmod {x^6}.$$ 這樣我們就求得了 $D^2+1$ 在 ${\mathbb{C}}_6[x]$ 上的逆為 $u(D)=D^4-D^2+1$, 這與我們前面用無窮冪級數取有限截斷得到的結果一致。 但從方法上講, 解同餘方程的方法更切中要害 (打蛇打七寸！), 從而也更容易理解。

7. 解常係數線性微分方程的新方法 : 非齊次項為多項式的情形

通過將微分方程轉化為多項式的同餘方程, 現在我們可以得到求解常係數線性微分方程 \begin{equation} P(D)y=f(x),\label{11} \end{equation} (其中 $P,f$ 是任意的複係數多項式) 的一個特解的一個新方法。

設 $f$ 的次數為 $m$, 則 $f\in \mathbb{C}_{m+1}[x]$ (其中 $\mathbb{C}_{m+1}[x]$ 表示次數小於 $m+1$ 的複係數多項式空間), 且 $D$ 在 $\mathbb{C}_{m+1}[x]$ 上 的極小多項式為 $x^{m+1}$。 因此(11)所對應的代數同餘方程為 \begin{equation} P(x)U(x)\equiv1\pmod {x^{m+1}}\label{12} \end{equation} 當且僅當 $P(x)$ 與 $x^{m+1}$ 互質------顯然, 這等價於 $P(0)\neq0$, 即 $P$ 的常數項非零------時, 我們可以求出 滿足方程(12)的 $U(x)$, 從而 $U(D)$ 是 $P(D)$ 的一個右逆 (實際上是真正的逆), 由此立即得到方程(11)的一個特解 $y=U(D)f(x)$, 注意, 這個解是多項式, 而且次數與 $f$ 的次數相等, 即為 $m$。

當 $P(0)=0$ 時怎麼辦呢？我們先看一個最特殊的例子, $P$ 以 $0$ 為唯一的根, 如 $P(x)=x^k$, 從而我們要解的方程即 \begin{equation*} D^ky=f(x). \end{equation*} 很明顯, 通過逐次積分, 可以得到這個方程的通解。

在一般情況, 若 $P(0)=0$, 我們對 $P(x)$ 作分解 $$P(x)=x^kQ(x),$$ 其中 $Q(0)\neq0$。 從而原方程(11)的一個特解, 可以通過以下兩步得到 : 第一, 用前述方法求出方程 $Q(D)z=f(x)$ 的一個特解 $z$, 這是一個次數等於 $m$ 的多項式； 第二, 求出方程 $D^ky=z$ 的一個特解 $y$, 通過每一次積分取平凡常數 (即 $0$), 我們 可以保證 $y=x^kw(x)$, 其中 $w(x)$ 是一個 $m$ 次多項式。事實上, 若寫 \begin{eqnarray*} z&=&a_m\frac{x^m}{m!}+a_{m-1}\frac{x^{m-1}}{(m-1)!}+\cdots+a_0,\\ {\hbox{則}} y&=&a_m\frac{x^{k+m}}{(k+m)!}+a_{m-1}\frac{x^{k+m-1}}{(k+m-1)!}+\cdots+a_0\frac{x^k}{k!}\\ &=&x^k\Big(a_m\frac{x^{m}}{(k+m)!}+a_{m-1}\frac{x^{m-1}}{(k+m-1)!}+\cdots+a_0\frac{1}{k!}\Big) \end{eqnarray*} 是方程 $D^ky=z$ 的一個特解。

評注 : 現在完全可以理解, 為什麼在待定係數法中, 我們可以假設具有所述形式的特解了。 然而, 說明這一點絕非我們的本意。 我們的目標, 說得宏大一些, 就是要把待定係數法取而代之！ 當然, 目前我們只解決了非齊次項為多項式的情形, 接下來我們會表明, 非齊次項為擬多項式的情形, 同樣可以「除之而後快」。

8. 解常係數線性微分方程的新方法 : 非齊次項為擬多項式的情形

現在我們轉向待定係數法所適用的更一般方程, 即我們來討論方程 \begin{equation} Ly=P(D)y=\mathrm{e}^{\lambda x} f(x),\label{13} \end{equation} 其中 $P,f$ 是任意的複係數多項式, $\lambda$ 是一個給定的 (複) 常數。 形如 $\mathrm{e}^{\lambda x} f(x)$ (其中 $f(x)$ 為多項式) 的函數, 稱為擬多項式。 當 $\lambda=0$ 時, 我們得到真正的多項式 $f(x)$。

我們的目標, 仍然是求方程 \eqref{13} 的一個特解。 注意, 當 $\lambda=0$ 時, 即回到上一節方程 \eqref{11}, 這是我們已經用新方法解決了的。 對於一般的 $\lambda$, 方程 \eqref{13} 怎麼解呢？

注意到, 形如 $\mathrm{e}^{\lambda x} g(x)$ (其中 $g(x)$ 為次數小於 $m$ 的複係數多項式) 的函數構成一個維數為 $m$ 的複向量空間 $V_m(\lambda)$, 而且它恰好可以看成 $V_m(0)=\mathbb{C}_m[x]$ (它顯然是 $m$ 維的) 在左乘變換 $g\mapsto \mathrm{e}^{\lambda x}g$ 之下的像, 而且 $V_m(\lambda)$ 也是微分運算元 $D$ 的不變子空間。 這個觀察允許我們將方程 \eqref{13} 的求解轉化為方程 \begin{equation} L_\lambda y=f(x),\label{14} \end{equation} 其中 $L_\lambda$ 是使得下述圖 ^{6 6} 更具體些, 其基礎是

可交換 (即 $\mathrm{e}^{\lambda x}L_\lambda=L \mathrm{e}^{\lambda x}$) 的唯一運算元, 即 $$L_\lambda=\mathrm{e}^{-\lambda x}L \mathrm{e}^{\lambda x}.$$ 將 $L=P(D)$ 代入並經過一個神奇的計算 (見最末一節的練習2), 上式可以化簡為 \begin{equation} L_\lambda=\mathrm{e}^{-\lambda x}L \mathrm{e}^{\lambda x}=\mathrm{e}^{-\lambda x}P(D) \mathrm{e}^{\lambda x}=P(D+\lambda).\label{15} \end{equation} 從而, 經過變數替換 $y=\mathrm{e}^{\lambda x}z$, 關於 $y$ 的方程(13)等價於關於 $z$ 的方程 \begin{equation} P(D+\lambda)z=f.\label{16} \end{equation} 這是一個我們在前一節就已經解決的問題。

評注1 : 不必囉嗦的是, 上述結果可以解釋為什麼在這一情形中待定係數法具有所述的特解形式。 也許有讀者要問, 此處的方法是否可以推廣到其它情形, 即非齊次項是具有其它形式的函數。斯特朗(Gilbert Strang) 的下述結果 表明, 本質上, 這個方法對擬多項式已經達到極致了。 換言之, 這個方法僅僅適用於擬多項式的線性組合。

定理1 : 設 $D$ 是作用在複值光滑函數空間 $C^\infty(\mathbb{R})$ 上的微分運算元, 則 $D$ 的任意一個有限維不變子空間具有形式 $$V=V_{m_1}(\lambda_1)\oplus\cdots\oplus V_{m_s}(\lambda_s),$$ 其中 $\lambda_1,\ldots,\lambda_s$ 是互異的 (複) 常數。

評注2 : 還有一種求解非齊次線性方程特解的線性代數方法, 即所謂的「零化子方法」 (annihilator method), 例如參見 第 8.14 節。 零化子方法可視為待定係數法的理論依據。 但不論是零化子方法還是待定係數法, 都不如這裡介绍的方法直接。 尤其是, 不論從理解還是算法的角度看, 此處介绍的方法都是最簡練的。

9. 一個應用 : 廣義特徵方程的求解

作為上一節結果的一個非凡應用, 我們來求解齊次方程 \begin{equation} (D-\lambda)^my=0,\label{17} \end{equation} 其中 $m$ 是正整數。

當 $m=1$ 時, 我們得到微分運算元 $D$ 的特徵方程, 即 \begin{equation} Dy=\lambda y.\label{18} \end{equation} 根據微積分的基本結果, 其通解為 $y=C\mathrm{e}^{\lambda x}$, 其中 $C$ 為任意常數。 (據此可以理解, 在微積分中何以指數函數如此重要, 因為它是微分運算元的特徵函數！)

對一般的 $m$, \eqref{17} 可謂 $D$ 的廣義特徵方程。 利用數學歸納法, 我們不難證明, 廣義特徵方程 \eqref{17} 的通解 (可謂之「廣義特徵函數」) 恰好是擬多項式空間, 即 $V_m(\lambda)$。 換言之, 我們有下述結果 :

定理2 : $(D-\lambda)^{m}y=0$ 的通解為 $$y=\mathrm{e}^{\lambda x}(C_1x^{m-1}+C_2x^{m-2}+\cdots+C_m),$$ 其中 $C_1,\ldots,C_m$ 為任意常數。

證明 : 用歸納法。$m=1$ 的情況, 已經考慮。

對 $m=2$, 我們由 $m=1$ 的結果首先推出 $$(D-\lambda)y=C_1\mathrm{e}^{\lambda x},$$ 從而得到非齊次的微分方程。 為求解上述方程, 引入 $y=\mathrm{e}^{\lambda x}z$, 從而上述方程變成 $$Dz=C_1,$$ 解得 $z=C_1x+C_2$, 從而 $y=\mathrm{e}^{\lambda x}(C_1x+C_2)$。

現在假設 $m=k$ 情況已證, 考慮 $m=k+1$ 的情況。 若 $y$ 滿足 $(D-\lambda)^{k+1} y=0$, 則 $(D-\lambda)y=\mathrm{e}^{\lambda x}(C_1x^{k-1}+C_2x^{k-2}+\cdots+C_k),$ 引入 $y=\mathrm{e}^{\lambda x}z$, 從而上述方程變成 $$Dz=C_1x^{k-1}+C_2x^{k-2}+\cdots+C_k,$$ 從而 $z=C_1'x^k+C_2'x^{k-1}+\cdots+C_k'+C_{k+1}',$ 進而 $y=\mathrm{e}^{\lambda x}(C_1'x^k+C_2'x^{k-1}+\cdots+C_k'+C_{k+1}')$。證畢！

評注 : 根據這一結果, 再結合線性代數中的下述基本結果 (例如, 見 [定理5.2.1]), 我們可以給出常係數線性齊次微分方程 $P(D)y=0$ 的完整討論, 但限於篇幅, 此處從略。

定理3 : 設 $A$ 是複向量空間 $V$ 的線性變換, 設 $f(x)$ 是複係數多項式, 假定 $f(x)=(x-\lambda_1)^{m_1}\cdots(x-\lambda_s)^{m_s}$, 則有以下對應的零化子空間分解 $$\ker \big(f(A)\big)=\ker \big((A-\lambda_1)^{m_1}\big)\oplus\cdots\oplus \ker\big((A-\lambda_s)^{m_s}\big),$$ 其中 $\ker B=\{v\in V, Bv=0\}$ 表示 $B$ 的零化子空間 (又稱核空間)。

10. 解常係數線性遞推關係的新方法

不難想見, 在離散情形, 我們可以 得到求解常係數線性遞推關係的新方法。

10.1. 齊次情形

設 $\ell=\{(x_0,x_1,\ldots,x_n,\ldots),\,x_i\in\mathbb{C}\}$ 為所有的複數值數列的集合。

令 $T$ 表示作用在 $\ell$ 上的右平移運算元 : \begin{equation} T\,:\, x=(x_0,x_1,\ldots,x_n,\ldots)\mapsto (x_1,\ldots,x_n,\ldots)=Tx.\label{19} \end{equation} 於是, 一個 $d$ 階常係數齊次遞推關係 \begin{equation}x_{n}=a_1x_{n-1}+a_2x_{n-2}+\cdots+a_dx_{n-d}=0,\qquad (n\geq d)\label{20} \end{equation} 可以寫成下述形式 : \begin{equation} P(T)x_n=0,\qquad (n\geq d)\label{21} \end{equation} 其中 $P(T)=T^d-a_1T^{d-1}-\cdots-a_d$。

根據定理 3, 為了求解方程 \eqref{21}, 假定我們有特徵多項式分解 $$p(t)=t^d-a_1t^{d-1}-\cdots-a_d=(t-\lambda_1)^{m_1}\cdots(t-\lambda_s)^{m_s},$$ 我們只需要求解對應的廣義特徵方程, 即 \begin{equation} (T-\lambda)^mx_n=0,\qquad (n\geq m)\label{22} \end{equation} 對此, 我們有下述平行于定理2的結果 :

定理$2'$ : 設 $\lambda$ 為常數, 則 $(T-\lambda)^{m}x_n=0\,(n\geq m)$ 的通解為 \begin{equation*} x_n=\lambda^n(C_1n^{m-1}+C_2n^{m-2}+\cdots+C_m),\qquad(n\geq m)\qquad\qquad \tag{*} \end{equation*} 其中 $C_1,\ldots,C_m$ 為任意常數。

證明 : 若 $\lambda=0$, 則此時我們要解的方程為 $T^mx=0$。 根據定義容易看出, $T^mx=0$ 當且僅當 $x_{m}=x_{m+1}=\cdots=0$。

若 $\lambda=1$, $T-\lambda=T-1=\Delta$ 是向前差分運算元 (即 $\Delta x_n=x_{n+1}-x_n$), 上述方程化簡為 $$\Delta^m x=0.$$ 根據朱世傑招差公式 (參見[定理3]), 我們推出 $x_n=x(n)$ 是一個次數不超過 $m-1$ 的多項式。

若 $\lambda\neq0,1$, 我們構造 $\ell$ 上的線性變換 $K_\lambda$ 如下 : $$K_\lambda:\,\,(x_0,x_1,\ldots,x_n,\ldots)\,\mapsto\,(x_0,\lambda^{-1}x_1,\ldots,\lambda^{-n}x_n,\ldots).$$ 那麼容易驗證, 我們有 $$K_\lambda(T-\lambda)=\lambda(T-1)K_\lambda,$$ 即 $$(T-\lambda)=\lambda K_\lambda^{-1}(T-1)K_\lambda=K_\lambda^{-1}(\lambda \Delta)K_\lambda,$$ 從而 $$(T-\lambda)^m=(K_\lambda^{-1}(\lambda \Delta)K_\lambda)^m=K_\lambda^{-1}(\lambda \Delta)^mK_\lambda.$$ 由此立即看出, $$(T-\lambda)^{m}x=0\iff (\lambda \Delta)^m y=0\iff \Delta^m y=0,$$ 其中 $y=K_\lambda x$。 由於 $ \Delta^m y=0$ 的解為次數不超過 $m-1$ 的多項式, 所以立即得出 $x=K_\lambda^{-1} y=K_{\lambda^{-1}} y$ 具有定理所述形式。

評注 : 此處我們的證明與定理2中的證明略有不同, 因為我們此時還未介紹非齊次遞推關係的求解。 待下面介紹這一方法後, 讀者可以自行補充一個完全平行的證明。 這同時建議讀者去思考定理 2 的一個平行的證明。

至此, 常係數齊次遞推關係的求解得到了圓滿的解決。

10.2. 非齊次情形

現在我們轉向非齊次情形, 即我們要求解方程 \begin{equation} P(T)x(n)=\lambda^nf(n),\qquad (n\geq d)\label{23} \end{equation} 其中 $\lambda$ 是一個非零常數, $P,f$ 是複係數多項式, 且 $P$ 的次數為 $d$。

先考慮 $\lambda=1$ 的特殊情況, 即方程 \begin{equation} P(T)x(n)=f(n),\qquad (n\geq d)\label{24} \end{equation} 其中 $P,f$ 是多項式, 次數分別為 $d,m$。

根據上一節的分析, 我們只要求出 $f(n)$ 生成的一個 $T$-不變子空間。由於 $T=\Delta+1$, 不難確定 $f(n)$ 的一個 $\Delta$-不變子空間從而也是 $T$-不變子空間不變子空間為 $${\mathbb{C}}_{m+1}=\{(x_0,\ldots,x_n,\ldots)\,|\,x_n=g(n),\,g\ \text{為次數不超過 $m$ 的複係數多項式}\},$$ 即以次數不超過 $m$ 的多項式為通項公式的數列。 由於差分運算元 $\Delta$ 在 $\mathbb{C}_{m+1}$ 上的極小多項式為 $t^{m+1}$, 所以右平移運算元 $T=\Delta+1$ 在 $\mathbb{C}_{m+1}$ 上的極小多項式為 $(t-1)^{m+1}$。 因此對方程 \eqref{24} 的討論分為以下兩種情況 :

1. 若 $P(1)\neq0$, 這意味著多項式 $P(t)$ 與 $(t-1)^{m+1}$ 互質, 從而用求一術可以解出一個多項式 $U(t)$, 使得 $$P(t)U(t)\equiv1 \pmod{(t-1)^{m+1}}.$$ 於是 \eqref{24} 的解為 $x=U(T)f$。 注意到, 對任意的正整數 $k$ 有, $T^kf(n)=f(n+k)$, 由此不難看出, 解 $x(n)=U(T)f(n)$ 是一個 $m$ 次多項式。
2. 不然有 $P(t)=(t-1)^kQ(t)$, 其中 $Q(1)\neq0$。 從而 $P(T)=(T-1)^kQ(T)=\Delta^k Q(T)$, 於是我們可以分兩步求解 \eqref{22}。 第一步, 用前述方法求出方程 $Q(T)z=f$ 的一個特解 $z=z(n)$, 這是一個 $m$ 次多項式； 第二步, 對第一步得到的 $z$, 求出方程 $\Delta^kx=z$ 的一個特解 $x$, 我們 可以保證 $x(n) =[n(n-1)\cdots(n-k+1)]w(n)$, 其中 $w(n)$ 是一個 $m$ 次多項式。事實上, 根據 反復應用楊輝恒等式的差分形式(參見 [p. 71] 第 \eqref{16} 式), 我們有 : 若 $z(n)$ 的朱世傑招差公式(參見 [p.67] 定理3) 為: \begin{eqnarray*} z(n)&=&a_m{n\choose m}+a_{m-1}{n\choose {m-1}}+\cdots+a_0,\\ {\hbox{則}} x(n)&=&a_m{n\choose {k+m}}+a_{m-1}{n\choose {k+m-1}}+\cdots+a_0{n\choose k}\\ &=&[n(n-1)\cdots(n-k+1)]\cdot \sum_{s=0}^m\left[\frac{a_s s!}{(k+s)!}{{n-k}\choose {s}}\right] \end{eqnarray*} 是方程 $\Delta^kx=z$ 的一個特解。

   

現在我們轉向一般情形 \eqref{23} 的討論。 一個簡單的分析指出 (這裡也有一個類似的交換圖, 如上), 借助變數替換 $x_n=\lambda^ny_n$, 方程 \eqref{21} 等價於關於 $y$ 的方程 \begin{equation} P(\lambda T)y=f,\label{25} \end{equation} 而這個方程我們在上面已經解決了！

10.3. 評注

評注1 : 不難想到, 我們也有一個平行于定理1的類似結果, 此處從略。 這樣的結果限定了我們這個方法的適用範圍, 即只能應用于非齊次項為形如(*)的函數的線性組合的情形。

評注2 : 在通常的教科書 (如 ) 中, 對於常係數線性遞推關係的求解, 作者一般優先介紹母函數法, 當然也有其它方法, 如華羅庚先生在 第四章講了包括母函數法在內的四個方法, 但正如華先生指出的 : 「雖然我們講了四個方法, 但實質上並沒有太多的差異, 講來講去繞不過分項分數 (注 : 即我們所謂的部分分式分解) 這一關。」

當然, 也有作者為了繞開「部分分式分解」, 就採取「待定係數法」, 如高德納(Donald E. Knuth)及其合作者就在書中 說[p.17] :

對處理所有的常係數線性遞推關係, 部分分式分解完全能夠勝任。 然而, 為簡單起見, 我們將介紹一個不同的方法, 這在許多舊一點的文獻中都可以找見。 這個方法是基於試探解, 類似於微分方程的求解。在某些例子中, 這個方法能夠最快提供答案, 但其規則看起來像黑魔法, 為了理解這種「掐指一算」為何奏效, 疑惑的讀者終究要回到作為其基礎的部分分式理論。}

我們要指出, 事實上本文介紹的新方法, 不僅可以避開待定係數法和母函數法, 甚至可以幫助我們重新理解作為其公共基礎的部分分式理論。限於篇幅, 此處我們不再展開。

10.4. 應用舉例

例1 (取自中譯本[p.189]) : 求遞推關係 $x_n=5x_{n-1}-6x_{n-2}+7^n,\,(n\geq2)$ 的一個特解。

解 : 令 $T$ 為右平移運算元, 則上述方程可寫成 (注意, 右邊非齊次項多出一個 $7^2$) $$P(T)x_n=7^27^n,$$ 其中 $P(T)=T^2-5T+6$。 令 $x_n=7^nz_n$, 則上述方程等價於 $$Q(T)z_n=7^2=49,$$ 其中 $Q(T)=P(7T)=49T^2-35T+6$。 由於 $Q(1)=49-35+6=20\neq0$, 所以 我們只要求出一個多項式 $U(T)$ 使得 $$Q(T)U(T)\equiv1\mod{(T-1)}.$$ 當然我們可以用求一術, 但此處有更簡單的方法。 因為上述條件等價於 $$Q(1)U(1)-1=0,$$ 從而 $U(1)=\dfrac{1}{P(1)}=\dfrac{1}{20}$, 因此 $U(T)=\dfrac{1}{20}$ 是同餘方程的一個特解, 從而 $z_n=U(T)49=\dfrac{49}{20}$ 是 $Q(T)z_n=49$ 的一個特解, $x_n=\dfrac{49}{20}7^n$ 是原方程的一個特解。

例2 (取自[p.341]) : 求遞推關係 $x_n=x_{n-1}+2x_{n-2}+(-1)^n,\,(n\geq2)$ 的通解。

解 : 分三步。

第一步, 求出一個特解。 令 $T$ 為右平移運算元, 則上述方程等價於 $$P(T)x_n=(-1)^{n+2}=(-1)^n,$$ 其中 $P(T)=T^2-T-2=(T+1)(T-2)$。 令 $x_n=(-1)^nz_n$, 則上述方程等價於 $z_n$ 的下述方程 $$Q(T)z_n=1,$$ 其中 $Q(T)=P(-T)=(T-1)(T+2)=\Delta(T+2)$。 因此我們先求出 $y_n$ 使得 $(T+2)y_n=1$, 只要取常數列 $y_n=\dfrac{1}{3}$。 進一步求解 $\Delta z_n=y_n=\dfrac{1}{3}$, 即得 $Q(T)z_n=1$ 的一個特解 $z_n=\dfrac{1}{3}n$。 從而原方程的一個特解為 $x_n=\dfrac{1}{3}n(-1)^n$。

第二步, 求出對應齊次方程的通解。因為特徵方程為 $P(T)=T^2-T-2=(T+1)(T-2)$, 根據定理3, 我們只要分別求出 $(T+1)x_n=0$ 與 $(T-2)x_n=0$ 的通解, 並作疊加, 根據定理 $2'$, 結果為 $$C_1(-1)^n+C_22^n,\qquad\text{其中 $C_1,C_2$ 為任意複常數}.$$

第三步, 根據線性原理, 得到非齊次方程的通解公式為 : $$\frac{1}{3}n(-1)^n+C_1(-1)^n+C_22^n,\qquad\text{其中 $C_1,C_2$ 為任意複常數}.$$

11. 練習

1. 證明以下結果, 並說明這是求一術解整數同餘方程的理論依據 (提示 : 左乘即行變換)。 命題 : 設二階整數矩陣 $$A=\begin{bmatrix} a&1\\ ~b~&~0~ \end{bmatrix}.$$
   1. 若二階整數矩陣 $$B=\begin{bmatrix} u&\ast\\ ~\ast~&~\ast~ \end{bmatrix}$$ 使得 $$BA=\begin{bmatrix} 1&v\\ ~\ast~&~\ast~ \end{bmatrix},$$ 則 $u=v$, 且 $x=v$ 滿足同餘方程 $ax\equiv 1\pmod {b}$。
   2. 若二階整數矩陣 $$C=\begin{bmatrix} ~\ast~&~\ast~\\ u&\ast \end{bmatrix}$$ 使得 $$CA=\begin{bmatrix} ~\ast~&~\ast~\\ 1&v \end{bmatrix},$$ 則 $u=v$, 且 $x=v$ 滿足同餘方程 $ax\equiv 1\pmod {b}$。
2. 設 $D$ 是作用在光滑函數空間上的微分運算元, 證明, 對任意的多項式 $P$ 與常數 $\lambda$, 有 $\mathrm{e}^{-\lambda x}P(D)\mathrm{e}^{\lambda x}=P(D+\lambda)$。 並據此給出定理 2 的一個簡單證明。
3. 設 $T$ 是作用在數列空間上的右平移運算元, 證明, 對任意的多項式 $P$ 與非零常數 $\lambda$, 有 $\mathrm{\lambda}^{-n}P(T)\lambda^n =P(\lambda T)$。 請據此化簡定理$2'$的證明。
4. 利用本文介紹的方法, 求出斐波拉契數列 (Fibonacci sequence) $$x_n=x_{n-1}+x_{n-2}(n\geq 2);\qquad x_0=x_1=1$$ 的通項公式。 (讀者可以參見[p.36])

致謝 : 三年前, 筆者第一次從史特格茲(Steven Strogatz)的科普著作 中見到 從線性代數的觀點來求解斐波拉契數列的通項公式, 令人耳目一新。 在本文的寫作過程中, 筆者從與天津大學理學院劉雲朋教授、 中央民族大學理學院王兢教授、 麻省理工學院(MIT)數學系斯特朗(Gilbert Strang)教授、 威斯康辛大學密爾沃基分校許光午教授、 上海大學李常品教授、 勞倫斯伯克利國家實驗室邵美悅博士的討論中受益良多, 特表感謝！ 筆者對當初追問待定係數法道理何在的那些學生深表感謝, 在某種程度上, 這篇文章首先是筆者對他們的一個答覆。 感謝本刊審稿人對初稿提出了批評指正和進一步的修改意見。

參考文獻

---本文作者任教於中國西北農林科技大學理學院---