Неявный (обратный) метод Эйлера и жёсткие задачи

Иногда явный метод Эйлера не просто неточен — он взрывается, выдавая числа, скачущие от плюс к минус бесконечности, хотя настоящее решение мирно затухает к нулю. Лечится это сменой одной точки, в которой берётся наклон: не в начале шага, а в его конце. Так рождается неявный метод Эйлера — главный инструмент для жёстких задач.

Неявный (обратный) метод Эйлера — численная схема $y_{n+1} = y_n + h\cdot f(t_{n+1}, y_{n+1})$, где наклон берётся в конце шага. Неизвестное $y_{n+1}$ входит в обе части, поэтому на каждом шаге его приходится находить, решая уравнение.

Одна буква разницы — и совсем другое поведение

Сравните две формулы. Явный Эйлер: $y_{n+1} = y_n + h\cdot f(t_n, y_n)$. Неявный Эйлер: $y_{n+1} = y_n + h\cdot f(t_{n+1}, y_{n+1})$. Отличие только в индексах внутри $f$: $n$ заменили на $n+1$. Геометрически это значит, что мы шагаем не по касательной в текущей точке, а по касательной, проведённой в той точке, куда хотим попасть. Звучит как замкнутый круг — ведь чтобы узнать наклон в конце, нужно уже знать конец. И верно: $y_{n+1}$ теперь стоит и слева, и справа (внутри $f$), это уравнение относительно неизвестного, а не готовая формула. Отсюда и слово «неявный».

Зачем терпеть такое усложнение: устойчивость

Расплата за неявность — необходимость на каждом шаге решать уравнение. Зачем оно того стоит? Ответ — устойчивость. Есть класс уравнений, на которых явные методы требуют абсурдно крошечного шага, иначе разваливаются. Их называют жёсткими (stiff). Типичный признак — в решении есть очень быстро затухающая компонента: что-то стремительно «садится» на равновесие, а дальше система живёт спокойно.

Возьмём модельную задачу $y' = -15y$, $y(0) = 1$. Точное решение — $e^{-15 t}$, оно стремительно падает к нулю: за $t = 0.2$ уже $0.05$, за $t = 0.4$ почти ноль. Решение скучное и спокойное. Но посмотрим, что с ним делает явный Эйлер при умеренном шаге $h = 0.2$.

Для линейного $f = -15y$ явный шаг — это $y_{n+1} = y_n + h\cdot(-15 y_n) = (1 + h\cdot(-15))\cdot y_n = (1 - 3)\cdot y_n = -2\cdot y_n$. Множитель равен $-2$! Каждый шаг численное решение умножается на $-2$: оно прыгает $1,\,-2,\,4,\,-8,\,16,\,\dots$ — знакопеременный взрыв, уносящийся в бесконечность, тогда как точное решение спокойно стремится к нулю. Катастрофа — и не из-за точности, а из-за неустойчивости.

Теперь неявный Эйлер на той же задаче. Уравнение $y_{n+1} = y_n + h\cdot(-15 y_{n+1})$ для линейного $f$ решается явно: переносим всё с $y_{n+1}$ влево, $y_{n+1}\cdot(1 + 15h) = y_n$, откуда $y_{n+1} = \frac{y_n}{1 + 15h}$. При $h = 0.2$ знаменатель $1 + 3 = 4$, множитель $1/4 = 0.25$ — по модулю меньше единицы. Решение затухает: $1,\,0.25,\,0.0625,\,\dots$ — качественно ведёт себя правильно, как и точное. Никакого взрыва при любом $h$: знаменатель $1 + 15h$ всегда больше единицы, неявная схема устойчива безусловно.

import math

LAM = -15.0          # жёсткий коэффициент: y' = -15 y
Y0  = 1.0
h   = 0.2
n   = 6

# Явный Эйлер: y_{n+1} = y_n + h*(-15*y_n)
def explicit(h, n):
    y, seq = Y0, [Y0]
    for _ in range(n):
        y = y + h * (LAM * y)
        seq.append(y)
    return seq

# Неявный Эйлер для линейного f: y_{n+1} = y_n / (1 - h*LAM)
def implicit(h, n):
    y, seq = Y0, [Y0]
    for _ in range(n):
        y = y / (1 - h * LAM)
        seq.append(y)
    return seq

expl = explicit(h, n)
impl = implicit(h, n)

print(f"mnozhitel yavnyy  = {1 + h * LAM:.2f}")
print(f"mnozhitel neyavnyy = {1 / (1 - h * LAM):.2f}")
print()
print(f"{'t':>4} {'yavnyy':>12} {'neyavnyy':>10} {'tochnoe':>10}")
for k in range(n + 1):
    t = k * h
    print(f"{t:4.1f} {expl[k]:12.5f} {impl[k]:10.5f} {math.exp(LAM * t):10.5f}")

