§ 1.8. Жесткие системы

Здесь описывается важный и часто встречающийся в приложениях класс так называемых жестких систем дифференциальных уравнений. Изучаются особенности численного решения жестких систем.

1.8.1. Пример. Рассмотрим систему дифференциальных уравнений

xў1= l1x1,   xў2= l2x2,

(1)

в которой l₁ = -1, а l₂ = -10⁶. Общее решение этой системы, как известно, имеет вид

x1(t) = c1e-t,   x2(t) = c2e-106t.

Допустим, нам нужно найти решение (1), удовлетворяющее начальным условиям

x1(0) = 1,   x2(0) = 1.

Воспользуемся, например, явным методом Эйлера. Поскольку его область устойчивости S (см. задачу 1.7.6) есть единичный круг в C с центром в точке -1, для устойчивого счета нужно, чтобы выполнялись неравенства

| l1t + 1| < 1,   | l2t + 1| < 1,

откуда следует, что t должно быть меньше 2·10^-5. С другой стороны, при t > 5·10^-5 компонента x₂(t) решения отличается от нулевой функции менее, чем на e = 10^-20 и, по-существу, второе уравнение после трех шагов длины t = 2·10^-5 в решении не нуждается. Первая же компонента системы (1) может интегрироваться с шагом t, выбранным из требования | l₁t + 1| < 1, т. е. с любым шагом t < 2. Таким образом, за выбор шага при интегрировании системы (1) отвечает та компонента, поведение которой за исключением малых интервалов на оси t несущественно.

1.8.2. Жесткие системы.

Описанная выше ситуация возникает из-за большого разброса собственных значений матрицы системы (1): l₂/l₁ = 10⁶. Компонента с бóльшим (по модулю) собственным значением вынуждает выбирать мелкий шаг и, одновременно, быстро перестает влиять на решение. Класс дифференциальных уравнений с таким поведением выделяется в теории численных методов понятием жестких уравнений. Точнее, система линейных автономных дифференциальных уравнений

xў = Ax

(2)

называется жесткой, если, во-первых, все собственные значения l_i (i = 1, ..., m) матрицы A имеют отрицательную вещественную часть (т. е. система (2) экспоненциально устойчива) и, во-вторых,

S =  min 1ЈiЈm Re li max 1ЈiЈm Re li  >> 1.

Число S при этом называют коэффициентом жесткости системы (2). Значок >> ("значительно превосходит") на практике обычно означает, что S > 100, хотя часто встречаются задачи (например, в теории электрических цепей, в химической кинетике и т. п.) с коэффициентом жесткости »10⁶ и более.

В случае общей системы дифференциальных уравнений (E) жесткие системы можно выделять по-разному. Наиболее простым определением является, по-видимому, следующее. Пусть j — решение уравнения (E) на отрезке [t₀, t₀ + T] и A(t) = ¶f(t, x)/¶x|_x=j(t). Говорят, что система (E) жесткая вдоль решения j на отрезке [t₀, t₀+T], если при всех s О [t₀, t₀ + T] (автономная) система

xў = A(s)x

жесткая.

Подчеркнем, что понятие жесткости относится к дифференциальному уравнению, а не к разностной схеме.

1.8.3. Еще один пример.

Попытаемся решить систему (1) неявным методом Эйлера. Как мы знаем (см. п. 1.7.7), его область устойчивости представляет собой внешность единичного круга в C с центром в точке +1: в частности, она содержит всю левую открытую полуплоскость {l О C: Re l < 0}. Поскольку l₁ и l₂ в (1) отрицательны, l₁t, l₂t О S при любом t > 0. Таким образом, при использовании неявного метода Эйлера для решения устойчивой системы нет ограничений на величину шага, вызванных требованиями устойчивости метода.

1.8.4. A-устойчивые методы.

Как упоминалось выше, если область абсолютной устойчивости метода содержит левую открытую полуплоскость {l О C: Re l < 0}, то для любой жесткой системы (2) метод абсолютно устойчив. В самом деле, поскольку для жестких систем Re l_i < 0, произведение tl_i лежит в S при любых t > 0. Методы, область абсолютной устойчивости которых содержит левую открытую полуплоскость, называются A-устойчивыми. Для A-устойчивых методов шаг t может выбираться только из соображений аппроксимации.

