Методи на Рунге-Кута

В математическия анализ, методите на Рунге-Кута са важна част от имплицитните и експлицитните итеративни методи за получаване на приближено решение на обикновени диференциални уравнения. Тези похвати са били въведени около 1900 г. от немските математици Карл Рунге и Мартин Вилхелм Кута.

Метод на Рунге-Кута от 4-ти ред

Методът на Рунге-Кута от 4-ти ред е най-често използваният, също известен като „РК“, „Класически метод на Рунге-Кута“ или просто „Метод на Рунге-Кута“.

Дадена е следната начална задача (още се нарича „задача на Коши“):

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

Казано с думи, това означава, че стойностите, които получава y са функции на y и на t(time). В началото времето е $t_{0}$, a y е $y_{0}$. В уравнението y може да е скаларна или векторна променлива.

Методът на РК от 4-ти ред за тази задача е получен от следните формули:

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

където $y_{n+1}$ приближението с РК4 на $y(t_{n+1})$, и

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

(Забележка: горните уравнения имат различни, но равнозначни дефиниции в различните текстове).

Така новата стойност ($y_{n+1}$) се пресмята с помощта на текущата стойност($y_{n}$) плюс стойност на четирите нараствания с определени тегла, където всяко нарастване е произведение от големината на интервала h и оценения наклон, зададен от функцията f от дясната страна на диференциалното уравнение.

• $k_{1}$ е нарастване, основаващо се на наклона от началото на интервала, използвайки $y_{n}$, Метод на (Ойлер);
• $k_{2}$ е нарастване, основаващо се на наклона от междинна точка от интервала, използвайки $y_{n}+{\tfrac {1}{2}}k_{1}$;
• $k_{3}$ е отново нарастване, основаващо се на налкона от междинна точка, но използвайки $y_{n}+{\tfrac {1}{2}}k_{2}$;
• $k_{4}$ е нарастване, основаващо се на наклона в края на интервала, използвайки $y_{n}+k_{3}$.

При усредняване на четирите нараствания се дава по-голямо тегло на нарастванията в междинните точки. Теглата са избрани така, че ако $f$ не зависи от $y$ и диференциалното уравнение има пръв интеграл, тогава РК4 e „правило на Симпсън“.

РК4 е метод от 4-ти ред, което означава, че грешката на всяка стъпка е от порядъка на $h^{5}$, докато крайната натрупана грешка има ред $h^{4}$.

Явни методи на Рунге-Кута

Семейството на експлицитните методи на Рунге-Кута са обобщение на РК4 метода, споменат преди малко. Те са представени чрез:

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

където

$k_{1}=hf(t_{n},y_{n}),\,$

$k_{2}=hf(t_{n}+c_{2}h,y_{n}+a_{21}k_{1}),\,$

$k_{3}=hf(t_{n}+c_{3}h,y_{n}+a_{31}k_{1}+a_{32}k_{2}),\,$

$\vdots$

$k_{s}=hf(t_{n}+c_{s}h,y_{n}+a_{s1}k_{1}+a_{s2}k_{2}+\cdots +a_{s,s-1}k_{s-1}).$

(Забележка: горните уравнения имат различни, но равнозначни дефиниции в различните текстове).

За да се специфицира определен метод, е нужно да се зададе цяло число s (брой етапи) и коефициентите a_ij (за 1 ≤ j < i ≤ s), b_i (за i = 1, 2, ..., s) и c_i (за i = 2, 3, ..., s). Матрицата [a_ij] се нарича матрица на Рунге-Кута, където b_i и c_i са познати като тегла и възли. Таблицата на Бутчер е следната:

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

Методът на Рунге-Кута е коректен, ако

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

Има и други съпътстващи изисквания, ако искаме метода да има определен ред p, което ознавача, че локалното грешката от закръгляване е O(h^p+1). Това може да бъде извлечено от дефиницията за грешка от закръгляване. Например, 2-етапен метод е от 2-ри ред ако b₁ + b₂ = 1, b₂c₂ = 1/2, и a₂₁ = c₂.

Примери

РК4 методът има следната структура. Неговата таблицата е както следва:

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

Най-простият метод на Рунге-Кута е методът на Ойлер, представен с формулата $y_{n+1}=y_{n}+hf(t_{n},y_{n})$. Това е единственият коректно поставен експлицитен метод на Рунге-Кута от един етап. Таблицата, която отговаря на тази формула е:

	0
		1

Методи от 2-ри ред

Методът на междинната точка представлява пример за метод от 2-ри ред с два етапа:

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

Таблицата е:

	0
	1/2	1/2
		0	1

TМетодът на междинната точка не е единствен метод на РК от 2-ри ред с два етапа. Всъщност, има семейство от такива методи, параметризирани с α и представени с формула

$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 )}.$

Таблицата на Бутчер е:

	0
	$\alpha$	$\alpha$
		$1-{\tfrac {1}{2\alpha }}$	${\tfrac {1}{2\alpha }}$

