1. 引言 : 亥維賽的妙招

解常係數線性微分方程和遞推關係的 新方法一 一秦九韶和亥維賽的遺產

林開亮

1. 引言 : 亥維賽的妙招

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

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

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

亥維賽的洞見是, 微分運算元通常可以如普通數一樣處理, 這個事實對求解某些類型 的微分方程非常有用。 我們來看一個例子, 考慮微分方程

y + d²y

dx² = x⁵. (1)

我們想要求出一個特解。 亥維賽的方法是, 像對待普通數那樣對待微分運算元 d/dx。

為讓它看起來更合 「情理」, 我們用一個單獨的字母 D 來表示這個運算元 D = d

dx. (2)

兩次作用 D 所得到的平方運算元 D² 表示兩次微分, 即二階求導運算元 d²/dx²; D³ 表示三階求導運算元 d³/dx³; 依次類推, 於是我們的方程成為 y + D²y = x⁵, 我們可以表示為

(1 + D²)y = x⁵. (3) 我們可以通過 「除以係數 1 + D²」 而形式地 「解」 方程, 從而得到解 y = (1 +

63

D²)⁻¹x⁵。 將 (1 + D²)⁻¹ 展開成 「D 的冪級數」:

(1 + D²)⁻¹ = 1− D²+ D⁴− D⁶+· · · , (4) 從而我們求出一個 (正確的!) 特解 :

y =(1 + D²)⁻¹x⁵

=(1− D² + D⁴− D⁶+· · · )x⁵

=(1− D² + D⁴)x⁵ (D⁶x⁵ = 0, 等等) (5)

=x⁵− (x⁵)^+ (x⁵)^

=x⁵− 20x³+ 120x.

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

2. 尋找關鍵點

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

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

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

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

3. 隱藏的本質

可以設想, 亥維賽的方法當初之所以會遭到反對, 最主要的一點, 是無窮和 (4) 不好接 受。¹正如彭羅斯所指出的, 這是可以嚴格化的。² 但實際上, 有一個更簡單的方法, 完全不涉及 無窮和!

1因此, 通常的微分方程教科書並不介紹這個妙招。 經典的北大教材 [2] 第一版倒是介紹了, 但是作為加星號的一節 (「運算元法和拉普拉 斯變換法簡介」)。 也許是反響不好, 在第二版中, 這一節又被刪掉了!

2對此, 包括維納 (Norbert Wiener) 在內的許多數學家都作出了貢獻, 其關鍵字是 operational calculus, 有興趣的讀者可以進一 步瞭解。

我們只要觀察到 :

(1 + D²)(1− D²+ D⁴) = 1 + D⁶.

這意味著, 在次數小於 6 的複係數多項式函數空間 C6[x]上 (有 D⁶ = 0), 其實有運算元乘積 (1 + D²)(1− D² + D⁴) = 1 + D⁶ = 1,

即, (1 − D²+ D⁴) 是 (1 + D²) 的右逆, 為簡單起見, 我們仍然記為 (1 + D²)⁻¹! 從而, 對 任意的 f(x) ∈ C6[x], 可以直接寫出方程 (1 + D²)y = f (x) 的解為

y = (1 + D²)⁻¹f (x) = (1− D²+ D⁴)f (x).

特別的, 對 f(x) = x⁵, 我們即得到原來方程 (1 + D²)y = x⁵ 的解 y = (1− D²+ D⁴)x⁵ = x⁵− 20x³+ 120x.

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

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

為了說明這個更一般的方法, 我們不妨假設所考慮的, 是一個更一般的方程

P (D)y = x⁵, (6)

其中 P = P (x) 是一個複係數多項式。 根據前面的分析容易看出, 為了求出方程 (6) 的一個特 解, 我們只要求出 P (D) 在 C6[x]上的一個右逆。 由於 D 在 C6[x]上的極小多項式是 x⁶, 所 以我們只要求出滿足同餘條件³

P (x)U (x)≡ 1 (mod x⁶) (7) 的一個複係數多項式 U = U(x), 即可得到 P (D) 在 C6[x] 上的一個右逆 U(D)。 這是因為, 在 P (x)U(x) − 1 = x⁶Q(x) 中令 x = D 即可得到

P (D)U (D)− I = D⁶Q(D) = 0· Q(D) = 0, 即 P (D)U(D) = I。

