Методы Рунге–Кутты

Семейство неявных и явных итерационных методов

Сравнение методов Рунге-Кутты для дифференциального уравнения (красный цвет — точное решение) ${\displaystyle y'=\sin(t)^{2}\cdot y}$

В численном анализе используются методы Рунге –Кутты ( англ.: / ˈ r ʊ ŋ ə ˈ k ʊ t ɑː / ^ⓘ RUUNG-ə-KUUT-tah^[1]) представляют собой семействонеявных и явныхитерационных методов, включающееметод Эйлера, используемый привременной дискретизациидля приближенных решенийодновременных нелинейных уравнений.^[2]Эти методы были разработаны около 1900 года немецкими математикамиКарлом РунгеиВильгельмом Куттой.

Метод Рунге–Кутта

Уклоны, используемые классическим методом Рунге-Кутта

Наиболее известный член семейства методов Рунге–Кутты обычно называют «RK4», «классическим методом Рунге–Кутты» или просто «методом Рунге–Кутты».

Пусть задача с начальными данными задана следующим образом:

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

Вот неизвестная функция (скалярная или векторная) времени , которую мы хотели бы аппроксимировать; нам говорят, что , скорость, с которой изменяется, является функцией и самого себя. В начальный момент времени соответствующее значение равно . Функция и начальные условия , даны. ${\displaystyle y}$ ${\displaystyle t}$ ${\displaystyle {\frac {dy}{dt}}}$ ${\displaystyle y}$ ${\displaystyle t}$ ${\displaystyle y}$ ${\displaystyle t_{0}}$ ${\displaystyle y}$ ${\displaystyle y_{0}}$ ${\displaystyle f}$ ${\displaystyle t_{0}}$ ${\displaystyle y_{0}}$

Теперь выберем размер шага h > 0 и определим:

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

для n = 0, 1, 2, 3, ..., используя ^[3]

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

( Примечание: приведенные выше уравнения имеют разные, но эквивалентные определения в разных текстах. ^[4] )

Здесь представлена ​​аппроксимация RK4 для , а следующее значение ( ) определяется текущим значением ( ) плюс средневзвешенное значение четырех приращений, где каждое приращение является произведением размера интервала h и предполагаемого наклона, заданного функцией f в правой части дифференциального уравнения. ${\displaystyle y_{n+1}}$ ${\displaystyle y(t_{n+1})}$ ${\displaystyle y_{n+1}}$ ${\displaystyle y_{n}}$

• ${\displaystyle k_{1}}$- наклон в начале интервала, с использованием ( метода Эйлера ); ${\displaystyle y}$
• ${\displaystyle k_{2}}$— наклон в средней точке интервала, с использованием и ; ${\displaystyle y}$ ${\displaystyle k_{1}}$
• ${\displaystyle k_{3}}$снова наклон в средней точке, но теперь с использованием и ; ${\displaystyle y}$ ${\displaystyle k_{2}}$
• ${\displaystyle k_{4}}$— наклон в конце интервала с использованием и . ${\displaystyle y}$ ${\displaystyle k_{3}}$

При усреднении четырех наклонов больший вес придается наклонам в средней точке. Если не зависит от , так что дифференциальное уравнение эквивалентно простому интегралу, то RK4 — это правило Симпсона . ^[5] ${\displaystyle f}$ ${\displaystyle y}$

Метод RK4 является методом четвертого порядка, что означает, что локальная ошибка усечения составляет порядка , тогда как общая накопленная ошибка составляет порядка . ${\displaystyle O(h^{5})}$ ${\displaystyle O(h^{4})}$

Во многих практических приложениях функция независима от (так называемая автономная система или система, не зависящая от времени, особенно в физике), а их приращения вообще не вычисляются и не передаются в функцию , а используется только окончательная формула для . ${\displaystyle f}$ ${\displaystyle t}$ ${\displaystyle f}$ ${\displaystyle t_{n+1}}$

Явные методы Рунге–Кутты

Семейство явных методов Рунге–Кутты является обобщением упомянутого выше метода RK4. Оно задается формулой

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

где ^[6]

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

( Примечание: приведенные выше уравнения могут иметь разные, но эквивалентные определения в некоторых текстах. ^[4] )

Чтобы указать конкретный метод, необходимо указать целое число 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 известны как веса и узлы . ^[7] Эти данные обычно располагаются в мнемоническом устройстве, известном как таблица Бутчера (в честь Джона К. Бутчера ):

Разложение в ряд Тейлора показывает , что метод Рунге–Кутты является последовательным тогда и только тогда, когда

${\displaystyle \sum _{i=1}^{s}b_{i}=1.}$

Существуют также сопутствующие требования, если требуется, чтобы метод имел определенный порядок p , что означает, что локальная ошибка усечения равна O( h ^{p +1} ). Они могут быть выведены из определения самой ошибки усечения. Например, двухэтапный метод имеет порядок 2, если b ₁ + b ₂ = 1, b ₂ c ₂ = 1/2 и b ₂ a ₂₁ = 1/2. ^[8] Обратите внимание, что популярным условием для определения коэффициентов является ^[8]

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

Однако одно это условие не является ни достаточным, ни необходимым для согласованности. ^[9]

В общем случае, если явный -этапный метод Рунге-Кутты имеет порядок , то можно доказать, что число этапов должно удовлетворять , а если , то . ^[10] Однако неизвестно, являются ли эти границы точными во всех случаях. В некоторых случаях доказано, что граница не может быть достигнута. Например, Бутчер доказал, что для , явного метода с этапами не существует . ^[11] Бутчер также доказал, что для , явного метода Рунге-Кутты с этапами не существует . ^[12] В общем случае, однако, остается открытой проблема, каково точное минимальное число этапов для явного метода Рунге-Кутты, чтобы иметь порядок . Некоторые известные значения: ^[13] ${\displaystyle s}$ ${\displaystyle p}$ ${\displaystyle s\geq p}$ ${\displaystyle p\geq 5}$ ${\displaystyle s\geq p+1}$ ${\displaystyle p>6}$ ${\displaystyle s=p+1}$ ${\displaystyle p>7}$ ${\displaystyle p+2}$ ${\displaystyle s}$ ${\displaystyle p}$

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

Доказуемая граница выше затем подразумевает, что мы не можем найти методы порядков , которые требуют меньше стадий, чем методы, которые мы уже знаем для этих порядков. Работа Бутчера также доказывает, что методы 7-го и 8-го порядка имеют минимум 9 и 11 стадий соответственно. ^{[11] [12]} Пример явного метода порядка 6 с 7 стадиями можно найти в ^[14] . Явные методы порядка 7 с 9 стадиями ^[11] и явные методы порядка 8 с 11 стадиями ^[15] также известны. См. ^{[16] [17]} для резюме. ${\displaystyle p=1,2,\ldots ,6}$

Примеры

Метод RK4 попадает в эту структуру. Его таблица ^[18]

Небольшая вариация метода «Рунге–Кутта» также принадлежит Кутте в 1901 году и называется правилом 3/8. ^[19] Основное преимущество этого метода в том, что почти все коэффициенты ошибок меньше, чем в популярном методе, но он требует немного больше FLOP (операций с плавающей точкой) на шаг времени. Его таблица Бутчера

Однако простейшим методом Рунге–Кутты является (прямой) метод Эйлера , заданный формулой . Это единственный последовательный явный метод Рунге–Кутты с одним этапом. Соответствующая таблица имеет вид ${\displaystyle y_{n+1}=y_{n}+hf(t_{n},y_{n})}$

Методы второго порядка с двумя этапами

Примером метода второго порядка с двумя этапами является явный метод средней точки :

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

Соответствующая таблица:

Метод средней точки — не единственный метод Рунге–Кутты второго порядка с двумя этапами; существует семейство таких методов, параметризованных α и заданных формулой ^[20]

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

Его картина Мясника - это

В этом семействе есть метод средней точки , есть метод Хойна , ^[5] и есть метод Ралстона. ${\displaystyle \alpha ={\tfrac {1}{2}}}$ ${\displaystyle \alpha =1}$ ${\displaystyle \alpha ={\tfrac {2}{3}}}$

Использовать

В качестве примера рассмотрим двухэтапный метод Рунге-Кутты второго порядка с α = 2/3, также известный как метод Ральстона . Он задается таблицей

с соответствующими уравнениями

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

Этот метод используется для решения задачи начального значения

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

с размером шага h = 0,025, поэтому метод должен выполнить четыре шага.

Метод осуществляется следующим образом:

Численные решения соответствуют подчеркнутым значениям.

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

Адаптивные методы предназначены для получения оценки локальной ошибки усечения одного шага Рунге-Кутты. Это делается с помощью двух методов, одного с порядком и одного с порядком . Эти методы переплетены, т. е. имеют общие промежуточные шаги. Благодаря этому оценка ошибки имеет небольшие или пренебрежимо малые вычислительные затраты по сравнению с шагом с методом более высокого порядка. ${\displaystyle p}$ ${\displaystyle p-1}$

Во время интегрирования размер шага адаптируется таким образом, чтобы предполагаемая ошибка оставалась ниже порогового значения, определенного пользователем: если ошибка слишком велика, шаг повторяется с меньшим размером шага; если ошибка намного меньше, размер шага увеличивается для экономии времени. Это приводит к (почти) оптимальному размеру шага, что экономит время вычислений. Более того, пользователю не нужно тратить время на поиск подходящего размера шага.

Шаг низшего порядка определяется как

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

где те же, что и для метода высшего порядка. Тогда ошибка равна ${\displaystyle k_{i}}$

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

что равно . Таблица Бутчера для этого вида метода расширена, чтобы дать значения : ${\displaystyle O(h^{p})}$ ${\displaystyle b_{i}^{*}}$

Метод Рунге –Кутты–Фельберга имеет два метода порядков 5 и 4. Его расширенная таблица Бутчера имеет вид:

Однако простейший адаптивный метод Рунге–Кутты включает в себя объединение метода Хойна , который имеет порядок 2, с методом Эйлера , который имеет порядок 1. Его расширенная таблица Бутчера имеет вид:

Другими адаптивными методами Рунге–Кутты являются метод Богацки–Шампайна (порядки 3 и 2), метод Кэша–Карпа и метод Дорманда–Принса (оба с порядками 5 и 4).

Неконфлюэнтные методы Рунге–Кутты

Метод Рунге–Кутты называется неконфлюэнтным ^[21], если все различны. ${\displaystyle c_{i},\,i=1,2,\ldots ,s}$

Методы Рунге–Кутты–Нистрема.

Методы Рунге–Кутты–Нистрёма являются специализированными методами Рунге–Кутты, оптимизированными для дифференциальных уравнений второго порядка. ^{[22] [23]} Общий метод Рунге–Кутты–Нистрёма для системы ОДУ второго порядка

${\displaystyle {\ddot {y}}_{i}=f_{i}(y_{1},y_{2},\cdots ,y_{n})}$

с заказом связано с формой ${\displaystyle s}$

${\displaystyle {\begin{cases}g_{i}=y_{m}+c_{i}h{\dot {y}}_{m}+h^{2}\sum _{j=1}^{s}a_{ij}f(g_{j}),&i=1,2,\cdots ,s\\y_{m+1}=y_{m}+h{\dot {y}}_{m}+h^{2}\sum _{j=1}^{s}{\bar {b}}_{j}f(g_{j})\\{\dot {y}}_{m+1}={\dot {y}}_{m}+h\sum _{j=1}^{s}b_{j}f(g_{j})\end{cases}}}$

который образует стол Мясника с формой

${\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 &{\bar {b}}_{1}&{\bar {b}}_{2}&\dots &{\bar {b}}_{s}\\&b_{1}&b_{2}&\dots &b_{s}\\\end{array}}={\begin{array}{c|c}\mathbf {c} &\mathbf {A} \\\hline &\mathbf {\bar {b}} ^{\top }\\&\mathbf {b} ^{\top }\end{array}}}$.

