有限差分法

有限差分法（finite-difference methods，FDM），是在数学中的一种微分方程数值方法，是通过有限差分來近似导數，从而寻求微分方程的近似解。

由泰勒展開式的推導
首先假設要近似函數的各級導數都有良好的性質，依照泰勒定理，可以形成以下的泰勒展開式：

:f(x_0 + h) = f(x_0) + \frac{f'(x_0)}{1!}h + \frac{f^{(2)}(x_0)}{2!}h^2 + \cdots + \frac{f^{(n)}(x_0)}{n!}h^n + R_n(x),

其中n!表示是n的階乘，Rn(x)為餘數，表示泰勒多項式和原函數之間的差。可以推導函數f一階導數的近似值：

:f(x_0 + h) = f(x_0) + f'(x_0)h + R_1(x),

設定x0=a，可得：

:f(a+h) = f(a) + f'(a)h + R_1(x),

除以h可得：

:{f(a+h)\over h} = {f(a)\over h} + f'(a)+{R_1(x)\over h}

求解f'(a)：

:f'(a) = {f(a+h)-f(a)\over h} - {R_1(x)\over h}

假設R_1(x)相當小，因此可以將"f"的一階導數近似為：

:f'(a)\approx {f(a+h)-f(a)\over h}.

準確度及誤差
近似解的誤差定義為近似解及解析解之間的差值。有限差分法的兩個誤差來源分別是捨入誤差及截尾誤差（或稱為離散化誤差），前者是因為電腦計算小數時四捨五入造成的誤差，後者則是用有限阶级数表示导数引起的误差。

在運用有限差分法求解一問題（或是說找到問題的近似解）時，第一步需要將問題的定義域離散化。一般會將問題的定義域用均勻的網格分割（可參考右圖）。因此有限差分法會制造一組導數的離散數值近似值。

一般會關注近似解的，會用大O符號表示，局部截尾誤差是指應用有限差分法一次後產生的誤差，因此為f'(x_i) - f'_i，此時f'(x_i)是實際值，而f'_i為近似值。泰勒多項式的餘數項有助於分析局部截尾誤差。利用f(x_0 + h)泰勒多項式的餘數項，也就是

:
R_n(x_0 + h) = \frac{f^{(n+1)}(\xi)}{(n+1)!} (h)^{n+1}
, 其中 x_0 ,

可以找到局部截尾誤差的主控項，例如用前項差分法計算一階導數，已知f(x_i)=f(x_0+i h),

: f(x_0 + i h) = f(x_0) + f'(x_0)i h + \frac{f*(\xi)}{2!} (i h)^{2},

利用一些代數的處理，可得

: \frac{f(x_0 + i h) - f(x_0)}{i h} = f'(x_0) + \frac{f*(\xi)}{2!} i h,

注意到左邊的量是有限差分法的近似，右邊的量是待求解的量再加上一個餘數，因此餘數就是局部截尾誤差。上述範例可以用下式表示：

: \frac{f(x_0 + i h) - f(x_0)}{i h} = f'(x_0) + O(h).

在此例中，局部截尾誤差和時間格點的大小成正比。

範例：常微分方程
例如考慮以下的常微分方程
: u'(x) = 3u(x) + 2. \,
利用數值方法中歐拉法求解，利用以下的有限差分式
:\frac{u(x+h) - u(x)}{h} \approx u'(x)
來近似導數，並配合一些代數處理（等號兩側同乘以h，再加上u(x)），可得
: u(x+h) = u(x) + h(3u(x)+2). \,
最後的方程式即為有限差分方程，求解此方程則可得到原方程的近似解。

範例：熱傳導方程
考慮正規化的一維熱傳導方程式，為齊次的狄利克雷邊界條件

: U_t=U_{xx} \,
: U(0,t)=U(1,t)=0 \, （邊界條件）
: U(x,0) =U_0(x) \, （初始條件）

對此問題求數值解的一種方式是用差分去近似所有的導數，可以將空間分割為 x_0, ..., x_J ，將時間也分割為 t_0, ...., t_N 。假設在時間及空間都是均勻的網格切割，空間中兩個連續位置的間隔為h，兩個連續時間之間的間隔為k。點

: u(x_j,t_n) = u_{j}^n

表示 u(x_j, t_n)的數值近似解。

顯式方法
]]
利用在時間 t_n 的前向差分，以及在位置 x_j 的二階中央差分（FTCS格式），可以得到以下的迭代方程：

: \frac{u_{j}^{n+1} - u_{j}^{n}}{k} = \frac{u_{j+1}^n - 2u_{j}^n + u_{j-1}^n}{h^2}. \,

這是用求解一維導熱傳導方程的顯式方法。

可以用以下的式子求解 u_j^{n+1}

: u_{j}^{n+1} = (1-2r)u_{j}^{n} + ru_{j-1}^{n} + ru_{j+1}^{n}

其中 r=k/h^2.

因此配合此迭代關係式，已知在時間n的數值，可以求得在時間n+1的數值。 u_0^n 及 u_J^n 的數值可以用邊界條件代入，在此例中為0。