Однако, как оказалось, A-устойчивых методов не так уж много. А именно, доказано, что явных линейных многошаговых A-устойчивых методов (как, впрочем, и явных A-устойчивых методов Рунге — Кутты) просто нет. Порядок же аппроксимации неявного линейного A-устойчивого многошагового метода не может превышать двух.

Задача 1.8.1. Покажите, что метод трапеций

xi - xi-1 =   t 2 [f(ti, xi) + f(ti-1, xi-1)]

является A-устойчивым методом второго порядка аппроксимации.

1.8.5. A(a)-устойчивые методы.

Из вышеизложенного следует, что требование A-устойчивости резко сужает класс доступных схем. Естественно попытаться снизить требования к области устойчивости метода. Один из возможных вариантов таков. Заметим, что если спектр s(A) матрицы A лежит в левой открытой полуплоскости, в силу своей конечности, он лежит в некотором клине W_a = {l О C: |arg l - p| < a}, где 0 Ј a < p/2 (рис. 10).

Рис. 10.

Поэтому, если область абсолютной устойчивости метода содержит клин W_a (такие методы называют A(a)-устойчивыми), то его применение к системе (2) с расположенным в W_a спектром не будет иметь ограничений на шаг из-за требований устойчивости.

Здесь доказано, что не бывает явных даже A(0)-устойчивых методов, но зато для любого 0 < a < p/2 существует неявный линейный многошаговый A(a)-устойчивый метод четвертого порядка аппроксимации.

Задача 1.8.2. Покажите, что неявный линейный четырехшаговый метод

xi –  48 25 xi–1 +  36 25 xi–2 –  16 25 xi–3 +  3 25 xi–4 = t 12 25 f(ti, xi)

(3)

является методом четвертого порядка аппроксимации

Задача 1.8.3. Докажите, что метод (3) A(a)-устойчив при некотором a > 0.

Отметим, что если заранее известно, что решаемая (устойчивая) система дифференциальных уравнений имеет только вещественные (отрицательные в силу устойчивости) собственные значения (например, если матрица линеаризованной системы симметрична), то для устойчивого счета достаточно уже A(0)-устойчивости метода.

Подчеркнем еще раз, что поскольку, как уже отмечалось, не существует явных даже A(0)-устойчивых методов, для решения жестких задач могут применяться только неявные методы.

В заключение параграфа опишем еще один класс неявных линейных многошаговых методов, активно применяющихся при решении жестких задач.

1.8.6. Методы, основанные на дифференцировании.

В этих методах коэффициенты a_s, b_s линейного многошагового метода подбираются из следующих соображений (ср. c п. 1.2.12). Построим интерполяционный полином U(t) по узлам t_i-p, ..., t_i и значениям x_i-p, ... , x_i и потребуем, чтобы он удовлетворял уравнению (E) в точке t_i-r сетки, т. е. чтобы выполнялось равенство

Uў(ti-r) = f(ti-r, xi-r).

(4)

Например, если p = 1 и r = 1, то

U(t) =  t - ti-1 ti - ti-1 xi +  t - ti ti-1 - ti xi-1

и после несложных преобразований равенство (4) превращается в явный метод Эйлера.

Если r = 1, то методы вида (4) явные, если же r = 0, то неявные. Эти методы при r = 0 называют еще формулами дифференцирования назад, поскольку в левой части равенства (4) стоит аппроксимация производной решения в точке t_i-r = t_i по узлам t_i-p, ..., t_i.

Задача 1.8.4. Выведите формулы дифференцирования назад для p = 1, 2, 3.

Задача 1.8.5. Исследуйте область абсолютной устойчивости построенных в предыдущей задаче методов.

При p = 1, 2 формулы дифференцирования назад являются A-устойчивыми, при p = 3, ..., 6 — A(a)-устойчивыми при некотором a > 0, при p > 6 они не являются даже A(0)-устойчивыми.

File based on translation from T_EX by T_TH, version 3.05.
Created 29 May 2002, 10: 24.