Четырёхстадийный метод Гаусса–Лежандра – это неявный метод Рунге–Кутты порядка точности p=8 . Он строится на квадратурной формуле Гаусса с четырьмя узлами и относится к числу наиболее точных симплектических A-устойчивых методов интегрирования ОДУ.

Задача

Рассматривается задача Коши:

y' = f(t,y), \:\:\: y(t_0) = y_0

(1.1)

Метод вычисляет:

  • четыре внутренние стадии,
  • после чего формирует y_{n+1} .

Узлы Гаусса

Узлы на интервале [0,1] :

\begin{aligned} c_{1,2} & = \frac{1}{2} \pm \frac{1}{2} \sqrt{\frac{3}{7}+\frac{2}{7}\sqrt{\frac{6}{5}}}, \:\:\:\: c_{3,4} & = \frac{1}{2} \pm \frac{1}{2} \sqrt{\frac{3}{7}-\frac{2}{7}\sqrt{\frac{6}{5}}} . \\ \end{aligned}

(1.2)

Таблица Батчера

Полная таблица достаточно громоздка, поскольку содержит вложенные радикалы. В компактной форме:

\renewcommand{\arraystretch}{1.5} \begin{array}{c|cccc} c_1 & a_{11} & a_{12} & a_{13} & a_{14}\\ c_2 & a_{21} & a_{22} & a_{23} & a_{24}\\ c_3 & a_{31} & a_{32} & a_{33} & a_{34}\\ c_4 & a_{41} & a_{42} & a_{43} & a_{44}\\[0.4em] \hline & b_1 & b_2 & b_3 & b_4 \end{array}\:\:\:\:\:\rightarrow\:\:\:\:\: \begin{matrix} \textbf{c} & \vline & \textbf{A} \\ \hline \, & \vline & \textbf{b}^T \end{matrix}

(1.3)

Веса квадратуры

\begin{aligned} b_{1} = b_{4} & = \frac{18-\sqrt{30}}{72}, \:\:\:\: b_{2} = b_{3} & = \frac{18+\sqrt{30}}{72} . \\ \end{aligned}

(1.4)

Формулы метода

Стадии:

\begin{aligned} Y_{i} & = y_{n} + h \sum^{\substack{3}}_{\substack{j=1}}a_{i}f(t_{n}+c_{j}h, Y_{j}), \\ \end{aligned}

(1.5)

где:

  • Y_{i} – внутренние стадии,
  • a_{ij} – коэффициенты таблицы Батчера.

Переход:

\begin{aligned} y_{n+1} & = y_{n} + h \sum^{\substack{3}}_{\substack{j=1}}b_{i}f(t_{n}+c_{i}h, Y_{i}). \\ \end{aligned}

(1.6)

Функция устойчивости

Для тестового уравнения:

\begin{aligned} y' & = \lambda y \\ \end{aligned}

(1.7)

получается рациональная функция устойчивости:

\begin{aligned} R(z) & = \frac{1+\frac{1}{2}z+\frac{3}{28}z^2+\frac{1}{84}z^3+\frac{1}{1680}z^4}{1-\frac{1}{2}z+\frac{3}{28}z^2-\frac{1}{84}z^3+\frac{1}{1680}z^4}, \:\:\:\:\;\: z = h\lambda . \\ \end{aligned}

(1.8)

Это Pade-аппроксимация:

\begin{aligned} R(z) & = P_{[4/4]}(e^z). \\ \end{aligned}

(1.9)

Динамические характеристики

Как мы говорили ранее, Динамической характеристикой \textbf{D} численного метода дискретизации ОДУ называется логарифм натуральный его функции устойчивости:

\textbf{D} = ln\:R(z) .

(1.7)

Так как, в общем случае, z – комплексная величина, то функцию \textbf{D} можно представить в виде комплексной динамической харатреристики и рассматривать вещественную:

\textbf{D}_{R} = Re\:\textbf{D},

(1.8)

и мнимую:

\textbf{D}_{I} = Im\:\textbf{D}

(1.9)

динамические характеристики.

Вещественная характеристика

На следующем рисунке представлена вещественная динамическая характеристика 4-стадийного метода Гаусса-Лежандра.

Fig.1. Вещественная динамическая характеристика метода Гаусса-Лежандра ( \textbf{D}_{R} ).

Если представить z в виде комплексного числа z = \mu + i\nu , то координатная ось, отрисованная синим цветом, будет вещественной, а ось, отрисованная зеленым – будет отражать мнимую координату точки z .

Ось, представленная красным цветом, показывает величину \textbf{D}_{R} .

Клавишей <P> (на клавиатуре) можно переключать перспективный и ортогональный режимы камеры, а колесом мыши – отдалять или приближать изображение.

Идеальная вещественная динамическая характеристика представляет собой плоскость, которой принадлежит мнимая ось Im(z) = \nu , и повернутая на 45 градусов вокруг оси Im(z) = \nu .

Фактически, отклонение вещественной динамической характеристики \textbf{D}_{R} от идеальной, представляет собой амплитудную погрешность численного метода (см. Рис. 3.1).

Мнимая характеристика

Мнимая динамическая характеристика 4-стадийного метода Гаусса-Лежандра представлена на рис.2:

Fig.2. Мнимая динамическая характеристика метода Гаусса-Лежандра ( \textbf{D}_{I} ).

Идеальной мнимой динамической характеристикой является плоскость, которой принадлежит вещественная ось Re(z) = \mu , и повернутая на 45 градусов вокруг оси Re(z) = \mu .

Фактически, отклонение мнимой динамической характеристики \textbf{D}_{I} от идеальной, представляет собой частотную погрешность численного метода (см. Рис. 3.1).

Области динамической согласованности

На рис.3 представлена область (зеленый контур), внутри которой относительные погрешности метода не превыщают одного процента. Вычислить границы данной области легко по формулам (3.25) и (3.25), представленными выше Определения 2, динамической теории погрешностей численных методов.

Fig.3. Область динамической согласованности метода Гаусса-Лежандра.

Кроме того, на рис.3 представлены области (красные контуры), гарантированно содержащие нули функции устойчивости метода, а также области (синие контуры), гарантированно содержащие полюса функции устойчивости метода.