Два явных метода РКН четвертого порядка задаются следующими таблицами Бутчера:

${\displaystyle {\begin{array}{c|ccc}c_{i}&&a_{ij}&\\{\frac {3+{\sqrt {3}}}{6}}&0&0&0\\{\frac {3-{\sqrt {3}}}{6}}&{\frac {2-{\sqrt {3}}}{12}}&0&0\\{\frac {3+{\sqrt {3}}}{6}}&0&{\frac {\sqrt {3}}{6}}&0\\\hline {\overline {b_{i}}}&{\frac {5-3{\sqrt {3}}}{24}}&{\frac {3+{\sqrt {3}}}{12}}&{\frac {1+{\sqrt {3}}}{24}}\\\hline b_{i}&{\frac {3-2{\sqrt {3}}}{12}}&{\frac {1}{2}}&{\frac {3+2{\sqrt {3}}}{12}}\end{array}}}$ ${\displaystyle {\begin{array}{c|ccc}c_{i}&&a_{ij}&\\{\frac {3-{\sqrt {3}}}{6}}&0&0&0\\{\frac {3+{\sqrt {3}}}{6}}&{\frac {2+{\sqrt {3}}}{12}}&0&0\\{\frac {3-{\sqrt {3}}}{6}}&0&-{\frac {\sqrt {3}}{6}}&0\\\hline {\overline {b_{i}}}&{\frac {5+3{\sqrt {3}}}{24}}&{\frac {3-{\sqrt {3}}}{12}}&{\frac {1-{\sqrt {3}}}{24}}\\\hline b_{i}&{\frac {3+2{\sqrt {3}}}{12}}&{\frac {1}{2}}&{\frac {3-2{\sqrt {3}}}{12}}\end{array}}}$