В това семейство при частния случай $\alpha ={\tfrac {1}{2}}$ получаваме метода на междинната точка, а при $\alpha =1$ метод на Хойн.

Използване

Като пример, да разгледаме метода на РК от 2-ри ред с 2 етапа с α=2/3. Той се представя с таблицата:

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

със съответните уравнения:

${\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}}$

Методът се използва за да реши началната задача:

$y'=\tan(y)+1,\quad y(1)=1,\ t\in [1,1.1]$

със стъпка h=0.025, следователно методът трябва да направи 4 стъпки.

Изпълнението на метода е както следва:

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

Числените решения отговарят на подчертаните стойности.

Адаптивни методи на Рунге-Кута

Адаптивните методи са съставени, така че да оценяват локалната грешката от закръгляване на всяка РК-стъпка. Това става използвайки два метода в таблици, единия от ред $p$, а другия от ред $p-1$.

Стъпката от по-нисък ред се получава от:

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

където $k_{i}$ са същите като за метод от по-висок ред. Тогава грешката е :

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

което е$O(h^{p})$.

Таблицата на Бутчер за тови вид метод е разширена за да се въведат стойностите на $b_{i}^{*}$:

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

Методът на Рунге-Кута-Фелберг използва два метода от 5-и и 4-ти ред. Неговата разширена таблица е:

	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

Най-простият адаптивен метод на РК включва комбинирането на метод на Хойн, който е от 2-ри ред, с метод на Ойлер, който е от 1-ви ред. Неговата разширена таблица е:

	0
	1	1
		1/2	1/2
		1	0

Оценката на грешката се използва за контролиране на големината на стъпката.

Неявни методи на Рунге-Кута

Всички споменати методи на Рунге-Кута са експлицитни. Експлицитните методи на Рунге-Кута са неподходящи за решаване на така наречените „stiff equations“ (твърди уравнения), защото обхвата на абсолютната им устойчивост е малък, по-точно ограничен. Този проблем е особено важен при решаването на частни диференциални уравнения.

Неустойчивостта на експлицитните методи подбужда създаването на имплицитните методи. Имплицитен метод на РК има формата:

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

където

$k_{i}=hf{\bigl (}t_{n}+c_{i}h,y_{n}+\sum _{j=1}^{s}a_{ij}k_{j}{\bigr )},\quad i=1,\ldots ,s.$

Разликата с експлицитния метод е, че при експлицитния метод сумата на j стига само до i – 1. Това се вижда и в таблицата. За имплицитен метод, матрицата на коефициенти $a_{ij}$ не е непременно долно-триъгълна:

${\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}\\\end{array}}={\begin{array}{c|c}\mathbf {c} &A\\\hline &\mathbf {b^{T}} \\\end{array}}$

Последиците от тази разлика са, че на всяка стъпка се решава система от алгебрични уравнения. Това значително повишава изчислителния разход. Ако метод със s етапа се използва за решаване на диференциално уравнение с m компонента, системата от алгебрични уравнения ще има ms компонента. За сравнение, при имплицитен s-стъпков линеен метод се решава система от алгебрични уравнения със само s компонента.

Примери

Най-простият имплицитен метод на РК е представен с метод на Ойлер назад:

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

Таблицата е просто:

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

Тази таблица отговаря на формулите:

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

които могат да бъдат пренаредени за да се получи формула за метод на Ойлер назад, описан по-горе.

Друг пример за имплицитен метод на Рунге-Кута е правилото на трапеца. Съответната таблица е:

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

Методите на Гаус-Льожандър формират семейство методи, базирани на квадратурата на Гаус. Метод на Гаус-Льожандър със s етапа е от ред 2s. Метод с два етапа (и 4-ти ред) има следната таблица:

${\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}}\end{array}}$

Устойчивост

Предимството на имплицитните методи на Рунге-Кута пред експлицитните е по-голямата им устойчивост, особено когато се прилагат върху твърди уравнения. Разглеждаме линейното тестово уравнение y' = λy. Методът на РК, приложен за това уравнение, свежда решаването му до изпълняване на итерациите $y_{n+1}=r(h\lambda )\,y_{n}$, където r е дадено чрез

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

където е е единичния вектор. Функцията r се нарича функция на устойчивостта. Експлицитните методи имат строга ниско-триъгълна матрица А, която предполага, че det(I − zA) = 1 и функцията на устойчивостта е полином.

Численото решение на линейното тестово уравнение се разпада до нула ако |r(z)| < 1 при z=hλ. Множеството от такива стойности на z се нарича област на абсолютна устойчивост. В частност, методът трябва да бъде A-стабилен, ако всяко z с Re(z)<0 се намира в областта на абсолютната устойчивост. Функцията на устойчивостта на експлицитния метод на РК е полином, следователно експлицитните методи на РК никога не могат да бъдат A-стабилни.

Извеждане на метода на Рунге-Кута от 4-ти ред

Методът на Рунге-Кута от ред $s$ може да бъде написан така:

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

където:

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

са нараствания получени чрез пресмятане на производните на $y_{t}$ в $i$-тия ред.