3這個條件等價於 x⁶| P (x)U(x) − 1, 即 x⁶ 整除 P (x)U(x) − 1。

因此, 我們要做的, 就是求解多項式同餘方程 (7)。 幸運的是, 這種類型的問題是古人早已 解決的。 特別是, 當我們用整數代替多項式時, 相應的問題即整數的同餘方程⁴

ax≡ 1 (mod b), (8)

其中 a, b 是給定的正整數, x 是要求的整數。

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

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

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

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

秦九韶–黃宗憲的方法可用矩陣表述。 首先寫出一個 2 行 2 列的陣

! a 1 b 0

"

,

其中第一列 a, b 都是源自方程 (8), 而第二列的兩個元素 1, 0 則是添加進來的 (用於探測未知 數 x)。

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

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

作為例子, 我們用秦九韶–黃宗憲的方法來求

5x≡ 1 (mod 7) (9)

的一個解。

4這個條件等價於 b | ax − 1, 即 b 整除 ax − 1。

解 : 求一術步驟如下 :

! 5 1 7 0

"

7=1·5+2

−−−−−−−−−−→

下行減去上行的1倍

! 5 1 2 −1

"

5=2·2+1

−−−−−−−−−−→

上行減去下行的2倍

! 1 3 2 −1

"

根據求一術, x 在 1 的右邊, 即 x = 3。 這是很容易驗證的 : 5· 3 = 15 ≡ 1 (mod 7).

當然, 你或許以為我是把問題搞複雜了, 你甚至在一開始就試出來 x = 3 是一個解。 然而, 正如吳文俊先生多次強調的 (例如, 見 [16]), 中國古代數學講的是一種演算法。 簡單的例子你 用技巧可以解決, 但如果換成一個稍微複雜的例子, 如解方程

250x≡ 1 (mod 2017), 你可能就無計可施了!⁵

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

整數的帶餘除法 : 設 a 和 b 是兩個整數, 其中 b > 0, 則存在唯一的整數 q 和 r 使得 a = qb + r,

其中 r 滿足 0 ≤ r < b。

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

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

多項式的帶餘除法 : 設 a(x) 和 b(x) 是兩個多項式, 其中 b(x) = 0, 則存在唯一的多項式 q(x) 和 r(x) 使得

a(x) = q(x)b(x) + r(x),

其中 r(x) 滿足 deg r(x) < deg b(x), 這裡 deg r(x), deg b(x) 分別表示多項式 r(x), b(x) 的 次數(degree, 縮寫為 deg).

因此, 同樣有求解一般多項式同餘方程 (假定其中 a(x), b(x) 互質, 這是方程有解的充分 必要條件)

a(x)u(x)≡ 1 (mod b(x)) (10) 的求一術 :

5這讓我們回想起著名數學家、 數學教育家波利亞 (George P´olya) 的一句名言 ([8]「傳統的數學教授」 一節) : 「方法與技巧的差別 何在? 方法乃是放之四海而皆準的技巧。 (What is the difference between method and device? A method is a device which you used twice.)」 中國古代數學, 注重的是方法 (「法」「術」 同義) 而非技巧。 讀者若想領教秦九韶–黃宗憲的求一術的威力, 不妨用上面的方程 250x ≡ 1 (mod 2017) 一試!

首先寫出一個 2 行 2 列的陣

!a(x) 1 b(x) 0

"

.

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

評注 : 這裡用非零常數 c 取代了整數同餘方程求一術中的目標“1”, 可以這麼理解, 非零常數 是多項式環中的單位 (可逆元), 而 ±1 也是整數環中的單位。 (在整數的情況, 由於我們限定了 餘數大於等於 0, 所以把 −1 的情況就排除了 。 如果我們限定餘數是關於 0 最對稱的一組剩餘 系, 求一術就成為了 「求正負一術」。) 所以求一術的目標 「得一」 (因而又名 「得一術」), 實則是

「得單位」。 因此, 最恰當的稱謂既不是 「求一術」, 也不是 「得一術」, 而是 「求逆術」。 不論是整 數還是多項式的情況, 我們所討論的方程 (8) (9) 實際上就是求一個給定元素 (在某個商環中) 的逆。

作為例子, 我們來看一開頭彭羅斯所給出的微分方程 (3) 所確定的多項式同餘方程 : (x²+ 1)u(x)≡ 1 (mod x⁶).