Эти две схемы также обладают симплектическими свойствами сохранения, когда исходное уравнение выводится из консервативной классической механической системы, т.е. когда

${\displaystyle f_{i}(x_{1},\cdots ,x_{n})={\frac {\partial V}{\partial x_{i}}}(x_{1},\cdots ,x_{n})}$

для некоторой скалярной функции . ^[24] ${\displaystyle V}$

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

Все методы Рунге–Кутты, упомянутые до сих пор, являются явными методами . Явные методы Рунге–Кутты, как правило, не подходят для решения жестких уравнений , поскольку их область абсолютной устойчивости мала; в частности, она ограничена. ^[25] Этот вопрос особенно важен при решении уравнений в частных производных .

Нестабильность явных методов Рунге–Кутты мотивирует разработку неявных методов. Неявный метод Рунге–Кутты имеет вид

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

где

${\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.}$ ^[26]

Разница с явным методом заключается в том, что в явном методе сумма по j увеличивается только до i − 1. Это также проявляется в таблице Бутчера: матрица коэффициентов явного метода является нижней треугольной. В неявном методе сумма по j увеличивается до s , а матрица коэффициентов не является треугольной, что дает таблицу Бутчера вида ^[18] ${\displaystyle a_{ij}}$

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

Объяснение этой строки см. в разделе «Адаптивные методы Рунге-Кутты» выше. ${\displaystyle b^{*}}$

