Введение в вычислительную физику - Федоренко Р.П.
ISBN 5-7417-0002-0
Скачать (прямая ссылка):
Вторая область — область устойчивости, где Itf1I < I, Itf2I < 1. Желательно, чтобы эта область покрывала значительную часть по-луплоскости Re % < 0. Легко проверить, что | tft | « | tf2| » 1/V211| при ЦІ »1. Можно найти такое значение R, что Itf1I < 1, Itf2I < 1 при ЦІ > R. Представляет интерес граница области устойчивости — линия шах {| tfj(|) I, tf2| И) = 1. Считается необходимым, чтобы зона устойчивости содержала какую-то достаточно широкую окрестность линии Im I = 0, Re I < 0. В частности, если область устойчивости есть полуплоскость Re I < 0, схему называют Л-устойчивой. Доказана теорема о том, что Л-устойчивыми могут быть только неявные схемы не выше второю порядка аппроксимации. Схему называют А(а)-устойчивой, если область устой-
§ 17] ЖЕСТКИЕ СИСТЕМЫ ОДУ 223
чивости содержит конус I Im И ^ sin а |Re || (Re | < 0). Схему называют L-устойчивой, если в области Re | < —а2 (а =*= 0) решения разностного уравнения убывают, как qn, ще q < 1 и не зависит от %.
Безытерационные схемы типа схемы Розенброка. В последние годы была предложена некоторая общая конструкция схем интегрирования жестких систем, в которых система нелинейных уравнений не решается. Рассмотрим пример подобной схемы. Стандартный шаг интегрирования состоит из следующих операций. Пусть имеется точка хп. Вычисляется матрица A = fx(xn), и хп+j находится из уравнения
(Е - ах А - Px2A2) Хп+\Х" = /(хп + ух /(*„)). (20)
Таким образом, шаг стоит двух вычислений /, вычисления А и решения системы линейных уравнений. Параметры а, (3, 7 подбираются так, чтобы обеспечить возможно более высокий порядок аппроксимации и необходимую устойчивость.
Проиллюстрируем характерную технику подбора параметров. Разложим решение в ряд Тейлора в точке tn:
x(tn + X) = x(tn) + Xfn + xJ (fxf)n + х~ (fxxff + fJJ)n. (21)
Здесь /„ = x(tn), (fxf)n = x(tn), (fxxff + fjxf)n = x(tn). Остальные члены ряда опущены. Для схемы можно написать аналогичное разложение. Из (20) следует, что
*„+i = хп +(Е- ах А - Px2^2)-1 f(xn + ух f(xn)) = хп + xfn +
+ х2(а + 7)Af + т3(р + а2 + ay)A2f + \у2хг (fxxff)n. (22)
Распорядимся параметрами так, чтобы все выписанные члены (21) и (22) совпали. Тем самым будет обеспечен третий порядок аппроксимации. В результате мы получаем систему уравнений
а + 7=1/2, §+а(а +7) = 1/6, 72 = 1/3,
которая легко решается. Приведем числовые значения: а = 1.077, P = —0.372, 7 = —0.577. (Решение с 7 = 0.577 неинтересно.)
Перейдем к анализу устойчивости, используя схему (20) для уравнения х = Xx. После несложных преобразований имеем
*„+i = \ = хХ,
224
ЧИСЛЕННЫЕ МЕТОДЫ МАТЕМАТИЧЕСКОЙ ФИЗИКИ
[Ч. II
где
q(%) = (1 - 0.0771 - 0.205%2)/(1 - 1,077% + 0.37|2).
Простой анализ показывает, что | д(%) \ < 1 при Im | = 0, Re Sj < 0. Легко проверить, что при IIJ >10 также q(%) < 1. Более точные сведения о границе области устойчивости можно получить численно.
Отметим, что схема «устойчива» и в большей части правой полуплоскости. В этом случае качественные поведения траекторий дифференциального и разностного уравнений принципиально различны. Это не относится, разумеется, к области точности.
Регулярные жесткие системы. Изучение жестких систем обнаружило их большое сходство с сингулярно-возмущенными. Сложилось впечатление, что жесткую систему можно получить из сингулярно-возмущенной, сделав гладкую замену переменных. При такой замене теряется четкое разделение переменных на быстрые и медленные (х и у в (12)), маскируется то многообразие Г, около устойчивых ветвей которого происходит медленное движение фазовой точки (из системы (12) мы сразу получаем для Г уравнение /(•*> У) = 0). Отсутствие асимптотической теории для общих жестких систем, аналогичной теории сингулярно-возмущенных систем, затрудняет разработку и оценку численных методов. Теория объясняет, какой должна быть траектория и чего мы вправе требовать от численного метода.
Опишем возможный вариант такой асимптотической теории. Основным ее объектом является многообразие Г, определяемое уравнениями
(/(*), Ф*(х))=0, /=1,2,...,7. (23)
Здесь Ф*(х) — собственные векторы матрицы /*(*), соответствующие точкам жесткого спектра. Хотя это определение не конструктивно (так как оно оперирует с трудно вычисляемыми объектами), тем не менее в теоретическом анализе его использовать удалось.
Было показано, что в малой окрестности Г движение происходит так же, как и в окрестности поверхности f(x, у) = 0 для системы (12): все траектории очень быстро (со скоростью O(L)) входят в 0(L~2)-окрестность Г и движутся в ней со скоростью 0(1). Удалось построить «обратную» замену переменных, сводящих общую жесткую систему к сингулярно-возмущенной, разложить X на быструю и медленную компоненты и выписать уравнения их эволюции.
Имея достаточно ясное представление о том, как устроена траектория, оказалось возможным дать достаточно аккуратное обоснование некоторых разностных схем и получить оценку погрешности приближенного решения в терминах двух малых параметров: т и
ЖЕСТКИЕ СИСТЕМЫ ОДУ
