Numerov算法和轉(zhuǎn)移矩陣法的分析和比較

本文對薛定諤方程進(jìn)行了初步的討論，在介紹了波函數(shù)的一些基本性質(zhì)以及基本解法的前提下，選取了較為實用的兩種求解微分方程的算法——Numerov算法和轉(zhuǎn)移矩陣法，以一維定態(tài)薛定諤方程求解為例進(jìn)行求解，并對它們的原理進(jìn)行了基本的闡述。根據(jù)波函數(shù)的實際性質(zhì)，如有界性和漸變性等，對兩種方法進(jìn)行了分析和比較，繼而得到一些較為實用的結(jié)論，以及對兩種方法使用過程中的建議。

其中）（x 被稱為定態(tài)波函數(shù)，簡稱波函數(shù)，這說明，微觀粒子在一個固定的勢場中運(yùn)動時，只需要知道其能量，便可以求出其波函數(shù)。理論上，對于任何的E，都有一個波函數(shù)與之對應(yīng)，但在物理還要求波函數(shù)滿足自然條件，例如邊界值為零等，故方程只對應(yīng)一些特定的能量才有解。需要指出，定態(tài)并不意味著與時間無關(guān)，只不過與時間相關(guān)的部分都?xì)w結(jié)到相因子上去而已。

即使是一維定態(tài)薛定諤方程這樣看似簡單的方程，也只有在勢能形式極為特殊的情況才有精確解。例如：當(dāng)勢能取諧振子勢時，波函數(shù)具有厄米多項式形式。所以對于一維薛定諤方程來說，仍需要求諸于數(shù)值手段。歷史上，人們針對薛定諤方程發(fā)展了很多優(yōu)秀的計算機(jī)算法，包括矩陣對角化算法，轉(zhuǎn)移矩陣法，虛時間法，快速傅里葉變換法等等。這些算法各有優(yōu)缺點，適用于不同勢的薛定諤方程求解。

本文中，主要討論、分析和和比較求解薛定諤方程N(yùn)ummerov算法和轉(zhuǎn)移矩陣法。所有使用的程序都是在Matlab環(huán)境中開發(fā)的。

算法分析

打靶法+Numerov +Simpson算法

Numerov可以求解波函數(shù)，打靶法可以搜索本征能量值。要求非零能量本征值。本征值小于零的解是束縛解，這時它表現(xiàn)的較像兩端被固定的，在振動琴弦。波函數(shù)在經(jīng)典力學(xué)允許的范圍內(nèi)是震蕩的，在經(jīng)典力學(xué)禁止的范圍內(nèi)為指數(shù)行為。由于波函數(shù)分為兩個有不同行為的區(qū)域，直接積分會是交界區(qū)域不和諧或者數(shù)值不穩(wěn)定。知道波函數(shù)的圖像及其導(dǎo)數(shù)必須連續(xù)，也就是通常所說的平滑的曲線，可以通過從節(jié)分交界處向左積分，再向右積分。因為一定有一種歸一化系數(shù)使兩函數(shù)在交界處的值相等，所以只需要讓其導(dǎo)函數(shù)連續(xù)便可。

其中 <ψ和>ψ分別代表向左積分向右積分得到的波函數(shù)。

于是，利用拼湊的方法將不同區(qū)域內(nèi)的波函數(shù)湊成完整的一個波函數(shù)，利用打靶法獲得本征能量時，需要從最小值一點一點向上搜索，若發(fā)現(xiàn)F變號，則說明剛剛錯過了一個本征能量值，或說明波函數(shù)多了一個零點，于是減小步長，在這個范圍內(nèi)，重新搜索，獲得一個更加精確的能量值，并把它記錄下來。然后繼續(xù)向上搜索，獲得更多的能量本征值，直到粒子跳出束縛態(tài)，即E>0為止。對于一個本征波函數(shù)，采取具有較高精度的Simpson算法進(jìn)行積分，對其進(jìn)行歸一化處理。

轉(zhuǎn)移矩陣法

式（2）可寫成如下形式

由于E表示某一特定能量，V（X）也可視為某一特定勢能，故可將求解區(qū)域均勻分為N個區(qū)域，每個區(qū)域內(nèi)的波函數(shù)都近似為常數(shù)。于是式（7）可寫成