Следствием этого различия является то, что на каждом шаге необходимо решить систему алгебраических уравнений. Это значительно увеличивает вычислительные затраты. Если метод с s этапами используется для решения дифференциального уравнения с m компонентами, то система алгебраических уравнений имеет ms компонентов. Это можно противопоставить неявным линейным многошаговым методам (другое большое семейство методов для ОДУ): неявный s -шаговый линейный многошаговый метод должен решить систему алгебраических уравнений только с m компонентами, поэтому размер системы не увеличивается с увеличением числа шагов. ^[27]

Примеры

Простейшим примером неявного метода Рунге–Кутты является обратный метод Эйлера :

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

Таблица Мясника для этого проста:

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

Эта таблица Бутчера соответствует формулам

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

которую можно переформулировать, чтобы получить формулу для обратного метода Эйлера, указанную выше.

Другим примером неявного метода Рунге–Кутты является правило трапеций . Его таблица Бутчера имеет вид:

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

Правило трапеции — это метод коллокации (как обсуждалось в этой статье). Все методы коллокации — это неявные методы Рунге–Кутты, но не все неявные методы Рунге–Кутты являются методами коллокации. ^[28]

Методы Гаусса –Лежандра образуют семейство методов коллокации, основанных на квадратуре Гаусса . Метод Гаусса–Лежандра с s этапами имеет порядок 2 s (таким образом, могут быть построены методы с произвольно высоким порядком). ^[29] Метод с двумя этапами (и, следовательно, с порядком четыре) имеет таблицу Бутчера:

${\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}}}$ ^[27]

Стабильность

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

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

где e обозначает вектор единиц. Функция r называется функцией устойчивости . ^[31] Из формулы следует, что r является частным двух полиномов степени s, если метод имеет s этапов. Явные методы имеют строго нижнюю треугольную матрицу A , что подразумевает, что det( I − zA ) = 1 и что функция устойчивости является полиномом. ^[32]

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

Если метод имеет порядок p , то функция устойчивости удовлетворяет условию . Таким образом, интересно изучить отношения полиномов заданных степеней, которые наилучшим образом приближают экспоненциальную функцию. Они известны как аппроксимации Паде . Аппроксимация Паде с числителем степени m и знаменателем степени n является A-устойчивой тогда и только тогда, когда m ≤ n ≤ m + 2. ^[33] ${\displaystyle r(z)={\textrm {e}}^{z}+O(z^{p+1})}$ ${\displaystyle z\to 0}$

Метод Гаусса–Лежандра с s этапами имеет порядок 2 s , поэтому его функция устойчивости является аппроксимацией Паде с m = n = s . Отсюда следует, что метод является A-устойчивым. ^[34] Это показывает, что A-устойчивый метод Рунге–Кутты может иметь произвольно высокий порядок. Напротив, порядок A-устойчивых линейных многошаговых методов не может превышать двух. ^[35]

B-стабильность

Концепция A-устойчивости для решения дифференциальных уравнений связана с линейным автономным уравнением . Далквист предложил исследование устойчивости численных схем применительно к нелинейным системам, удовлетворяющим условию монотонности. Соответствующие концепции были определены как G-устойчивость для многошаговых методов (и связанных с ними одноногих методов) и B-устойчивость (Butcher, 1975) для методов Рунге–Кутты. Метод Рунге–Кутты, примененный к нелинейной системе , который проверяет , называется B-устойчивым , если это условие подразумевает для двух численных решений. ${\displaystyle y'=\lambda y}$ ${\displaystyle y'=f(y)}$ ${\displaystyle \langle f(y)-f(z),\ y-z\rangle \leq 0}$ ${\displaystyle \|y_{n+1}-z_{n+1}\|\leq \|y_{n}-z_{n}\|}$

Пусть , и — три матрицы, определяемые с помощью Метод Рунге–Кутты называется алгебраически устойчивым ^[36], если матрицы и обе неотрицательно определены. Достаточным условием B-устойчивости ^[37] является: и неотрицательно определены. ${\displaystyle B}$ ${\displaystyle M}$ ${\displaystyle Q}$ ${\displaystyle s\times s}$ $${\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}.}$$ ${\displaystyle B}$ ${\displaystyle M}$ ${\displaystyle B}$ ${\displaystyle Q}$

Вывод метода Рунге–Кутты четвертого порядка

В общем случае метод порядка Рунге–Кутты можно записать так: ${\displaystyle s}$

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

где:

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

являются приращениями, полученными путем оценки производных в -ом порядке. ${\displaystyle y_{t}}$ ${\displaystyle i}$

Мы разрабатываем вывод ^[38] для метода Рунге–Кутты четвертого порядка, используя общую формулу с оценками, как объяснено выше, в начальной точке, средней точке и конечной точке любого интервала ; таким образом, мы выбираем: ${\displaystyle s=4}$ ${\displaystyle (t,\ t+h)}$

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

и в противном случае. Начнем с определения следующих величин: ${\displaystyle \beta _{ij}=0}$

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

где и . Если мы определим: ${\displaystyle y_{t+h/2}^{1}={\dfrac {y_{t}+y_{t+h}^{1}}{2}}}$ ${\displaystyle y_{t+h/2}^{2}={\dfrac {y_{t}+y_{t+h}^{2}}{2}}}$

${\displaystyle {\begin{aligned}k_{1}&=f(y_{t},\ t)\\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)\\k_{3}&=f\left(y_{t+h/2}^{2},\ t+{\frac {h}{2}}\right)=f\left(y_{t}+{\frac {h}{2}}k_{2},\ t+{\frac {h}{2}}\right)\\k_{4}&=f\left(y_{t+h}^{3},\ t+h\right)=f\left(y_{t}+hk_{3},\ t+h\right)\end{aligned}}}$

