Математика / Начала анализа
Неявная трапециевидная схема Кранка — Николсон для ОДУ
Неявная трапециевидная схема усредняет производные на концах шага, включая неизвестный новый уровень. Она имеет второй порядок и A-устойчивость, но при крупном шаге затухание может сопровождаться сменой знака.
Формула
x — время t от 0 до 0,25; y от 0,5 до 1. Линия exp(−2t), точка — шаг 0,6.
Обозначения
- $y_n,y_{n+1}$
- старое и новое состояния, единица y
- $t_n,t_{n+1}$
- концы шага, единица t
- $h$
- шаг, единица t
- $f$
- правая часть ОДУ, единица y/единица t
Условия применения
Требуется начальное состояние и возможность решить алгебраическое уравнение для y_{n+1}. Для нелинейного f применяют Ньютон или другой внутренний решатель с контролем невязки. Гладкость решения нужна для второго порядка; шаг не должен пересекать неразрешённый скачок.
Ограничения
Безусловная линейная устойчивость не означает монотонность: при больших шагах на жёстком затухании возможны слабые колебания. Схема не является автоматически точной для нелинейностей и требует контроля сходимости внутреннего решения. При разрыве коэффициентов нужно согласовать узел с событием.
Подробное объяснение
Трапециевидная схема следует из приближения интеграла производной площадью трапеции: среднее значение f в начале и конце шага умножается на h. Конечная производная зависит от неизвестного y_{n+1}, поэтому формула не даёт число одним присваиванием. Для каждого шага возникает нелинейное или линейное алгебраическое уравнение.
На линейной модели y'=λy схема имеет множитель R(z)=(1+z/2)/(1-z/2), z=hλ. При отрицательном λ модуль множителя не превосходит единицы, но при большом |z| он стремится к -1. Это объясняет одновременно устойчивость и возможное чередование знака у быстро затухающей компоненты.
Как и у явного Хойна, глобальный порядок для гладкой задачи равен двум, но происхождение и устойчивость различаются. Внутренний решатель должен решить уравнение с точностью, не скрывающей дискретизационную ошибку; иначе дорогое неявное вычисление не приносит пользы.
Для линейной системы шаг превращается в решение матричной системы. Если матрица постоянна, факторизацию можно переиспользовать. Для нелинейной системы Якобиан и начальное приближение влияют на стоимость и сходимость метода Ньютона.
Название Кранка—Николсон особенно распространено для параболических уравнений после пространственной дискретизации. В статье об ОДУ корректнее говорить о неявном трапециевидном методе; связь с Кранком—Николсон — это тот же временной оператор, а не новая формула для любой задачи.
Как пользоваться формулой
- Запишите остаток G(z)=z-y_n-h[f(t_n,y_n)+f(t_{n+1},z)]/2 и задайте начальное приближение z.
- Для линейной f решите полученное равенство напрямую; для нелинейной примените Ньютон или другой проверенный решатель.
- Проверьте невязку G(z), а затем примите z как y_{n+1}; не смешивайте ошибку решателя с ошибкой сетки.
- Перенесите время на h и повторите шаг, контролируя знак, масштаб и возможные события модели.
- Сравните сетки h и h/2; отдельно проверьте устойчивость на тесте y'=λy при отрицательном λ.
Историческая справка
Трапециевидное правило интегрирования относится к классической квадратуре и задолго предшествует компьютерным решателям ОДУ. Когда его применяют к интегральной форме y(t_{n+1})-y(t_n)=∫fdt, получается неявная трапециевидная схема. В анализе временной дискретизации уравнений диффузии эта же идея получила имя Кранка—Николсон после работ Джона Кранка и Филлиса Николсона в середине XX века. Для ОДУ метод рассматривается как двухточечный линейный многошаговый или одношаговый неявный метод с вторым порядком. Позднейшая теория устойчивости, жёсткости и неявных решателей уточнила, почему он полезен и где даёт колебания. Современные библиотеки используют трапециевидные шаги вместе с контролем невязки, адаптацией и обработкой событий; одна формула без этих механизмов не является готовым промышленным решателем.
Историческая линия формулы
Неявная трапециевидная формула выросла из общей квадратуры; имя Кранка—Николсон закрепилось прежде всего за временной схемой для диффузии. Атрибуция не означает, что авторы создали всю теорию неявных ОДУ или единственную реализацию метода.
Пример
Решим y'=-2y, y(0)=1 одним шагом h=0,25. Неявная трапеция даёт y1=1+0,125[-2-2y1]. После раскрытия скобок y1=1-0,25-0,25y1, поэтому 1,25y1=0,75 и y1=0,6. Точное значение e^-0,5≈0,606530660, ошибка около 0,006530660. Множитель схемы здесь (1-0,25)/(1+0,25)=0,6. Ответ показывает, что устойчивость не равна точности: при более крупном шаге решение всё ещё затухает, но может заметно отличаться от экспоненты. При h=0,125 два шага дают примерно 0,604938, поэтому сравнение сеток отделяет дискретизацию от ошибки решателя.
Частая ошибка
Не подставляйте в правую часть старое y вместо неизвестного нового. Не объявляйте метод явным только потому, что формула короткая. Не останавливайте Ньютон по одному малому приращению без проверки невязки. Не путайте отсутствие роста с высокой точностью. Не переносите вывод о монотонности на большие шаги без проверки множителя.
Практика
Задачи с решением
Линейное затухание
Условие. Для y'=-y, y(0)=1 и h=0,5 найдите один шаг.
Решение. y1=1+0,25(-1-y1), поэтому 1,25y1=0,75.
Ответ. y1=0,6.
Постоянная скорость
Условие. Для y'=3, y(0)=2 и h=0,4 найдите y1.
Решение. Оба наклона равны 3, поэтому y1=2+0,4·3.
Ответ. y1=3,2.
Дополнительные источники
- https://fncbook.com/multistep/
- https://fncbook.com/implicit/
Связанные формулы
Математика
Неявный метод Эйлера: шаг и устойчивость затухания
Неявный метод Эйлера использует производную в конце шага. Новое значение входит в обе части равенства, поэтому каждый переход требует решения уравнения, зато схема хорошо подавляет быстро затухающие компоненты.
Математика
Метод Хойна: явное усреднение двух наклонов
В этой двухстадийной схеме сначала делают прогноз Эйлера на полный шаг, затем усредняют начальный и прогнозный конечный наклоны. Метод имеет второй порядок, но названия его вариантов различаются в литературе.
Математика
Метод Адамса — Башфорта второго порядка для ОДУ
Двухшаговый метод Адамса — Башфорта использует два уже вычисленных наклона и получает явное приближение второго порядка с одним новым вызовом правой части.