Phương pháp Runge-Kutta

Trong giải tích số, các phương pháp Runge-Kutta là một họ của các phương pháp lặp ẩn (implicit) và hiện (explicit), trong đó bao gồm thường trình nổi tiếng được gọi là các phương pháp Euler, được sử dụng trong việc rời rạc hóa thời gian để tìm lời giải gần đúng cho các phương trình vi phân thường. Những phương pháp này được phát triển vào khoảng năm 1900 bởi hai nhà toán học người Đức C. Runge và M. W. Kutta.

Xem bài viết về các phương pháp số áp dụng cho việc tìm lời giải của các phương trình vi phân thường để hiểu hơn về nền tảng cơ sở và các phương pháp số khác. Xem thêm Danh sách các phương pháp Runge-Kutta.

Phương pháp Runge-Kutta

Thành viên được biết đến rộng rãi nhất của họ Runge-Kutta là "RK4", "phương pháp Runge-Kutta cổ điển" hoặc đơn giản là "phương pháp Runge-Kutta".

Cho một bài toán giá trị ban đầu được chỉ rõ như sau:

${\displaystyle {\dot {y}}=f(t,y),\quad y(t_{0})=y_{0}.}$

Ở đây y là một hàm chưa biết (vô hướng hoặc vector) của thời gian t, mà chúng ta muốn tìm lời giải gần đúng; chúng ta biết rằng ${\displaystyle {\dot {y}}}$, tốc độ thay đổi của y, là một hàm của t và y. Tại thời điểm ban đầu ${\displaystyle t_{0}}$ giá trị y tương ứng là ${\displaystyle y_{0}}$. Hàm f và các dữ liệu ${\displaystyle t_{0}}$, ${\displaystyle y_{0}}$ được cho trước.

Bây giờ chọn một bước kích thước h> 0 và định nghĩa

${\displaystyle {\begin{aligned}y_{n+1}&=y_{n}+{\tfrac {h}{6}}\left(k_{1}+2k_{2}+2k_{3}+k_{4}\right),\\t_{n+1}&=t_{n}+h\\\end{aligned}}}$

với n = 0, 1, 2, 3,..., sử dụng^[1]

${\displaystyle {\begin{aligned}k_{1}&=f(t_{n},y_{n}),\\k_{3}&=f(t_{n}+{\frac {h}{2}},y_{n}+{\frac {h}{2}}k_{2}),\\k_{4}&=f(t_{n}+h,y_{n}+hk_{3}).\end{aligned}}}$

(Lưu ý: các phương trình trên có những định nghĩa khác nhau nhưng tương đương trong các văn bản khác nhau).^[2]

Ở đây ${\displaystyle y_{n+1}}$ là giá trị gần đúng (xấp xỉ) RK4 của ${\displaystyle y(t_{n+1})}$, và giá trị tiếp theo (${\displaystyle y_{n+1}}$) được xác định bằng giá trị hiện tại (${\displaystyle y_{n}}$) cộng với trung bình trọng lượng của bốn số gia, trong đó mỗi số gia là sản phẩm của kích cỡ của khoảng thời gian, h, và một độ dốc ước tính được chỉ rõ bởi hàm f ở phía bên tay phải của phương trình vi phân.

• ${\displaystyle k_{1}}$ là số gia dựa trên độ dốc tại điểm bắt đầu của khoảng thời gian, sử dụng ${\displaystyle y}$ (phương pháp Euler);
• ${\displaystyle k_{2}}$ là số gia dựa trên độ dốc tại điểm giữa của khoảng thời gian, sử dụng ${\displaystyle y+{\tfrac {h}{2}}k_{1}}$;
• ${\displaystyle k_{3}}$ cũng là số gia dựa trên độ dốc ở điểm giữa, nhưng sử dụng ${\displaystyle y+{\tfrac {h}{2}}k_{2}}$;
• ${\displaystyle k_{4}}$ là số gia dựa trên độ dốc tại điểm cuối của khoảng thời gian, sử dụng ${\displaystyle y+hk_{3}}$.

Trong việc trung bình của bốn số gia, trọng lượng lớn hơn được trao cho các số gia tại điểm giữa. Nếu ${\displaystyle f}$ độc lập với ${\displaystyle y}$, để các phương trình vi phân tương đương với một tích phân đơn giản, thì RK4 là quy tắc Simpson.^[3]

Phương pháp RK4 là một phương pháp bậc 4, có nghĩa rằng sai số cắt cụt cục bộ bậc ${\displaystyle O(h^{5})}$, trong khi tổng lỗi tích lũy là bậc ${\displaystyle O(h^{4})}$.

Các phương pháp Runge-Kutta hiện

Họ các phương pháp Runge-Kutta hiện (explicit) là một sự tổng quát hóa của phương pháp RK4 đã đề cập bên trên. Nó được cho bởi công thức:

${\displaystyle y_{n+1}=y_{n}+h\sum _{i=1}^{s}b_{i}k_{i},}$

trong đó ^[4]