и для предыдущих соотношений можно показать, что справедливы следующие равенства : где: — полная производная по времени. ${\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}}}$$ $${\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}}_{t}+f_{t}(y_{t},\ t):={\ddot {y}}_{t}}$$ ${\displaystyle f}$

Если теперь выразить общую формулу, используя то, что мы только что вывели, то получим: $${\displaystyle {\begin{aligned}y_{t+h}={}&y_{t}+h\left\lbrace a\cdot f(y_{t},\ t)+b\cdot \left[f(y_{t},\ t)+{\frac {h}{2}}{\frac {d}{dt}}f(y_{t},\ t)\right]\right.+\\&{}+c\cdot \left[f(y_{t},\ t)+{\frac {h}{2}}{\frac {d}{dt}}\left[f\left(y_{t},\ t\right)+{\frac {h}{2}}{\frac {d}{dt}}f(y_{t},\ t)\right]\right]+\\&{}+d\cdot \left[f(y_{t},\ t)+h{\frac {d}{dt}}\left[f(y_{t},\ t)+{\frac {h}{2}}{\frac {d}{dt}}\left[f(y_{t},\ t)+\left.{\frac {h}{2}}{\frac {d}{dt}}f(y_{t},\ t)\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}}}$$

и сравнивая это с рядом Тейлора примерно : ${\displaystyle y_{t+h}}$ ${\displaystyle 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}}}$$

получаем систему ограничений на коэффициенты:

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

что при решении дает то, что указано выше. ${\displaystyle a={\frac {1}{6}},b={\frac {1}{3}},c={\frac {1}{3}},d={\frac {1}{6}}}$

Смотрите также

• Метод Эйлера
• Список методов Рунге–Кутты
• Численные методы для обыкновенных дифференциальных уравнений
• Метод Рунге–Кутты (SDE)
• Общие линейные методы
• Интегратор группы Ли

Примечания

1. "Метод Рунге-Кутты". Dictionary.com . Получено 4 апреля 2021 г. .
2. DEVRIES, Paul L.; HASBUN, Javier E. Первый курс вычислительной физики. Второе издание. Jones and Bartlett Publishers: 2011. стр. 215.
3. Пресс и др. 2007, с. 908; Сюли и Майерс 2003, с. 328
4. Atkinson (1989, стр. 423), Hairer, Nørsett & Wanner (1993, стр. 134), Kaw & Kalu (2008, §8.4) и Stoer & Bulirsch (2002, стр. 476) не учитывают фактор h в определении стадий. Ascher & Petzold (1998, стр. 81), Butcher (2008, стр. 93) и Iserles (1996, стр. 38) используют значения y в качестве стадий.
5. Süli & Mayers 2003, с. 328
6. Пресс и др. 2007, стр. 907
7. Исерлес 1996, стр. 38
8. Iserles 1996, стр. 39
9. В качестве контрпримера рассмотрим любую явную 2-этапную схему Рунге-Кутты с и и выбранными случайным образом. Этот метод является последовательным и (в общем случае) сходится в первом порядке. С другой стороны, 1-этапный метод с является непоследовательным и не сходится, хотя он тривиально утверждает, что . ${\displaystyle b_{1}=b_{2}=1/2}$ ${\displaystyle c_{1}}$ ${\displaystyle a_{21}}$ ${\displaystyle b_{1}=1/2}$ ${\displaystyle \sum _{j=1}^{i-1}a_{ij}=c_{i}{\text{ for }}i=2,\ldots ,s.}$
10. Мясник 2008, стр. 187
11. Butcher 1965, стр. 408
12. Butcher 1985
13. Мясник 2008, стр. 187–196.
14. Мясник 1964
15. Кертис 1970, стр. 268
16. Хайрер, Норсетт и Ваннер 1993, стр. 179
17. Мясник 1996, стр. 247
18. Süli & Mayers 2003, стр. 352
19. Хайрер, Норсетт и Ваннер (1993, стр. 138) относятся к Кутте (1901).
20. Süli & Mayers 2003, стр. 327
21. Ламберт 1991, стр. 278
22. Дорманд, Дж. Р.; Принс, П. Дж. (октябрь 1978 г.). «Новые алгоритмы Рунге–Кутты для численного моделирования в динамической астрономии». Небесная механика . 18 (3): 223–232. Bibcode : 1978CeMec..18..223D. doi : 10.1007/BF01230162. S2CID  120974351.
23. Fehlberg, E. (октябрь 1974 г.). Классические формулы Рунге–Кутты–Нюстрёма седьмого, шестого и пятого порядка с контролем размера шага для общих дифференциальных уравнений второго порядка (Отчет) (NASA TR R-432 ред.). Marshall Space Flight Center, AL: National Aeronautics and Space Administration.
24. Цинь, Мэн-Чжао; Чжу, Вэнь-Цзе (1991-01-01). «Канонические методы Рунге-Кутты-Нистрёма (RKN) для обыкновенных дифференциальных уравнений второго порядка». Компьютеры и математика с приложениями . 22 (9): 85–95. doi :10.1016/0898-1221(91)90209-M. ISSN  0898-1221.
25. Süli & Mayers 2003, стр. 349–351.
26. Изерлес 1996, с. 41; Сули и Майерс 2003, стр. 351–352.
27. Süli & Mayers 2003, с. 353
28. Исерлес 1996, стр. 43–44
29. Исерлес 1996, стр. 47
30. Hairer & Wanner 1996, стр. 40–41.
31. Хайрер и Ваннер 1996, стр. 40
32. Iserles 1996, стр. 60
33. Исерлес 1996, стр. 62–63
34. Исерлес 1996, стр. 63
35. Этот результат принадлежит Далквисту (1963).
36. Ламберт 1991, стр. 275
37. Ламберт 1991, стр. 274
38. Лю, Лин-Сяо (август 2016 г.). "Приложение C. Вывод формул численного интегрирования" (PDF) . Численное моделирование космической плазмы (I) Lecture Notes . Институт космических наук, Национальный центральный университет . Получено 17 апреля 2022 г. .