此顯式方法在 r\le 1/2 時，為數值稳定且收斂。其數值誤差和時間間隔成正比，和位置間隔的平方成正比：
: \Delta u = O(k)+O(h^2) \,

隱式方法
若使用時間 t_{n+1} 的後向差分，及位置 x_j 的二階中央差分（BTCS 格式），可以得到以下的迭代方程：

: \frac{u_{j}^{n+1} - u_{j}^{n}}{k} =\frac{u_{j+1}^{n+1} - 2u_{j}^{n+1} + u_{j-1}^{n+1}}{h^2}. \,

這是用求解一維導熱傳導方程的隱式方法。

在求解線性聯立方程後可以得到 u_j^{n+1} ：

: (1+2r)u_j^{n+1} - ru_{j-1}^{n+1} - ru_{j+1}^{n+1}= u_{j}^{n}

此方法不論r的大小，都數值稳定且收斂，但在計算量會較顯式方法要大，因為每前進一個時間間隔，就需要求解一個聯立的數值方程組。其數值誤差和時間間隔成正比，和位置間隔的平方成正比：

: \Delta u = O(k)+O(h^2) \,

克兰克－尼科尔森方法
若使用時間 t_{n+1/2} 的中間差分，及位置 x_j 的二階中央差分（CTCS 格式），可以得到以下的迭代方程：

: \frac{u_j^{n+1} - u_j^{n}}{k} = \frac{1}{2} \left(\frac{u_{j+1}^{n+1} - 2u_j^{n+1} + u_{j-1}^{n+1}}{h^2}+\frac{u_{j+1}^{n} - 2u_j^{n} + u_{j-1}^{n}}{h^2}\right).\,

此公式為克兰克－尼科尔森方法（Crank–Nicolson method）。

在求解線性聯立方程後可以得到 u_j^{n+1} ：

: (2+2r)u_j^{n+1} - ru_{j-1}^{n+1} - ru_{j+1}^{n+1}= (2-2r)u_j^n + ru_{j-1}^n + ru_{j+1}^n

此方法不論r的大小，都數值稳定且收斂，但在計算量會較顯式方法要大，因為每前進一個時間間隔，就需要求解一個聯立的數值方程組。其數值誤差和時間間隔的平方成正比，和位置間隔的平方成正比：

: \Delta u = O(k^2)+O(h^2). \,

若時間刻度較小時，克兰克－尼科尔森方法是最精確的，而顯式方法是最不精確的，而且可能會不穩定，但是是最容易計算的，其數值計算量也最少。若時間刻度較大時，隱式方法的效果最好。

相關條目
*有限元分析
*差分
*時域有限差分
*模版 (數值分析)
*有限差分係數
*五點模版
*Lax等價定理
*

參考資料

• K.W. Morton and D.F. Mayers, Numerical Solution of Partial Differential Equations, An Introduction. Cambridge University Press, 2005.
• Oliver Rübenkönig, [https://web.archive.org/web/20100111163529/http://www.imtek.de/simulation/mathematica/IMSweb/imsTOC/Lectures%20and%20Tips/Simulation%20I/FDM_introDocu.html The Finite Difference Method (FDM) - An introduction], (2006) Albert Ludwigs University of Freiburg
• Autar Kaw and E. Eric Kalu, Numerical Methods with Applications, (2008) [http://www.autarkaw.com/books/numericalmethods/index.html]

外部連結

• [https://web.archive.org/web/20140302111344/http://emlab.utep.edu/ee5390cem.htm Finite-Difference Method (see and listen to lecture 9)]

*[https://web.archive.org/web/20060818132006/http://math.fullerton.edu/mathews/n2003/finitediffpde/FiniteDifferencePDEBib/Links/FiniteDifferencePDEBib_lnk_1.html List of Internet Resources for the Finite Difference Method for PDEs]

• [http://numericalmethods.eng.usf.edu/topics/finite_difference_method.html Finite Difference Method of Solving ODEs (Boundary Value Problems) Notes, PPT, Maple, Mathcad, Matlab, Mathematica]

*[https://web.archive.org/web/20090206071754/http://ltl.iams.sinica.edu.tw/document/training_lectures/2006/SH_Chen/Finite_Difference_Methods.pdf Lecture Notes] Shih-Hung Chen, National Central University

• , [http://faculty.washington.edu/rjl/fdmbook/ Finite Difference Methods for Ordinary and Partial Differential Equations] , SIAM, 2007.

*[https://web.archive.org/web/20120226144927/http://www.adeptscience.co.uk/products/mathsim/maple/powertools/des/unit18.html#MapleAutoBookmark4 Finite Difference Method]
*[http://math.fullerton.edu/mathews/n2003/FiniteDifferenceMod.html Finite Difference Method for Boundary Value Problems]
*[https://web.archive.org/web/20120207024810/http://www.composite-agency.com/product.htm Finite Difference Methodology in Materials Science]
*[http://jwhaverkort.net23.net/documents/NumMethPDEs.pdf Numerical Methods for time-dependent Partial Differential Equations]