${\displaystyle {\begin{aligned}k_{1}&=f(t_{n},y_{n}),\\k_{2}&=f(t_{n}+c_{2}h,y_{n}+h(a_{21}k_{1})),\\k_{3}&=f(t_{n}+c_{3}h,y_{n}+h(a_{31}k_{1}+a_{32}k_{2})),\\&\ \ \vdots \\k_{s}&=f(t_{n}+c_{s}h,y_{n}+h(a_{s1}k_{1}+a_{s2}k_{2}+\cdots +a_{s,s-1}k_{s-1})).\end{aligned}}}$

(Lưu ý: các phương trình bên trên có các định nghĩa khác nhau nhưng tương đương trong các văn bản khác nhau).^[2]

Để chỉ rõ một phương pháp cụ thể, cần cung cấp số nguyên s (số giai đoạn), và các hệ số a_ij (đối với 1 ≤ j < i ≤ s), b_i (đối với i = 1, 2,..., s) và c_i (đối với i = 2, 3,..., s). Ma trận [a_ij] được gọi là ma trận Runge–Kutta, trong khi b_i và c_i được biết đến như là các trọng lượng và các điểm nút.^[5] Những dữ liệu này thường được sắp xếp trong một thiết bị ghi nhớ, được biết đến như là một hoạt cảnh (tableau) Butcher (đặt theo tên của John C. Butcher):

${\displaystyle 0}$
${\displaystyle c_{2}}$	${\displaystyle a_{21}}$
${\displaystyle c_{3}}$	${\displaystyle a_{31}}$	${\displaystyle a_{32}}$
${\displaystyle \vdots }$	${\displaystyle \vdots }$		${\displaystyle \ddots }$
${\displaystyle c_{s}}$	${\displaystyle a_{s1}}$	${\displaystyle a_{s2}}$	${\displaystyle \cdots }$	${\displaystyle a_{s,s-1}}$
	${\displaystyle b_{1}}$	${\displaystyle b_{2}}$	${\displaystyle \cdots }$	${\displaystyle b_{s-1}}$	${\displaystyle b_{s}}$

Phương pháp Runge-Kutta là nhất quán nếu

${\displaystyle \sum _{j=1}^{i-1}a_{ij}=c_{i}{\text{ for }}i=2,\ldots ,s.}$

Ngoài ra còn có các yêu cầu kèm theo nếu yêu cầu phương pháp phải có bậc p nào đó, có nghĩa rằng sai số cắt cụt cục bộ là O(h^p+1). Những điều này có thể được rút ra từ định nghĩa của chính sai số cắt cụt. Cho ví dụ, một phương pháp hai giai đoạn có bậc 2 nếu b₁ + b₂ = 1, b₂c₂ = 1/2, và a₂₁ = c₂.^[6]

Nhìn chung, nếu một phương pháp Runge-Kutta hiện ${\displaystyle s}$-giai đoạn có ${\displaystyle s\geq p}$, và nếu ${\displaystyle p\geq 5}$, thì ${\displaystyle s>p}$.^[7] ${\displaystyle s}$ nhỏ nhất cần thiết cho một phương pháp Runge-Kutta hiện ${\displaystyle s}$-giai đoạn để có bậc ${\displaystyle p}$ là một vấn đề mở. Một số giá trị đã biết:^[8]

${\displaystyle {\begin{array}{c|cccccccc}p&1&2&3&4&5&6&7&8\\\hline \min s&1&2&3&4&6&7&9&11\end{array}}}$

Ví dụ

Phương pháp RK4 trong khuôn khổ này. Hoạt cảnh của nó là ^[9]

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

Một sự thay đổi đôi chút của phương pháp Runge-Kutta cũng là do Kutta tạo ra năm 1901 và được gọi là quy tắc-3/8.^[10] Ưu điểm chính phương pháp này là hầu như tất cả các hệ số lỗi là nhỏ hơn so với phương pháp thông thường, nhưng nó đòi hỏi nhiều FLOPS hơn đôi chút (các thao tác điểm-nổi) cho mỗi bước thời gian. Hoạt cảnh Butcher của nó là:

0
1/3	1/3
2/3	−1/3	1
1	1	−1	1
	1/8	3/8	3/8	1/8

Tuy nhiên, phương pháp Runge-Kutta đơn giản nhất là phương pháp Euler (tiếp tới = forward), được cho bởi công thức ${\displaystyle y_{n+1}=y_{n}+hf(t_{n},y_{n})}$. Đây là phương pháp Runge-Kutte hiện nhất quán duy nhất với một giai đoạn. Hoạt cảnh tương ứng là

0
	1

Các phương pháp bậc-hai với hai giai đoạn

Một ví dụ của một phương pháp bậc hai với hai giai đoạn là phương pháp điểm giữa:

${\displaystyle y_{n+1}=y_{n}+hf\left(t_{n}+{\frac {1}{2}}h,y_{n}+{\frac {1}{2}}hf(t_{n},\ y_{n})\right).}$

Hoạt cảnh tương ứng là

0
1/2	1/2
	0	1

Phương pháp điểm giữa không phải là phương pháp Runge-Kutta bậc-hai duy nhất với hai giai đoạn; có một họ những phương pháp như vậy, được tham số hóa bởi α và cho bởi công thức^[11]

${\displaystyle y_{n+1}=y_{n}+h{\bigl (}(1-{\tfrac {1}{2\alpha }})f(t_{n},y_{n})+{\tfrac {1}{2\alpha }}f(t_{n}+\alpha h,y_{n}+\alpha hf(t_{n},y_{n})){\bigr )}.}$