Ссылки

• Рунге, Карл Давид Толме (1895), «Über die numerische Auflösung von Differentialgleichungen», Mathematische Annalen , 46 (2), Springer : 167–178, doi : 10.1007/BF01446807, S2CID  119924854.
• Кутта, Вильгельм (1901), «Beitrag zur näherungsweisen Integration Totaler Differentialgleichungen», Zeitschrift für Mathematik und Physik , 46 : 435–453.
• Ашер, Ури М.; Петцольд, Линда Р. (1998), Компьютерные методы решения обыкновенных дифференциальных уравнений и дифференциально-алгебраических уравнений , Филадельфия: Общество промышленной и прикладной математики , ISBN 978-0-89871-412-8.
• Аткинсон, Кендалл А. (1989), Введение в численный анализ (2-е изд.), Нью-Йорк: John Wiley & Sons , ISBN 978-0-471-50023-0.
• Бутчер, Джон К. (май 1963 г.), «Коэффициенты для изучения процессов интегрирования Рунге-Кутты», Журнал Австралийского математического общества , 3 (2): 185–201, doi : 10.1017/S1446788700027932.
• Бутчер, Джон К. (май 1964 г.), «О процессах Рунге-Кутты высокого порядка», Журнал Австралийского математического общества , 4 (2): 179–194, doi : 10.1017/S1446788700023387
• Бутчер, Джон К. (1975), «Свойство устойчивости неявных методов Рунге-Кутты», BIT , 15 (4): 358–361, doi :10.1007/bf01931672, S2CID  120854166.
• Бутчер, Джон К. (2000), «Численные методы для обыкновенных дифференциальных уравнений в 20 веке», J. Comput. Appl. Math. , 125 (1–2): 1–29, Bibcode : 2000JCoAM.125....1B, doi : 10.1016/S0377-0427(00)00455-6.
• Бутчер, Джон К. (2008), Численные методы решения обыкновенных дифференциальных уравнений , Нью-Йорк: John Wiley & Sons , ISBN 978-0-470-72335-7.
• Селье, Ф.; Кофман, Э. (2006), Непрерывное системное моделирование , Springer Verlag , ISBN 0-387-26102-8.
• Дальквист, Джермунд (1963), «Специальная проблема устойчивости для линейных многошаговых методов», BIT , 3 : 27–43, doi : 10.1007/BF01963532, hdl : 10338.dmlcz/103497 , ISSN  0006-3835, S2CID  120241743.
• Форсайт, Джордж Э.; Малкольм, Майкл А.; Молер, Клив Б. (1977), Компьютерные методы математических вычислений , Prentice-Hall(см. главу 6).
• Хайрер, Эрнст; Норсетт, Сиверт Пол; Ваннер, Герхард (1993), Решение обыкновенных дифференциальных уравнений I: Нежесткие задачи , Берлин, Нью-Йорк: Springer-Verlag , ISBN 978-3-540-56670-0.
• Хайрер, Эрнст; Ваннер, Герхард (1996), Решение обыкновенных дифференциальных уравнений II: Жесткие и дифференциально-алгебраические проблемы (2-е изд.), Берлин, Нью-Йорк: Springer-Verlag , ISBN 978-3-540-60452-5.
• Исерлес, Арье (1996), Первый курс численного анализа дифференциальных уравнений , Cambridge University Press , ISBN 978-0-521-55655-2.
• Ламберт, Дж. Д. (1991), Численные методы для обыкновенных дифференциальных систем. Задача начального значения , John Wiley & Sons , ISBN 0-471-92990-5
• Кау, Аутар; Калу, Эгву (2008), Численные методы и их приложения (1-е изд.), autarkaw.com.
• Press, William H.; Teukolsky, Saul A .; Vetterling, William T.; Flannery, Brian P. (2007), «Раздел 17.1 Метод Рунге-Кутты», Numerical Recipes: The Art of Scientific Computing (3-е изд.), Cambridge University Press , ISBN 978-0-521-88068-8. Также, Раздел 17.2. Адаптивное управление размером шага для Рунге-Кутты.
• Стер, Йозеф; Булирш, Роланд (2002), Введение в численный анализ (3-е изд.), Берлин, Нью-Йорк: Springer-Verlag , ISBN 978-0-387-95452-3.
• Сюли, Эндре; Майерс, Дэвид (2003), Введение в численный анализ , Cambridge University Press , ISBN 0-521-00794-1.
• Тан, Делин; Чэнь, Чжэн (2012), «Об общей формуле метода Рунге-Кутты четвертого порядка» (PDF) , Журнал математической науки и математического образования , 7 (2): 1–10.
• Справочник по дискретной математике ignou (код mcs033)
• Джон К. Бутчер: «B-Series: Algebraic Analysis of Numerical Methods», Springer (SSCM, том 55), ISBN 978-3030709556 (апрель 2021 г.). 
• Бутчер, Дж. К. (1985), «Несуществование десятиэтапных явных методов Рунге-Кутты восьмого порядка», BIT Numerical Mathematics , 25 (3): 521–540, doi :10.1007/BF01935372.
• Бутчер, Дж. К. (1965), «О достижимом порядке методов Рунге-Кутты», Математика вычислений , 19 (91): 408–417, doi :10.1090/S0025-5718-1965-0179943-X.
• Кертис, AR (1970), «Процесс Рунге-Кутты восьмого порядка с одиннадцатью оценками функций на шаг», Numerische Mathematik , 16 (3): 268–277, doi :10.1007/BF02219778.
• Купер, Г. Дж.; Вернер, Дж. Х. (1972), «Некоторые явные методы Рунге–Кутты высокого порядка», Журнал SIAM по численному анализу , 9 (3): 389–405, doi :10.1137/0709037.
• Бутчер, Дж. К. (1996), «История методов Рунге-Кутты», Прикладная численная математика , 20 (3): 247–260, doi :10.1016/0168-9274(95)00108-5.

Внешние ссылки

• «Метод Рунге-Кутта», Энциклопедия математики , EMS Press , 2001 [1994]
• Метод Рунге–Кутты 4-го порядка.
• Реализация библиотеки компонентов Tracker в Matlab — реализует 32 встроенных алгоритма Рунге-Кутты в RungeKStep, 24 встроенных алгоритма Рунге-Кутты-Нистрёма в RungeKNystroemSStepи 4 общих алгоритма Рунге-Кутты-Нистрёма в RungeKNystroemGStep.