將式（12）表示成矩陣形式，并帶入邊界值，即可求出式（10）中的系數(shù)，這樣就求出了波函數(shù)。

算法優(yōu)缺點分析

上文是對于兩種算法的分析和理解。理論上兩種算法都可以對任意勢陷的薛定諤方程求以解。

在Numerov算法計算中，效率較高，可一次性將所有的束縛態(tài)的本征能量值求出。但在運(yùn)行Numerov程序時，發(fā)現(xiàn)，本征能量值共有36個，而波函數(shù)卻有50個。這里排除了多個不同的波函數(shù)對應(yīng)同一個本征能量值的簡并情況，因為求解的是一維薛定諤方程，無其他維度，故一個本征能量值，僅僅對應(yīng)一個波函數(shù)。故多出來的波函數(shù)為假態(tài)，即非物理的數(shù)學(xué)解，在此函數(shù)和其導(dǎo)函數(shù)中，有突變的情況出現(xiàn)，如圖1。

所以要利用量子力學(xué)中對波函數(shù)的要求（如波函數(shù)連續(xù)性、波函數(shù)的導(dǎo)數(shù)連續(xù)性要求等等）把這些假態(tài)剔除出去，這是此算法的缺點之一：在計算的過程中會出現(xiàn)假態(tài)。

對比圖2的左圖和右圖，不難發(fā)現(xiàn)：在波函數(shù)圖像仿真度上，轉(zhuǎn)移矩陣法不及Numerov算法。原因是若波函數(shù)在N個格點上表示，對于Numerov算法來說波函數(shù)就是一個N維數(shù)組，但對于轉(zhuǎn)移矩陣法則一個N×N維矩陣。所以轉(zhuǎn)移矩陣法會消耗很多的內(nèi)存，這就大大限制了轉(zhuǎn)移矩陣法中波函數(shù)所在的格點數(shù)目。也因此轉(zhuǎn)移矩陣法給出的圖像顯得比Numerov算法要粗糙。

從結(jié)構(gòu)上說，轉(zhuǎn)移矩陣法有很多優(yōu)點，例如：1）它能夠很好的考慮邊界條件，可以從它的推導(dǎo)過程中看到：充分考慮了函數(shù)在邊界處的取值）（aψ和）（bψ；2）通過構(gòu)造矩陣的方法來求解，可以使用很多現(xiàn)成的庫函數(shù)，如 LAPACK（Linear Algebra Package）或者擅長矩陣運(yùn)算的軟件，如MATLAB。

但是，轉(zhuǎn)移矩陣法也有一些固有缺點：因涉及到矩陣求逆的過程，在求解過程中出現(xiàn)奇異矩陣，則該辦法失效。且隨著所取的格點數(shù)目的增多，矩陣維數(shù)上升很快，會影響計算速度和運(yùn)算效率。在其程序中，可以隨意規(guī)定勢阱，但必須手動輸入所滿足的本征能量，并且其與打靶法尋找本征能量的配合很差（轉(zhuǎn)移矩陣法本身沒有判斷能量是否為本征能量的能力）。所以對于像諧振子勢這種本征能量值有特定函數(shù)表示的勢能比較合適，對于本征能量并不清楚的情況則是采用Numerov算法收效更大。

結(jié)論

根據(jù)以上分析，知道：打靶法+Numerov算法更適合于本征能量事先不知道的薛定諤方程，并且計算代價小。缺點是在求解過程中容易出現(xiàn)假態(tài)，需要額外的機(jī)制作為刪選條件。轉(zhuǎn)移矩陣法則更適用于本征能量已經(jīng)給定的情況，且改法能夠很好地處理波函數(shù)的邊界條件。涉及矩陣運(yùn)算，可以調(diào)用專門的線性代數(shù)庫函數(shù)或軟件進(jìn)行求解。

參考文獻(xiàn)

[1]大猩猩Hψ=Eψ.一維定態(tài)薛定諤方程本征值的數(shù)值解法[EB/OL].豆瓣網(wǎng)，[2011-07-20].https：//www. douban.com/note/162335796/.

[2]王憶鋒，唐利斌.利用轉(zhuǎn)移矩陣和MATLAB求解一維薛定諤方程的一種簡潔方法[J].紅外技術(shù)，2010，32（3）：177-180.

（作者簡介：張存良，東平明湖中學(xué)。）endprint