解 : 求一術步驟如下 :

!x²+ 1 1 x⁶ 0

"

x⁶=(x⁴−x²+1)·(x²+1)−1

−−−−−−−−−−−−−−−−−→

下行減去上行的(x⁴− x²+ 1)倍

!x²+ 1 1

−1 −(x⁴− x²+ 1)

"

. 根據求一術, u(x) 等於 −1 的右邊的多項式除以 −1, 即

u(x) = −(x⁴− x²+ 1)

−1 = x⁴− x²+ 1;

這是很容易驗證的 :

(x²+ 1)· (x⁴− x²+ 1) = x⁶+ 1≡ 1 (mod x⁶).

這樣我們就求得了 D²+ 1 在 C6[x] 上的逆為 u(D) = D⁴− D²+ 1, 這與我們前面用無窮 冪級數取有限截斷得到的結果一致。 但從方法上講, 解同餘方程的方法更切中要害 (打蛇打七 寸!), 從而也更容易理解。

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

通過將微分方程轉化為多項式的同餘方程, 現在我們可以得到求解常係數線性微分方程

P (D)y = f (x), (11)

(其中 P, f 是任意的複係數多項式) 的一個特解的一個新方法。

設 f 的次數為 m, 則 f ∈ C_m+1[x] (其中 C_m+1[x] 表示次數小於 m + 1 的複係數多項 式空間), 且 D 在 C_m+1[x] 上的極小多項式為 x^m+1。 因此 (11) 所對應的代數同餘方程為

P (x)U (x)≡ 1 (mod x^m+1) (12) 當且僅當 P (x) 與 x^m+1 互質——顯然, 這等價於 P (0) = 0, 即 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, 從而我們要解的方程即

D^ky = f (x).

很明顯, 通過逐次積分, 可以得到這個方程的通解。

在一般情況, 若 P (0) = 0, 我們對 P (x) 作分解 P (x) = x^kQ(x),

其中 Q(0) = 0。 從而原方程 (11) 的一個特解, 可以通過以下兩步得到 : 第一, 用前述方法 求出方程 Q(D)z = f(x) 的一個特解 z, 這是一個次數等於 m 的多項式; 第二, 求出方程 D^ky = z 的一個特解 y, 通過每一次積分取平凡常數 (即 0), 我們可以保證 y = x^kw(x), 其 中 w(x) 是一個 m 次多項式。 事實上, 若寫

z = a_mx^m

m! + a_m−1 x^m−1

(m− 1)! +· · · + a0, 則

y = a_m x^k+m

(k + m)! + a_m−1 x^k+m−1

(k + m− 1)!+· · · + a0

x^k k!

= x^k 

a_m x^m

(k + m)! + a_m−1 x^m−1

(k + m− 1)!+· · · + a0

1 k!

是方程 D^ky = z 的一個特解。

評注 : 現在完全可以理解, 為什麼在待定係數法中, 我們可以假設具有所述形式的特解了。 然 而, 說明這一點絕非我們的本意。 我們的目標, 說得宏大一些, 就是要把待定係數法取而代之!

當然, 目前我們只解決了非齊次項為多項式的情形, 接下來我們會表明, 非齊次項為擬多項式的 情形, 同樣可以 「除之而後快」。

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

現在我們轉向待定係數法所適用的更一般方程, 即我們來討論方程

Ly = P (D)y = e^λxf (x), (13) 其中 P, f 是任意的複係數多項式, λ 是一個給定的 (複) 常數。 形如 e^λxf (x) (其中 f(x) 為 多項式) 的函數, 稱為擬多項式。 當 λ = 0 時, 我們得到真正的多項式 f(x)。

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

注意到, 形如 e^λxg(x) (其中 g(x) 為次數小於 m 的複係數多項式) 的函數構成一個維數 為 m 的複向量空間 Vm(λ), 而且它恰好可以看成 Vm(0) =Cm[x] (它顯然是 m 維的) 在左 乘變換 g → e^λxg 之下的像, 而且 Vm(λ) 也是微分運算元 D 的不變子空間。 這個觀察允許我 們將方程 (13) 的求解轉化為方程

L_λy = f (x), (14)

其中 Lλ 是使得下述圖⁶