Извеждането на метода на Рунге-Кута от 4-ти ред се постига чрез използване на общата формула с $s=4$, както е описано по-горе, в началната точка, междинна точка и последната точка на всеки интервал $(t,t+h)$, затова избираме:

$\alpha _{i}$	$\beta _{i}$
$\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$

и $\beta _{ij}=0$. Започваме с определянето на следните величини:

${\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}}$

където $y_{t+h/2}^{1}={\dfrac {y_{t}+y_{t+h}^{1}}{2}}$ and $y_{t+h/2}^{2}={\dfrac {y_{t}+y_{t+h}^{2}}{2}}$ Ако задедем:

${\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}}$

и имайки предвид горните изрази можем да покажем, че следните равенства имат точност ${\mathcal {O}}(h^{2})$:

${\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}}$

където:

${\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}$

е пълната производна на $f$ по отношение на времето t.

Ако изразим общата формула, използвайки това, което изведохме по-горе, получаваме:

${\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}}$

И като сравним това с реда на Тейлър за $y_{t+h}$ около $y_{t}$:

${\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}}$

получаваме система от изисквания за коефициентите:

$\left\{{\begin{aligned}&a+b+c+d=1\\&{\frac {1}{2}}b+{\frac {1}{2}}c+d={\frac {1}{2}}\\&{\frac {1}{4}}c+{\frac {1}{2}}d={\frac {1}{6}}\\&{\frac {1}{4}}d={\frac {1}{24}}\end{aligned}}\right.$

която решена дава $a={\frac {1}{6}},b={\frac {1}{3}},c={\frac {1}{3}},d={\frac {1}{6}}$.

Литература

• Ascher, Uri M., Petzold, Linda R. Computer Methods for Ordinary Differential Equations and Differential-Algebraic Equations. Philadelphia, Society for Industrial and Applied Mathematics, 1998. ISBN 978-0-89871-412-8..
• Atkinson, Kendall A. An Introduction to Numerical Analysis. 2nd. New York, John Wiley & Sons, 1989. ISBN 978-0-471-50023-0..
• Butcher, John C. Numerical Methods for Ordinary Differential Equations. New York, John Wiley & Sons, 2003. ISBN 978-0-471-96758-3..
• Cellier, F., Kofman, E. Continuous System Simulation. Springer Verlag, 2006. ISBN 0-387-26102-8..
• Dahlquist, Germund. A special stability problem for linear multistep methods. BIT. Т. 3. 1963. DOI:10.1007/BF01963532. с. 27 – 43..
• Forsythe, George E., Malcolm, Michael A., Moler, Cleve B. Computer Methods for Mathematical Computations. Prentice-Hall, 1977. (see Chapter 6).
• Hairer, Ernst, Nørsett, Syvert Paul, Wanner, Gerhard. Solving ordinary differential equations I: Nonstiff problems. Berlin, New York, Springer-Verlag, 1993. ISBN 978-3-540-56670-0..
• Hairer, Ernst, Wanner, Gerhard. Solving ordinary differential equations II: Stiff and differential-algebraic problems. 2nd. Berlin, New York, Springer-Verlag, 1996. ISBN 978-3-540-60452-5..
• Iserles, Arieh. A First Course in the Numerical Analysis of Differential Equations. Cambridge University Press, 1996. ISBN 978-0-521-55655-2..
• Kaw, Autar, Kalu, Egwu. Numerical Methods with Applications. 1st. autarkaw.com, 2008..
• Section 17.1 Runge-Kutta Method // Numerical Recipes: The Art of Scientific Computing. 3rd. Cambridge University Press, 2007. ISBN 978-0-521-88068-8.. Also, Section 17.2. Adaptive Stepsize Control for Runge-Kutta Архив на оригинала от 2011-08-11 в Wayback Machine..
• Stoer, Josef, Bulirsch, Roland. Introduction to Numerical Analysis. 3rd. Berlin, New York, Springer-Verlag, 2002. ISBN 978-0-387-95452-3..
• Süli, Endre, Mayers, David. An Introduction to Numerical Analysis. Cambridge University Press, 2003. ISBN 0-521-00794-1..

Външни препратки

• 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).
• GafferOnGames Архив на оригинала от 2012-09-05 в Wayback Machine. – A physics resource for computer programmers
• PottersWheel – Parameter calibration in ODE systems using implicit Runge–Kutta integration

Тази страница частично или изцяло представлява превод на страницата Runge–Kutta methods в Уикипедия на английски. Оригиналният текст, както и този превод, са защитени от Лиценза „Криейтив Комънс – Признание – Споделяне на споделеното“, а за съдържание, създадено преди юни 2009 година – от Лиценза за свободна документация на ГНУ. Прегледайте историята на редакциите на оригиналната страница, както и на преводната страница, за да видите списъка на съавторите. ​ ВАЖНО: Този шаблон се отнася единствено до авторските права върху съдържанието на статията. Добавянето му не отменя изискването да се посочват конкретни източници на твърденията, които да бъдат благонадеждни.​