Вывод:

mnozhitel yavnyy  = -2.00
mnozhitel neyavnyy = 0.25

   t       yavnyy   neyavnyy    tochnoe
 0.0      1.00000    1.00000    1.00000
 0.2     -2.00000    0.25000    0.04979
 0.4      4.00000    0.06250    0.00248
 0.6     -8.00000    0.01562    0.00012
 0.8     16.00000    0.00391    0.00001
 1.0    -32.00000    0.00098    0.00000
 1.2     64.00000    0.00024    0.00000

Контраст разительный. Явный столбец скачет с нарастающим размахом в разные стороны — полная чушь. Неявный аккуратно затухает к нулю вслед за точным решением. Неявный не идеально точен (на $t = 0.2$ он даёт 0.25 против истинных 0.05 — шаг великоват для точности), но он качественно верен: затухает, а не взрывается. А это можно исправить, чуть уменьшив шаг, тогда как явному методу пришлось бы брать $h \lt 2/15 \approx 0.133$ только чтобы не развалиться.

Как работает под капотом: что делать, когда f нелинейна

Нам повезло: $f = -15y$ линейна, и уравнение $y_{n+1} = y_n + h\cdot f(t_{n+1}, y_{n+1})$ решилось алгеброй за одну строчку. В общем случае $f$ нелинейна, и явно выразить $y_{n+1}$ не получится. Тогда на каждом шаге решают нелинейное уравнение численно. Два классических подхода:

Prostaya iteraciya (metod nepodvizhnoy tochki):
  vzyat nachalnoe priblizhenie  Y := y_n   (ili yavnyy shag Eylera)
  povtoryat:
      Y := y_n + h * f(t_{n+1}, Y)
  poka sosednie Y ne sblizyatsya
  zatem  y_{n+1} := Y

Metod Nyutona (bystree, nuzhna proizvodnaya):
  ishchem koren  g(Y) = Y - y_n - h * f(t_{n+1}, Y) = 0
  iteraciya:    Y := Y - g(Y) / g'(Y),
  gde          g'(Y) = 1 - h * df/dy (t_{n+1}, Y)
  shoditsya kvadratichno za 2-3 shaga

Простая итерация проще, но сходится только при достаточно малом $h$ (когда $|h\cdot \partial f/\partial y| \lt 1$) — а это ровно тот режим, ради ухода из которого мы и брали неявный метод. Поэтому в серьёзных решателях жёстких задач (тот же scipy.integrate.solve_ivp с методом 'BDF' или 'Radau') внутри сидит метод Ньютона: он сходится почти при любом шаге, и именно он позволяет неявной схеме реализовать своё преимущество — большой устойчивый шаг.

Частые ошибки

Знак в знаменателе. Неявная формула для $y' = \lambda\,y$ — это $y_{n+1} = \frac{y_n}{1 - h\lambda}$. С $\lambda = -15$ знаменатель равен $1 - h\cdot(-15) = 1 + 15h$. Легко перепутать знак и получить $1 - 15h$ — тогда при $h = 0.2$ вышло бы $-2$, и неявный метод «взорвался» бы на ровном месте. Аккуратно подставляйте знак $\lambda$.

«Неявный всегда точнее». Нет. Неявный Эйлер — тоже метод первого порядка, его точность не выше явного. Его преимущество — устойчивость, а не точность. На нежёстких задачах с маленьким шагом явный метод и проще, и не хуже.

Применять неявный метод везде. Решать уравнение на каждом шаге дорого. На обычных, нежёстких задачах это пустая трата времени — там хватает явных методов. Неявные схемы оправданы именно на жёстких задачах, где они позволяют брать большой шаг.

Считать жёсткость синонимом «сложности». Жёсткость — это не про нелинейность или хаос. Наша задача $y' = -15y$ линейна и проста, но жёстка: в ней есть быстрая мода, заставляющая явный метод семенить мелким шагом.

• Неявный Эйлер берёт наклон в конце шага: $y_{n+1} = y_n + h\cdot f(t_{n+1}, y_{n+1})$; неизвестное входит в обе части.
• Главное достоинство — безусловная устойчивость: на жёстких задачах он не взрывается при любом шаге.
• На $y' = -15y$, $h = 0.2$ явный множитель $-2$ (знакопеременный взрыв), неявный $0.25$ (корректное затухание).
• Для линейного $f$ уравнение шага решается алгеброй; для нелинейного — простой итерацией или методом Ньютона.
• Цена устойчивости — решение уравнения на каждом шаге, поэтому неявные методы применяют именно к жёстким задачам, а не повсюду.
• Неявный Эйлер по точности тоже первого порядка — он берёт устойчивостью, а не точностью.