龍格-庫塔法

數值分析中,龍格-庫塔法(英文:Runge-Kutta methods)是用於非線性常微分方程的解的重要的一類隱式或顯式迭代法。這些技術由數學家卡爾·龍格馬丁·威爾海姆·庫塔於1900年左右發明。

四階龍格-庫塔法

在各種龍格-庫塔法當中有一個方法十分常用,以至於經常被稱為「RK4」或者就是「龍格-庫塔法」。該方法主要是在已知方程導數和初始值時,利用計算機的仿真應用,省去求解微分方程的複雜過程。

初值問題表述如下。

 

則,對於該問題的RK4由如下方程給出:

 

其中

 
 
 
 

這樣,下一個值(yn+1)由現在的值(yn)加上時間間隔(h)和一個估算的斜率的乘積所決定。該斜率是以下斜率的加權平均:

  • k1是時間段開始時的斜率;
  • k2是時間段中點的斜率,通過歐拉法採用斜率k1來決定y在點tn + h/2的值;
  • k3也是中點的斜率,但是這次採用斜率k2決定y值;
  • k4是時間段終點的斜率,其y值用k3決定。

當四個斜率取平均時,中點的斜率有更大的權值:

 

RK4法是四階方法,也就是說每步的誤差是h5,而總積累誤差為h4階。

注意上述公式對於純量或者向量函數(y可以是向量)都適用。

顯式龍格-庫塔法

顯式龍格-庫塔法是上述RK4法的一個推廣。它由下式給出

 

其中

 
 
 
 
 

(注意:上述方程在不同著述中有不同但等價的定義)。

要給定一個特定的方法,必須提供整數s(級數),以及係數 aij(對於1 ≤ j < is), bi(對於i = 1, 2, ..., s)和ci(對於i = 2, 3, ..., s)。這些數據通常排列在一個助記工具中,稱為Butcher tableau(得名自約翰·C·布徹英語John C. Butcher):

0
   
     
     
         
         

龍格-庫塔法是自洽的,如果

 

如果要求方法的精度為p階,即截斷誤差為O(hp+1)的,則還有相應的條件。這些可以從截斷誤差本身的定義中導出。例如,一個2級2階方法要求b1 + b2 = 1, b2c2 = 1/2, 以及b2a21 = 1/2。

例子

RK4法處於這個框架之內。其表為:

0
1/2 1/2
1/2 0 1/2
1 0 0 1
1/6 1/3 1/3 1/6

然而,最簡單的龍格-庫塔法是(更早發現的)歐拉方法,其公式為 。這是唯一自洽的一級顯式龍格-庫塔法。相應的表為:

0
1

隱式龍格-庫塔法

以上提及的顯式龍格-庫塔法一般來講不適用於求解剛性方程。這是因為顯式龍格-庫塔法的穩定區域被局限在一個特定的區域裏。顯式龍格-庫塔法的這種缺陷使得人們開始研究隱式龍格-庫塔法,一般而言,隱式龍格-庫塔法具有以下形式:

 

其中

 

在顯式龍格-庫塔法的框架里,定義參數 的矩陣是一個下三角矩陣,而隱式龍格-庫塔法並沒有這個性質,這是兩個方法最直觀的區別:

 

需要注意的是,與顯式龍格-庫塔法不同,隱式龍格-庫塔法在每一步的計算里需要求解一個線性方程組,這相應的增加了計算的成本。

參考

  • George E. Forsythe, Michael A. Malcolm, and Cleve B. Moler. Computer Methods for Mathematical Computations. Englewood Cliffs, NJ: Prentice-Hall, 1977. (See Chapter 6.)
  • Ernst Hairer, Syvert Paul Nørsett, and Gerhard Wanner. Solving ordinary differential equations I: Nonstiff problems, second edition. Berlin: Springer Verlag, 1993. ISBN 3-540-56670-8.
  • William H. Press, Brian P. Flannery, Saul A. Teukolsky, William T. Vetterling. Numerical Recipes in C. Cambridge, UK: Cambridge University Press, 1988. (See Sections 16.1 and 16.2.)