Hoạt cảnh Butcher của nó là

0
${\displaystyle \alpha }$	${\displaystyle \alpha }$
	${\displaystyle (1-{\tfrac {1}{2\alpha }})}$	${\displaystyle {\tfrac {1}{2\alpha }}}$

Trong họ này, ${\displaystyle \alpha ={\tfrac {1}{2}}}$ đối với phương pháp điểm giữa, và ${\displaystyle \alpha =1}$ đối với phương pháp Heun.^[3]

Ứng dụng

Như một ví dụ, xem xét phương pháp Runge-Kutta bậc-hai hai-giai đoạn với α = 2/3, cũng được biết như là phương pháp Ralston. Nó được cho bởi hoạt cảnh

	0
	2/3	2/3
		1/4	3/4

với các phương trình tương ứng

${\displaystyle {\begin{aligned}k_{1}&=f(t_{n},\ y_{n}),\\k_{2}&=f(t_{n}+{\tfrac {2}{3}}h,\ y_{n}+{\tfrac {2}{3}}hk_{1}),\\y_{n+1}&=y_{n}+h\left({\tfrac {1}{4}}k_{1}+{\tfrac {3}{4}}k_{2}\right).\end{aligned}}}$

Phương pháp này được sử dụng để giải các bài toán giá trị-ban đầu