V_m(0) ^{eλx //}

Lλ



V_m(λ)

L

V_m(0) ^{eλx //}V_m(λ) 可交換 (即 e^λxL_λ = Le^λx) 的唯一運算元, 即

L_λ = e^−λxLe^λx.

將 L = P (D) 代入並經過一個神奇的計算 (見最末一節的練習2), 上式可以化簡為

L_λ = e^−λxLe^λx = e^−λxP (D)e^λx= P (D + λ). (15) 從而, 經過變數替換 y = e^λxz, 關於 y 的方程 (13) 等價於關於 z 的方程

P (D + λ)z = f. (16)

6更具體些, 其基礎是

f ^eλx //

D+λ

e^λxf

D

(Df + λf) ^eλx//e^λx(λf + Df )

這是一個我們在前一節就已經解決的問題。

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

定理1 : 設 D 是作用在複值光滑函數空間 C^∞(R) 上的微分運算元, 則 D 的任意一個有限維 不變子空間具有形式

V = V_m1(λ₁)⊕ · · · ⊕ Vms(λ_s), 其中 λ1, . . . , λ_s 是互異的 (複) 常數。

評注2 : 還有一種求解非齊次線性方程特解的線性代數方法, 即所謂的 「零化子方法」

(annihilator method), 例如參見 [1] 第 8.14 節。 零化子方法可視為待定係數法的理論依據。

但不論是零化子方法還是待定係數法, 都不如這裡介绍的方法直接。 尤其是, 不論從理解還是算 法的角度看, 此處介绍的方法都是最簡練的。

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

作為上一節結果的一個非凡應用, 我們來求解齊次方程

(D− λ)^my = 0, (17) 其中 m 是正整數。

當 m = 1 時, 我們得到微分運算元 D 的特徵方程, 即

Dy = λy. (18)

根據微積分的基本結果, 其通解為 y = Ce^λx, 其中 C 為任意常數。 (據此可以理解, 在微積分 中何以指數函數如此重要, 因為它是微分運算元的特徵函數!)

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

定理2 : (D − λ)^my = 0的通解為

y = e^λx(C₁x^m−1+ C₂x^m−2 +· · · + Cm),

其中 C1, . . . , C_m 為任意常數。

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

對 m = 2, 我們由 m = 1 的結果首先推出

(D− λ)y = C1e^λx,

從而得到非齊次的微分方程。 為求解上述方程, 引入 y = e^λxz,從而上述方程變成 Dz = C₁,

解得 z = C1x + C₂,從而 y = e^λx(C₁x + C₂)。

現在假設 m = k 情況已證, 考慮 m = k + 1 的情況。 若 y 滿足 (D − λ)^k+1y = 0, 則 (D− λ)y = e^λx(C₁x^k−1+ C₂x^k−2+· · · + Ck),引入 y = e^λxz, 從而上述方程變成

Dz = C₁x^k−1+ C₂x^k−2+· · · + Ck,

從而 z = C1^x^k+C₂^x^k−1+· · ·+C_k^+C_k+1^ ,進而 y = e^λx(C₁^x^k+C₂^x^k−1+· · ·+C_k^+C_k+1^ )。 證畢!

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

定理3 : 設 A 是複向量空間 V 的線性變換, 設 f(x) 是複係數多項式, 假定 f(x) = (x − λ₁)^m1· · · (x − λs)^ms, 則有以下對應的零化子空間分解

ker f (A)

= ker

(A− λ1)^m1

⊕ · · · ⊕ ker

(A− λs)^ms , 其中 ker B = {v ∈ V, Bv = 0} 表示 B 的零化子空間 (又稱核空間)。

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

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

10.1. _齊次情形

設  = {(x0, x₁, . . . , x_n, . . .), x_i ∈ C} 為所有的複數值數列的集合。

令 T 表示作用在  上的右平移運算元 :

T : x = (x₀, x₁, . . . , x_n, . . .)→ (x1, . . . , x_n, . . .) = T x. (19)

於是, 一個 d 階常係數齊次遞推關係

x_n= a₁x_n−1+ a₂x_n−2+· · · + adx_n−d = 0, (n ≥ d) (20) 可以寫成下述形式 :

P (T )x_n = 0, (n≥ d) (21) 其中 P (T ) = T^d− a1T^d−1− · · · − ad。

根據定理 3, 為了求解方程 (21), 假定我們有特徵多項式分解

p(t) = t^d− a1t^d−1− · · · − ad= (t− λ1)^m1· · · (t − λs)^ms, 我們只需要求解對應的廣義特徵方程, 即

(T − λ)^mx_n = 0, (n≥ m) (22) 對此, 我們有下述平行于定理2的結果 :

定理2^ : 設 λ 為常數, 則 (T − λ)^mx_n= 0 (n≥ m) 的通解為

x_n= λⁿ(C₁n^m−1+ C₂n^m−2 +· · · + C_m), (n≥ m) (*) 其中 C1, . . . , C_m 為任意常數。

證明 : 若 λ = 0, 則此時我們要解的方程為 T^mx = 0。 根據定義容易看出, T^mx = 0 當且僅 當 xm = x_m+1 =· · · = 0。

若 λ = 1, T − λ = T − 1 = Δ 是向前差分運算元 (即 Δxn = x_n+1− xn), 上述方程 化簡為

Δ^mx = 0.

根據朱世傑招差公式 (參見[7, 定理3]), 我們推出 xn = x(n) 是一個次數不超過 m − 1 的多 項式。

若 λ = 0, 1, 我們構造  上的線性變換 Kλ 如下 :

K_λ : (x₀, x₁, . . . , x_n, . . .) → (x0, λ⁻¹x₁, . . . , λ⁻ⁿx_n, . . .).

那麼容易驗證, 我們有

K_λ(T − λ) = λ(T − 1)K_λ, 即

(T − λ) = λK_λ⁻¹(T − 1)Kλ = K_λ⁻¹(λΔ)K_λ,

從而

(T − λ)^m = (K_λ⁻¹(λΔ)K_λ)^m = K_λ⁻¹(λΔ)^mK_λ. 由此立即看出,

(T − λ)^mx = 0 ⇐⇒ (λΔ)^my = 0 ⇐⇒ Δ^my = 0,

其中 y = Kλx。 由於 Δ^my = 0 的解為次數不超過 m − 1 的多項式, 所以立即得出 x = K_λ⁻¹y = K_λ−1y 具有定理所述形式。

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

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

10.2. _{非齊次情形}

現在我們轉向非齊次情形, 即我們要求解方程

P (T )x(n) = λⁿf (n), (n≥ d) (23) 其中 λ 是一個非零常數, P, f 是複係數多項式, 且 P 的次數為 d。

先考慮 λ = 1 的特殊情況, 即方程

P (T )x(n) = f (n), (n≥ d) (24) 其中 P, f 是多項式, 次數分別為 d, m。

根據上一節的分析, 我們只要求出 f(n) 生成的一個 T -不變子空間。 由於 T = Δ + 1, 不 難確定 f(n) 的一個 Δ-不變子空間從而也是 T -不變子空間不變子空間為

Cm+1 ={(x0, . . . , x_n, . . .)| xn= g(n), g 為次數不超過 m 的複係數多項式}, 即以次數不超過 m 的多項式為通項公式的數列。 由於差分運算元 Δ 在 Cm+1 上的極小多項 式為 t^m+1, 所以右平移運算元 T = Δ + 1 在 Cm+1 上的極小多項式為 (t − 1)^m+1。 因此對 方程 (24) 的討論分為以下兩種情況 :

(i) 若 P (1) = 0, 這意味著多項式 P (t) 與 (t − 1)^m+1 互質, 從而用求一術可以解出一個多 項式 U(t), 使得

P (t)U (t)≡ 1 (mod (t − 1)^m+1).

於是 (24) 的解為 x = U(T )f。 注意到, 對任意的正整數 k 有, T^kf (n) = f (n + k), 由 此不難看出, 解 x(n) = U(T )f(n) 是一個 m 次多項式。

(ii) 不然有 P (t) = (t−1)^kQ(t),其中 Q(1) = 0。 從而 P (T ) = (T −1)^kQ(T ) = Δ^kQ(T ), 於是我們可以分兩步求解 (22)。 第一步, 用前述方法求出方程 Q(T )z = f 的一個特解 z = z(n), 這是一個 m 次多項式; 第二步, 對第一步得到的 z, 求出方程 Δ^kx = z 的一 個特解 x, 我們可以保證 x(n) = [n(n − 1) · · · (n − k + 1)]w(n), 其中 w(n) 是一個 m 次多項式。 事實上, 根據反復應用楊輝粧等式的差分形式 (參見 [7, p. 71] 第 (16) 式), 我們有 : 若 z(n) 的朱世傑招差公式(參見 [7, p.67] 定理3) 為:

z(n) = a_m

n m



+ a_m−1

 n m− 1



+· · · + a0, 則

x(n) = a_m

 n k + m



+ a_m−1

 n k + m− 1



+· · · + a0

n k



= [n(n− 1) · · · (n − k + 1)] ·

m s=0

 a_ss!

(k + s)!

n− k s



是方程 Δ^kx = z 的一個特解。

x_n ^{λn //}

λT 

λⁿx_n

T

λx_n+1 ^λn//λⁿ⁺¹x_n+1

現在我們轉向一般情形 (23) 的討論。 一個簡單的分析指出 (這裡也有一個類似的交換圖, 如上), 借助變數替換 xn = λⁿy_n, 方程 (21) 等價於關於 y 的方程

P (λT )y = f, (25)

而這個方程我們在上面已經解決了!

10.3. _評注

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

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

當然, 也有作者為了繞開 「部分分式分解」, 就採取 「待定係數法」, 如高德納 (Donald E.

Knuth) 及其合作者就在書中說[4, p. 17] :

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

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

10.4. _應用舉例

例1 (取自中譯本[3, p. 189]) : 求遞推關係 xn = 5x_n−1− 6xn−2+ 7ⁿ, (n ≥ 2) 的一個特解。

解 : 令 T 為右平移運算元, 則上述方程可寫成 (注意, 右邊非齊次項多出一個 7²) P (T )x_n= 7²7ⁿ,

其中 P (T ) = T²− 5T + 6。 令 xn = 7ⁿz_n, 則上述方程等價於 Q(T )z_n= 7² = 49,

其中 Q(T ) = P (7T ) = 49T²− 35T + 6。 由於 Q(1) = 49 − 35 + 6 = 20 = 0, 所以我們 只要求出一個多項式 U(T ) 使得

Q(T )U (T ) ≡ 1 mod (T − 1).

當然我們可以用求一術, 但此處有更簡單的方法。 因為上述條件等價於 Q(1)U (1)− 1 = 0,

從而 U(1) = 1

P (1) = 1

20, 因此 U(T ) = 1

20 是同餘方程的一個特解, 從而 zn = U (T )49 = 49

20 是 Q(T )zn= 49 的一個特解, xn= 49

207ⁿ 是原方程的一個特解。

例2 (取自[3, p. 341]) : 求遞推關係 xn= x_n−1+ 2x_n−2+ (−1)ⁿ, (n≥ 2) 的通解。

解 : 分三步。

第一步, 求出一個特解。 令 T 為右平移運算元, 則上述方程等價於 P (T )x_n = (−1)ⁿ⁺²= (−1)ⁿ,

其中 P (T ) = T² − T − 2 = (T + 1)(T − 2)。 令 xn= (−1)ⁿz_n, 則上述方程等價於 zn 的 下述方程

Q(T )z_n= 1,

其中 Q(T ) = P (−T ) = (T −1)(T +2) = Δ(T +2)。 因此我們先求出 y_n使得 (T +2)y_n = 1, 只要取常數列 yn = 1

3。 進一步求解 Δzn = y_n = 1

3, 即得 Q(T )zn = 1 的一個特解 z_n= 1

3n。 從而原方程的一個特解為 xn = 1

3n(−1)ⁿ。

第二步, 求出對應齊次方程的通解。 因為特徵方程為 P (T ) = T²−T −2 = (T +1)(T − 2), 根據定理3, 我們只要分別求出 (T + 1)xn = 0 與 (T − 2)xn = 0 的通解, 並作疊加, 根 據定理 2^, 結果為

C₁(−1)ⁿ+ C₂2ⁿ, 其中 C1, C₂ 為任意複常數.

第三步, 根據線性原理, 得到非齊次方程的通解公式為 : 1

3n(−1)ⁿ+ C₁(−1)ⁿ+ C₂2ⁿ, 其中 C1, C₂ 為任意複常數.

11. 練習

1. 證明以下結果, 並說明這是求一術解整數同餘方程的理論依據 (提示 : 左乘即行變換)。

命題 : 設二階整數矩陣

A =

! a 1 b 0

"

.

(i) 若二階整數矩陣

B =

! u ∗

∗ ∗

"

使得

BA =

! 1 v

∗ ∗

"

, 則 u = v, 且 x = v 滿足同餘方程 ax ≡ 1 (mod b)。

(ii) 若二階整數矩陣

C =

! ∗ ∗ u ∗

"

使得

CA =

! ∗ ∗ 1 v

"

, 則 u = v, 且 x = v 滿足同餘方程 ax ≡ 1 (mod b)。

2. 設 D 是作用在光滑函數空間上的微分運算元, 證明, 對任意的多項式 P 與常數 λ, 有 e^−λxP (D)e^λx= P (D + λ)。 並據此給出定理 2 的一個簡單證明。

3. 設 T 是作用在數列空間上的右平移運算元, 證明, 對任意的多項式 P 與非零常數 λ, 有 λ⁻ⁿP (T )λⁿ= P (λT )。 請據此化簡定理 2^ 的證明。

4. 利用本文介紹的方法, 求出斐波拉契數列 (Fibonacci sequence) x_n = x_n−1+ x_n−2(n ≥ 2); x₀ = x₁ = 1 的通項公式。 (讀者可以參見[13, p. 36])

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

參考文獻

1. Tom M. Apostol, Linear Algebra: A First Course with Applications to Differential Equations, John Wiley & Sons, Inc. 1997. 沈灏, 沈佳辰 (譯)。 線性代數及其應用導論。 北京 : 人民郵電出版社, 2010。

2. 丁同仁, 李承治。 常微分方程教程 【第2版】。 北京 : 北京大學出版社。 2004。

3. Ronald Graham, Donald Knuth, and Oren Patashnik, Concrete Mathematics (2nd ed.).

Reading, MA: Addison-Wesley Professional. 1994. 張明堯、 張凡 (譯)。 具體數學。 北京 : 人 民郵電出版社, 2013。

4. Daniel H. Greene, Donald E. Knuth, Mathematics for the Analysis of Algorithms (3rd ed.), Birkh¨auser, 1990.

5. 龔升, 張德健。 線性代數五講 — 第五講 : 向量空間在線性運算元下的分解。 數學傳播季刊, 32(2), 34–53, 2008。

6. 華羅庚。 高等數學引論 (第四冊)。 北京 : 科學出版社。 1998。

7. 林開亮。 微積分之前奏 (或變奏) : 高階等差數列的求和。數學傳播季刊, 41(1), 61–79, 2017。

8. George P´olya, How to Solve It (with foreword by John Horton Conway and added ex- ercises), Princeton University Press, 2004. 塗泓、 馮承天 (譯)。 怎樣解題。 上海 : 上海科技 教育出版社, 2011。

9. Rogers Penrose, The Road to Reality, New York: Vintage Books, 2004.

10. 錢寶琮。 中國數學史。 北京 : 科學出版社。 1964。

11. Kenneth H. Rosen, Discrete Mathematics and its Applications. (Seventh Edition), McGraw-Hill. 徐六通, 楊娟, 吳斌 (譯), 陳瓊 (改編)。 離散數學及其應用 [本科教學版]。 北京 : 機械工業出版社, 2017。

12. Gilbert Strang, Nice function, manuscript. 見其個人主頁 http://math.mit.edu/~gs/dela/nice_functions.pdf.

13. Steven Strogatz, The Calculus of Friendship: What a Teacher and a Student Learned about Life While Corresponding about Math, Princeton University Press, 2009. 有兩個 中譯本 : 李曉東 (譯)。 心中有數的人生。 瀋陽, 萬卷出版公司, 2010; 蔡承志 (譯)。 學微積分, 也學 人生。 台北, 遠流出版公司, 2011。

14. 許光午, 李寶。 大衍求一術的演算法意義與分析。 https://arxiv.org/abs/1610.01175.2016.

15. 吳文俊。 鄧若鴻、 吳天驕 (訪問整理)。 走自己的路 — 吳文俊口述自傳。 長沙 : 湖南教育出版社, 2015。

16. 伍鴻熙。 陳省身的伯克萊歲月。 收入紀念陳省身先生文集, 丘成桐等主編, 杭州, 浙江大學出版社, 2005。

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