Лекция 10: Численные методы решения жестких систем обыкновенных дифференциальных уравнений
9.3. Решение линейных ЖС ОДУ и вычисление матричной экспоненты
Рассмотрим один из наиболее простых численных методов решения жестких линейных систем ОДУ [9.10], основанный на представлении решения в явном виде
для линейных систем жестких ОДУ вида
Его численная реализация связана с вычислением матричной экспоненты. О свойствах матричной экспоненты (и в целом функций от матриц) можно прочитать в [9.11]. Использование для этого разложения в ряд Тейлора
представляется непригодным, так как простые оценки показывают, что по вычислительным затратам этот алгоритм сопоставим с методом Эйлера с шагом :
Следовательно, для того, чтобы k - й член разложения ряда Тейлора был, по крайней мере, порядка O(1), необходимо выполнение условия
что соответствует условию численного интегрирования с шагом Количество членов ряда при этом что неприемлемо для решения.
Представим экспоненту в следующем эквивалентном виде:
При этом значение параметра p выбирают таким, чтобы
и можно было использовать ряд Тейлора с небольшим количеством членов. Действительно, в этом случае
что вполне приемлемо при соответствующем выборе параметра p.
В этом алгоритме сначала вычисляют матрицу
затем путем последовательных перемножений. При таком способе вычисления матрицы также имеется опасность, связанная с тем, что при некоторых p, слагаемые, соответствующие мягкой части спектра, окажутся много меньше единицы (и слагаемых, соответствующих жесткой части). В этом случае предпочтение отдается представлению
а последовательность вычислений имеет вид
При этом на мягкой части спектра, поскольку имеем
на жесткой части, поскольку теперь p небольшие и получим оценку
что уже приемлемо для вычислений.