${\displaystyle y'=\tan(y)+1,\quad y_{0}=1,\ t\in [1,1.1]}$

với kích cỡ bước h = 0.025, vì vậy phương pháp này cần phải thực hiện bốn bước

Phương pháp được tiến hành như sau:

${\displaystyle t_{0}=1\colon }$
	${\displaystyle y_{0}=1}$
${\displaystyle t_{1}=1.025\colon }$
	${\displaystyle y_{0}=1}$	${\displaystyle k_{1}=2.557407725}$	${\displaystyle k_{2}=f(t_{0}+{\tfrac {2}{3}}h,\ y_{0}+{\tfrac {2}{3}}hk_{1})=2.7138981184}$
	${\displaystyle y_{1}=y_{0}+h({\tfrac {1}{4}}k_{1}+{\tfrac {3}{4}}k_{2})={\underline {1.066869388}}}$
${\displaystyle t_{2}=1.05\colon }$
	${\displaystyle y_{1}=1.066869388}$	${\displaystyle k_{1}=2.813524695}$	${\displaystyle k_{2}=f(t_{1}+{\tfrac {2}{3}}h,\ y_{1}+{\tfrac {2}{3}}hk_{1})}$
	${\displaystyle y_{2}=y_{1}+h({\tfrac {1}{4}}k_{1}+{\tfrac {3}{4}}k_{2})={\underline {1.141332181}}}$
${\displaystyle t_{3}=1.075\colon }$
	${\displaystyle y_{2}=1.141332181}$	${\displaystyle k_{1}=3.183536647}$	${\displaystyle k_{2}=f(t_{2}+{\tfrac {2}{3}}h,\ y_{2}+{\tfrac {2}{3}}hk_{1})}$
	${\displaystyle y_{3}=y_{2}+h({\tfrac {1}{4}}k_{1}+{\tfrac {3}{4}}k_{2})={\underline {1.227417567}}}$
${\displaystyle t_{4}=1.1\colon }$
	${\displaystyle y_{3}=1.227417567}$	${\displaystyle k_{1}=3.796866512}$	${\displaystyle k_{2}=f(t_{3}+{\tfrac {2}{3}}h,\ y_{3}+{\tfrac {2}{3}}hk_{1})}$
	${\displaystyle y_{4}=y_{3}+h({\tfrac {1}{4}}k_{1}+{\tfrac {3}{4}}k_{2})={\underline {1.335079087}}.}$

Lời giải số tương ứng với giá trị được gạch chân.

Các phương pháp Runge-Kutta thích ứng

Các phương pháp thích ứng được thiết kế để tạo ra ước tính của sai số cắt cụt cục bộ của một bước Runge-Kutta đơn nhất. Điều này được thực hiện bằng cách dùng hai phương pháp trong hoạt cảnh, một với bậc ${\displaystyle p}$ và một với bậc ${\displaystyle p-1}$.

Bước bậc-thấp hơn được cho bởi

${\displaystyle y_{n+1}^{*}=y_{n}+h\sum _{i=1}^{s}b_{i}^{*}k_{i},}$

trong đó ${\displaystyle k_{i}}$ là giống hệt như đối với phương pháp bậc-cao hơn. Và lỗi là

${\displaystyle e_{n+1}=y_{n+1}-y_{n+1}^{*}=h\sum _{i=1}^{s}(b_{i}-b_{i}^{*})k_{i},}$

với bậc ${\displaystyle O(h^{p})}$. Hoạt cảnh Butcher cho loại phương pháp này được mở rộng để cho ra các giá trị của ${\displaystyle b_{i}^{*}}$:

	0
	${\displaystyle c_{2}}$	${\displaystyle a_{21}}$
	${\displaystyle c_{3}}$	${\displaystyle a_{31}}$	${\displaystyle a_{32}}$
	${\displaystyle \vdots }$	${\displaystyle \vdots }$		${\displaystyle \ddots }$
	${\displaystyle c_{s}}$	${\displaystyle a_{s1}}$	${\displaystyle a_{s2}}$	${\displaystyle \cdots }$	${\displaystyle a_{s,s-1}}$
		${\displaystyle b_{1}}$	${\displaystyle b_{2}}$	${\displaystyle \cdots }$	${\displaystyle b_{s-1}}$	${\displaystyle b_{s}}$
		${\displaystyle b_{1}^{*}}$	${\displaystyle b_{2}^{*}}$	${\displaystyle \cdots }$	${\displaystyle b_{s-1}^{*}}$	${\displaystyle b_{s}^{*}}$

Phương pháp Runge-Kutta-Fehlberg có hai phương pháp bậc 5 và 4. Hoạt cảnh Butcher mở rộng của nó là:

	0
	1/4	1/4
	3/8	3/32	9/32
	12/13	1932/2197	−7200/2197	7296/2197
	1	439/216	−8	3680/513	-845/4104
	1/2	−8/27	2	−3544/2565	1859/4104	−11/40
		16/135	0	6656/12825	28561/56430	−9/50	2/55
		25/216	0	1408/2565	2197/4104	−1/5	0

Tuy nhiên, phương pháp Runge-Kutta thích ứng đơn giản nhất liên quan đến việc kết hợp phương pháp Heun (có bậc 2), với phương pháp Euler (có bậc 1). Hoạt cảnh Butcher mở rộng là:

	0
	1	1
		1/2	1/2
		1	0

Ước tính lỗi được sử dụng để kiểm soát kích thước bước.

Các phương pháp Runge-Kutta thích ứng khác là phương pháp Bogacki-Shampine (bậc 3 và bậc 2), phương pháp Cash-Karp và phương pháp Dormand-Prince (cả hai đều có bậc 5 và 4).

Các phương pháp Runge-Kutta phi hợp lưu

Một phương pháp Runge-Kutta được gọi là phi hợp lưu ^[12] nếu tất cả ${\displaystyle c_{i},\,i=1,2,\ldots ,s}$ đều khác nhau.

Các phương pháp Runge-Kutta ẩn

Tất cả các phương pháp Runge-Kutta đã đề cập cho đến bây giờ đều là các phương pháp hiện (explicit). Các phương pháp Runge-Kutta hiện nói chung không phù hợp cho việc tìm lời giải của các phương trình cứng vì vùng ổn định tuyệt đối của chúng rất nhỏ; cụ thể là, bị bao bọc.^[13] Vấn đề này đặc biệt quan trọng trong việc tìm lời giải cho các phương trình vi phân từng phần.

Sự bất ổn định của phương pháp Runge-Kutta hiện thúc đẩy sự phát triển của các phương pháp ẩn (implicit). Một phương pháp Runge-Kutta ẩn có dạng

${\displaystyle y_{n+1}=y_{n}+h\sum _{i=1}^{s}b_{i}k_{i},}$

trong đó

${\displaystyle k_{i}=f\left(t_{n}+c_{i}h,\ y_{n}+h\sum _{j=1}^{s}a_{ij}k_{j}\right),\quad i=1,\ldots ,s.}$ ^[14]

Sự khác nhau so với phương pháp hiện là trong một phương pháp hiện, phép tính tổng đối với chỉ số j chỉ lên đến i − 1. Điều này cũng thấy rõ trong các hoạt cảnh Butcher: ma trận hệ số ${\displaystyle a_{ij}}$ của một phương pháp hiện là ma trận tam giác dưới Trong một phương pháp ẩn, phép tổng đối với j tính đến s và ma trận hệ số không phải là ma trận tam giác, cho kết quả là bảng Butcher có dạng^[9]

${\displaystyle {\begin{array}{c|cccc}c_{1}&a_{11}&a_{12}&\dots &a_{1s}\\c_{2}&a_{21}&a_{22}&\dots &a_{2s}\\\vdots &\vdots &\vdots &\ddots &\vdots \\c_{s}&a_{s1}&a_{s2}&\dots &a_{ss}\\\hline &b_{1}&b_{2}&\dots &b_{s}\\&b_{1}^{*}&b_{2}^{*}&\dots &b_{s}^{*}\\\end{array}}={\begin{array}{c|c}\mathbf {c} &A\\\hline &\mathbf {b^{T}} \\\end{array}}}$

Kết quả của sự khác biệt này là tại mỗi bước, một hệ các phương trình đại số phải được giải. Điều này làm tăng chi phí tính toán một cách đáng kể. Nếu một phương pháp với s giai đoạn được sử dụng để giải một phương trình vi phân với m thành phần, thì hệ các phương trình đại số có ms thành phần. Điều này có thể tương phản với các phương pháp đa bước tuyến tính ẩn (một họ lớn khác của các phương pháp sử dụng cho việc tìm lời giải của các phương trình vi phân thường ODEs): một phương pháp đa bước tuyến tính s-bước ẩn cần giải một hệ các phương trình đại số với chỉ m thành phần, vì vậy kích thước của hệ không tăng khi số lượng các bước tăng lên.^[15]

Ví dụ

Ví dụ đơn giản của một phương pháp Runge-Kutta implicit là phương pháp Euler lùi (backward):

${\displaystyle y_{n+1}=y_{n}+hf(t_{n}+h,\ y_{n+1}).\,}$

Hoạt cảnh Butcher chỉ đơn giản là:

${\displaystyle {\begin{array}{c|c}1&1\\\hline &1\\\end{array}}}$

Hoạt cảnh Butcher này tương ứng với các công thức

${\displaystyle k_{1}=f(t_{n}+h,\ y_{n}+hk_{1})\quad {\text{and}}\quad y_{n+1}=y_{n}+hk_{1},}$

có thể được sắp xếp lại để có được công thức cho phương pháp Euler lùi được liệt kê ở trên.

Một ví dụ nữa của phương pháp Runge-Kutta ẩn là quy tắc hình thang. Hoạt cảnh Butcher của nó là:

${\displaystyle {\begin{array}{c|cc}0&0&0\\1&{\frac {1}{2}}&{\frac {1}{2}}\\\hline &{\frac {1}{2}}&{\frac {1}{2}}\\&1&0\\\end{array}}}$

Quy tắc hình thang là một phương pháp sắp đặt theo thứ tự (collocation) (như thảo luận trong bài viết về phương pháp này). Tất cả các phương pháp sắp đặt theo thứ tự đều là các phương pháp Runge-Kutta ẩn, nhưng không phải tất cả các phương pháp Runge-Kutta ẩn là các phương pháp sắp đặt theo thứ tự.^[16]

Các phương pháp Gauss-Legendre tạo thành một họ của các phương pháp sắp đặt theo thứ tự dựa trên phép cầu phương Gauss quadrature. Một phương pháp Gauss-Legendre với s giai đoạn s có bậc 2s (như vậy, các phương pháp với bậc cao tùy ý có thể được xây dựng).^[17] Phương pháp với hai giai đoạn (và do đó có bậc bốn) có hoạt cảnh Butcher:

${\displaystyle {\begin{array}{c|cc}{\frac {1}{2}}-{\frac {1}{6}}{\sqrt {3}}&{\frac {1}{4}}&{\frac {1}{4}}-{\frac {1}{6}}{\sqrt {3}}\\{\frac {1}{2}}+{\frac {1}{6}}{\sqrt {3}}&{\frac {1}{4}}+{\frac {1}{6}}{\sqrt {3}}&{\frac {1}{4}}\\\hline &{\frac {1}{2}}&{\frac {1}{2}}\\&{\frac {1}{2}}+{\frac {1}{2}}{\sqrt {3}}&{\frac {1}{2}}-{\frac {1}{2}}{\sqrt {3}}\end{array}}}$ ^[15]

Độ ổn định

Ưu điểm của các phương pháp Runge-Kutta ẩn so với các phương pháp hiện là độ ổn định của chúng lớn hơn, đặc biệt khi áp dụng cho các phương trình cứng. Xét phương trình thử nghiệm tuyến tính y' = λy. Một phương pháp Runge-Kutta được áp dụng cho phương trình này, để có phép lặp ${\displaystyle y_{n+1}=r(h\lambda )\,y_{n}}$, với r được cho bởi

${\displaystyle r(z)=1+zb^{T}(I-zA)^{-1}e={\frac {\det(I-zA+zeb^{T})}{\det(I-zA)}},}$ ^[18]

trong đó e là viết tắt của vector với các thành phần đều bằng một. Hàm r được gọi là hàm tính ổn định.^[19] Từ công thức, ta có r là thương số của hai đa thức bậc s nếu phương pháp có s giai đoạn. Các phương pháp hiện có ma trận A là ma trận tam giác thấp với các thành phần bằng không dọc đường chéo, tức là det(I − zA) = 1 và rằng hàm tính ổn định là một đa thức.^[20]

Lời giải số cho phương trình thử nghiệm tuyến tính phân rã về không nếu | r(z) | < 1 với z = hλ. Tập hợp của các z đó được gọi là miền ổn định tuyệt đối. Đặc biệt, phương pháp được gọi là ổn định-A nếu tất cả z với Re(z) < 0 đều nằm trong miền ổn định tuyệt đối. Hàm ổn định của một phương pháp Runge-Kutta hiện là một đa thức, vì vậy các phương pháp Runge-Kutta hiện có thể không bao giờ là ổn định-A.^[20]

Nếu phương pháp có bậc p, thì hàm tính ổn định thỏa mãn ${\displaystyle r(z)={\textrm {e}}^{z}+O(z^{p+1})}$ khi ${\displaystyle z\to 0}$. Vì vậy, quan tâm nghiên cứu các thương số của các đa thức với bậc cho trước gần giống so với hàm mũ nhất. Chúng được gọi là các Padé approximants. Một Padé approximant với tử số bậc m và mẫu số bậc n là ổn định-A khi và chỉ khi m ≤ n ≤ m + 2.^[21]

Phương pháp Gauss-Legendre với s giai đoạn có bậc 2s, vì vậy hàm tính ổn định của nó là Padé approximant với m = n = s. Do đó phương pháp này là ổn định-A.^[22] Điều này cho thấy rằng Runge-Kutta ổn định-A có thể có bậc cao tùy ý. Ngược lại, bậc của các phương pháp đa bước tuyến tính ổn định-A không thể vượt quá hai.^[23]

Độ ổn định-B

Khái niệm độ ổn định-A cho lời giải của các phương trình vi phân, có liên quan đến phương trình tự trị tuyến tính ${\displaystyle y'=\lambda y}$. Dahlquist đã đề xuất việc nghiên cứu độ ổn định của các lược đồ số (numerical schemes) khi áp dụng cho các hệ phi tuyến thỏa mãn một điều kiện đơn điệu. Các khái niệm tương ứng được xác định như là độ ổn định-G đối với các phương pháp đa bước (và các phương pháp một-chân liên quan) và độ ổn định-B (Butcher, 1975) cho các phương pháp Runge-Kutta. Một phương pháp Runge-Kutta được áp dụng cho hệ phi tuyến tính ${\displaystyle y'=f(y)}$, với ${\displaystyle \langle f(y)-f(z),\ y-z\rangle <0}$, được gọi là ổn định-B, nếu điều kiện này dẫn đến ${\displaystyle \|y_{n+1}-z_{n+1}\|\leq \|y_{n}-z_{n}\|}$.

Cho ${\displaystyle B}$, ${\displaystyle M}$ và ${\displaystyle Q}$ là ba ma trận ${\displaystyle s\times s}$ được xác định bởi

${\displaystyle B=\operatorname {diag} (b_{1},b_{2},\ldots ,b_{s}),\,M=BA+A^{T}B-bb^{T},\,Q=BA^{-1}+A^{-T}B-A^{-T}bb^{T}A^{-1}.}$

Một phương pháp Runge-Kutta được gọi là ổn định về mặt đại số ^[24] nếu các ma trận ${\displaystyle B}$ and ${\displaystyle M}$ đều xác định không âm. Điều kiện đủ cho độ ổn định-B ^[25] là: B and Q là xác định không âm.

Nguồn gốc của phương pháp Runge-Kutta bậc bốn

Nói chung một phương pháp Runge-Kutta bậc 4 có thể được viết như sau:

${\displaystyle y_{t+h}=y_{t}+h\cdot \sum _{i=1}^{s}a_{i}k_{i}+{\mathcal {O}}(h^{s+1}),}$

trong đó:

${\displaystyle k_{i}=f\left(y_{t}+h\cdot \sum _{j=1}^{s}\beta _{ij}k_{j},\ t_{n}+\alpha _{i}h\right)}$

là các số gia đạt được đánh giá các đạo hàm của ${\displaystyle y_{t}}$ tại bậc thứ ${\displaystyle i}$.

Chúng ta rút ra^[26] phương pháp Runge-Kutta bậc bốn bằng cách áp dụng công thức tổng quát với trường hợp ${\displaystyle s=4}$, như đã giải thích ở trên, tại điểm khởi đầu, điểm giữa và điểm cuối của bất kỳ khoảng ${\displaystyle (t,\ t+h)}$, do đó chúng ta chọn:

${\displaystyle {\begin{aligned}&\alpha _{i}&&\beta _{ij}\\\alpha _{1}&=0&\beta _{21}&={\frac {1}{2}}\\\alpha _{2}&={\frac {1}{2}}&\beta _{32}&={\frac {1}{2}}\\\alpha _{3}&={\frac {1}{2}}&\beta _{43}&=1\\\alpha _{4}&=1&&\\\end{aligned}}}$

và mặt khác ${\displaystyle \beta _{ij}=0}$. Chúng ta bắt đầu bằng cách định nghĩa các đại lượng sau đây:

${\displaystyle {\begin{aligned}y_{t+h}^{1}&=y_{t}+hf\left(y_{t},\ t\right)\\y_{t+h}^{2}&=y_{t}+hf\left(y_{t+h/2}^{1},\ t+{\frac {h}{2}}\right)\\y_{t+h}^{3}&=y_{t}+hf\left(y_{t+h/2}^{2},\ t+{\frac {h}{2}}\right)\end{aligned}}}$

trong đó ${\displaystyle y_{t+h/2}^{1}={\dfrac {y_{t}+y_{t+h}^{1}}{2}}}$ và ${\displaystyle y_{t+h/2}^{2}={\dfrac {y_{t}+y_{t+h}^{2}}{2}}}$ Nếu chúng ta định nghĩa:

${\displaystyle {\begin{aligned}k_{1}&=f(y_{t},\ t)\\k_{2}&=f\left(y_{t+h/2}^{1},\ t+{\frac {h}{2}}\right)\\k_{3}&=f\left(y_{t+h/2}^{2},\ t+{\frac {h}{2}}\right)\\k_{4}&=f\left(y_{t+h}^{3},\ t+h\right)\end{aligned}}}$

và với các mối quan hệ trước đó chúng ta có thể thấy rằng các đẳng thức sau đúng đến bậc ${\displaystyle {\mathcal {O}}(h^{2})}$:

${\displaystyle {\begin{aligned}k_{2}&=f\left(y_{t+h/2}^{1},\ t+{\frac {h}{2}}\right)=f\left(y_{t}+{\frac {h}{2}}k_{1},\ t+{\frac {h}{2}}\right)\\&=f\left(y_{t},\ t\right)+{\frac {h}{2}}{\frac {d}{dt}}f\left(y_{t},\ t\right)\\k_{3}&=f\left(y_{t+h/2}^{2},\ t+{\frac {h}{2}}\right)=f\left(y_{t}+{\frac {h}{2}}f\left(y_{t}+{\frac {h}{2}}k_{1},\ t+{\frac {h}{2}}\right),\ t+{\frac {h}{2}}\right)\\&=f\left(y_{t},\ t\right)+{\frac {h}{2}}{\frac {d}{dt}}\left[f\left(y_{t},\ t\right)+{\frac {h}{2}}{\frac {d}{dt}}f\left(y_{t},\ t\right)\right]\\k_{4}&=f\left(y_{t+h}^{3},\ t+h\right)=f\left(y_{t}+hf\left(y_{t}+{\frac {h}{2}}k_{2},\ t+{\frac {h}{2}}\right),\ t+h\right)\\&=f\left(y_{t}+hf\left(y_{t}+{\frac {h}{2}}f\left(y_{t}+{\frac {h}{2}}f\left(y_{t},\ t\right),\ t+{\frac {h}{2}}\right),\ t+{\frac {h}{2}}\right),\ t+h\right)\\&=f\left(y_{t},\ t\right)+h{\frac {d}{dt}}\left[f\left(y_{t},\ t\right)+{\frac {h}{2}}{\frac {d}{dt}}\left[f\left(y_{t},\ t\right)+{\frac {h}{2}}{\frac {d}{dt}}f\left(y_{t},\ t\right)\right]\right]\end{aligned}}}$

trong đó:

${\displaystyle {\frac {d}{dt}}f(y_{t},\ t)={\frac {\partial }{\partial y}}f(y_{t},\ t){\dot {y}}_{t}+{\frac {\partial }{\partial t}}f(y_{t},\ t)=f_{y}(y_{t},\ t){\dot {y}}+f_{t}(y_{t},\ t):={\ddot {y}}_{t}}$

là đạo hàm toàn phần của ${\displaystyle f}$ theo thời gian.

Nếu bây giờ chúng ta biểu diễn công thức tổng quát bằng những gì chúng ta vừa rút ra, chúng ta có được:

${\displaystyle {\begin{aligned}y_{t+h}&=y_{t}+h\left\lbrace a\cdot f(y_{t},\ t)+b\cdot \left[f\left(y_{t},\ t\right)+{\frac {h}{2}}{\frac {d}{dt}}f\left(y_{t},\ t\right)\right]\right.+\\&{}+c\cdot \left[f\left(y_{t},\ t\right)+{\frac {h}{2}}{\frac {d}{dt}}\left[f\left(y_{t},\ t\right)+{\frac {h}{2}}{\frac {d}{dt}}f\left(y_{t},\ t\right)\right]\right]+\\&{}+d\cdot \left[f\left(y_{t},\ t\right)+h{\frac {d}{dt}}\left[f\left(y_{t},\ t\right)+{\frac {h}{2}}{\frac {d}{dt}}\left[f\left(y_{t},\ t\right)+\left.{\frac {h}{2}}{\frac {d}{dt}}f\left(y_{t},\ t\right)\right]\right]\right]\right\rbrace +{\mathcal {O}}(h^{5})\\&=y_{t}+a\cdot hf_{t}+b\cdot hf_{t}+b\cdot {\frac {h^{2}}{2}}{\frac {df_{t}}{dt}}+c\cdot hf_{t}+c\cdot {\frac {h^{2}}{2}}{\frac {df_{t}}{dt}}+\\&{}+c\cdot {\frac {h^{3}}{4}}{\frac {d^{2}f_{t}}{dt^{2}}}+d\cdot hf_{t}+d\cdot h^{2}{\frac {df_{t}}{dt}}+d\cdot {\frac {h^{3}}{2}}{\frac {d^{2}f_{t}}{dt^{2}}}+d\cdot {\frac {h^{4}}{4}}{\frac {d^{3}f_{t}}{dt^{3}}}+{\mathcal {O}}(h^{5})\end{aligned}}}$

và so sánh với chuỗi Taylor của ${\displaystyle y_{t+h}}$ xung quanh ${\displaystyle y_{t}}$:

${\displaystyle {\begin{aligned}y_{t+h}&=y_{t}+h{\dot {y}}_{t}+{\frac {h^{2}}{2}}{\ddot {y}}_{t}+{\frac {h^{3}}{6}}y_{t}^{(3)}+{\frac {h^{4}}{24}}y_{t}^{(4)}+{\mathcal {O}}(h^{5})=\\&=y_{t}+hf(y_{t},\ t)+{\frac {h^{2}}{2}}{\frac {d}{dt}}f(y_{t},\ t)+{\frac {h^{3}}{6}}{\frac {d^{2}}{dt^{2}}}f(y_{t},\ t)+{\frac {h^{4}}{24}}{\frac {d^{3}}{dt^{3}}}f(y_{t},\ t)\end{aligned}}}$

chúng ta có được một hệ phương trình của các hệ số:

${\displaystyle {\begin{cases}&a+b+c+d=1\\[6pt]&{\frac {1}{2}}b+{\frac {1}{2}}c+d={\frac {1}{2}}\\[6pt]&{\frac {1}{4}}c+{\frac {1}{2}}d={\frac {1}{6}}\\[6pt]&{\frac {1}{4}}d={\frac {1}{24}}\end{cases}}}$

giải hệ trên ta có ${\displaystyle a={\frac {1}{6}},b={\frac {1}{3}},c={\frac {1}{3}},d={\frac {1}{6}}}$ như đã nói trên.

Xem thêm

• Phương pháp Euler
• Danh sách các phương pháp Runge–Kutta
• Phân tích phương trình vi phân thường bằng giải tích số
• PottersWheel – Hiệu chỉnh tham số trong các hệ phương trình vi phân thường (ODE) bằng cách sử dụng tích phân Runge-Kutta ẩn
• Phương pháp Runge–Kutta (SDE)
• Các phương pháp tuyến tính tổng quát

Tham khảo

Chú thích

Sách

• Ascher, Uri M.; Petzold, Linda R. (1998), Computer Methods for Ordinary Differential Equations and Differential-Algebraic Equations, Philadelphia: Society for Industrial and Applied Mathematics, ISBN 978-0-89871-412-8.
• Atkinson, Kendall A. (1989), An Introduction to Numerical Analysis (ấn bản 2), New York: John Wiley & Sons, ISBN 978-0-471-50023-0.
• Butcher, John C. (tháng 5 năm 1963), Coefficients for the study of Runge-Kutta integration processes, 3, tr. 185–201, doi:10.1017/S1446788700027932.
• Butcher, John C. (1975), “A stability property of implicit Runge-Kutta methods”, BIT, 15: 358–361, doi:10.1007/bf01931672.
• Butcher, John C. (2008), Numerical Methods for Ordinary Differential Equations, New York: John Wiley & Sons, ISBN 978-0-470-72335-7.
• Cellier, F.; Kofman, E. (2006), Continuous System Simulation, Springer Verlag, ISBN 0-387-26102-8.
• Dahlquist, Germund (1963), “A special stability problem for linear multistep methods”, BIT, 3: 27–43, doi:10.1007/BF01963532, ISSN 0006-3835.
• Forsythe, George E.; Malcolm, Michael A.; Moler, Cleve B. (1977), Computer Methods for Mathematical Computations, Prentice-Hall (see Chapter 6).
• Hairer, Ernst; Nørsett, Syvert Paul; Wanner, Gerhard (1993), Solving ordinary differential equations I: Nonstiff problems, Berlin, New York: Springer-Verlag, ISBN 978-3-540-56670-0.
• Hairer, Ernst; Wanner, Gerhard (1996), Solving ordinary differential equations II: Stiff and differential-algebraic problems (ấn bản 2), Berlin, New York: Springer-Verlag, ISBN 978-3-540-60452-5.
• Iserles, Arieh (1996), A First Course in the Numerical Analysis of Differential Equations, Cambridge University Press, ISBN 978-0-521-55655-2.
• Lambert, J.D (1991), Numerical Methods for Ordinary Differential Systems. The Initial Value Problem, John Wiley & Sons, ISBN 0-471-92990-5
• Kaw, Autar; Kalu, Egwu (2008), Numerical Methods with Applications (ấn bản 1), autarkaw.com.
• Kutta, Martin Wilhelm (1901), “Beitrag zur näherungsweisen Integration totaler Differentialgleichungen”, Zeitschrift für Mathematik und Physik, 46: 435–453.
• Press, William H.; Flannery, Brian P.; Teukolsky, Saul A.; Vetterling, William T. (2007), “Section 17.1 Runge-Kutta Method”, Numerical Recipes: The Art of Scientific Computing (ấn bản 3), Cambridge University Press, ISBN 978-0-521-88068-8, Bản gốc lưu trữ ngày 11 tháng 8 năm 2011, truy cập ngày 29 tháng 9 năm 2016. Also, Section 17.2. Adaptive Stepsize Control for Runge-Kutta Lưu trữ 2011-08-11 tại Wayback Machine.
• Stoer, Josef; Bulirsch, Roland (2002), Introduction to Numerical Analysis (ấn bản 3), Berlin, New York: Springer-Verlag, ISBN 978-0-387-95452-3.
• Süli, Endre; Mayers, David (2003), An Introduction to Numerical Analysis, Cambridge University Press, ISBN 0-521-00794-1.
• Tan, Delin; Chen, Zheng (2012), “On A General Formula of Fourth Order Runge-Kutta Method” (PDF), Journal of Mathematical Science & Mathematics Education, 7.2: 1–10.

Liên kết ngoài

• Hazewinkel, Michiel biên tập (2001), “Runge-Kutta method”, Bách khoa toàn thư Toán học, Springer, ISBN 978-1-55608-010-4
• Runge–Kutta 4th-Order Method
• Runge Kutta Method for O.D.E.'s
• DotNumerics: Ordinary Differential Equations for C# and VB.NET — Initial-value problem for nonstiff and stiff ordinary differential equations (explicit Runge–Kutta, implicit Runge–Kutta, Gear's BDF and Adams–Moulton).

Bản mẫu:Numerical integrators