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

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

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

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

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

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

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

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

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

Теперь неявный Эйлер на той же задаче. Уравнение y_{n+1} = y_n + h·(-15 y_{n+1}) для линейного f решается явно: переносим всё с y_{n+1} влево, y_{n+1}·(1 + 15h) = y_n, откуда y_{n+1} = y_n / (1 + 15h). При h = 0.2 знаменатель 1 + 3 = 4, множитель 1/4 = 0.25 — по модулю меньше единицы. Решение затухает: 1, 0.25, 0.0625, ... — качественно ведёт себя правильно, как и точное. Никакого взрыва при любом 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 < 2/15 ≈ 0.133 только чтобы не развалиться.

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

Нам повезло: f = -15y линейна, и уравнение y_{n+1} = y_n + h·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·df/dy| < 1) — а это ровно тот режим, ради ухода из которого мы и брали неявный метод. Поэтому в серьёзных решателях жёстких задач (тот же scipy.integrate.solve_ivp с методом 'BDF' или 'Radau') внутри сидит метод Ньютона: он сходится почти при любом шаге, и именно он позволяет неявной схеме реализовать своё преимущество — большой устойчивый шаг.

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

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

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

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

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

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