Федоренко Р.П. Введение в вычислительную физику (1185915), страница 46
Текст из файла (страница 46)
Представляет интерес граница области устойчивости — линия |пах (!д,Д)!, д !Ц) = 1. Считается необходимым, чтобы зона устойчивости содержала какую-то достаточно широкую окрестность линии 1ш» =О, йе 8 с О. В частности, если область устойчивости есть полуплоскость Ве 8 < О, схему называют А-устойчивой. Доказана теорема о том, что А-устойчивыми могут быть только неявные схемы не выше второго порядка аппроксимации. Схему называют А(п)-усгойчивой, если область устой- ЖЕСТКИЕ СИСТЕМЫ ОДУ ззз 6 171 чивости содержит конус ~1ш Ц < з(п а )ВеЦ (Ве г, < О). Схему называют ь-устойчивой, если в области Ве с < — аз (а ~ О) решения разностного уравнения убывают, как в", где д < 1 и не зависит от $.
Безытерационные схемы типа схемы Розенброна. В последние годы была предложена некоторая общая конструкция схем интегрирования жестких систем, в которых система нелинейных уравнений не решается. Рассмотрим пример подобной схемы. Стандартный шаг интегрирования состоит из следующих операций. Пусть имеется точка х„. Вычисляется матрица А =/ (х„), и х„ , находится нз уравнения (Š— атА — ртзА2) "+' " = /(х„ + тт /(х„)). (20) Таким образом, шаг стоит двух вычислений /', вычисления А и решения системы линейных уравнений. Параметры а, р, 7 подбираются так, чтобы обеспечить возможно более высокий порядок аппроксимации и необходимую устойчивость.
Проиллюстрируем характерную технику подбора параметров. Разложим решение в ряд Тейлора в точке 1„: 2 з х(1„+ т) — х(Т„) + т/„+ — (/„/)„+ — (/„„О+ /,/„/)„. (21) Здесь /„= х(1„), (/„/)„= х(/„), (/ „//+ /,/„/)„= х(г„). Остальные члены ряда опущены. Для схемы можно написать аналогичное разложение. Из (20) следует, что ' х„+, —— х„+ (Š— атА — рт2А2) ' /(х„+ тт /(х„)) = х„+ т/„+ .1- 22(а -1- 7)А/+ т~(р+ аз+ ат)А2/+ — 72т~ (/„,//)„. (22) Распорядимся параметрами так, чтобы все выписанные члены (21) и (22) совпали. Тем самым будет обеспечен третий порядок аппроксимации. В результате мы получаем систему уравнений а + 7 = 1/2, р + а(а + 7) = 1/б, уз = 1/3, которая легко решается.
Приведем числовые значения: а= 1.077, р = — 0.372, у = — 0.577. (Решение с у = 0.577 неинтересно.) Перейдем к анализу устойчивости, используя схему (20) для уравнения х= Хх. После несложных преобразований имеем х„ь, = д(~)х„, ф = Т1, гг« численные методы млтемхтической ьизики !ч. и где !гД) (1 0 077~ 0 205гг)/(! 1 077~ + 0 37!г) Простой анализ показывает, что !!!(Ц)~ < 1 при 1ш Ц= О, Ке К < О. Легко проверить, что при ! 1,! > ГО также д(г) <!. Более точные сведения о границе области устойчивости можно получить численно. Отметим, что схема «устойчива» и в большей части правой полуплоскостн.
В этом случае качественные поведения траекторий дифференциального н разностного уравнений принципиально различны. Это не относится, разумеется, к области точности. Регулярные жесткие системы. Изучение жестких систем обнаружило нх большое сходство с сннгулярно-возмущеннымн. Сложилось впечатление, что жесткую систему можно получить из сингулярно-возмущенной, сделав гладкую замену переменных. При такой замене теряется четкое разделение переменных на быстрые и медленные (х и у в (12)), маскируется то многообразие Г, около устойчивых ветвей которого происходит медленное движение фазовой точки (из системы (12) мы сразу получаем для Г уравнение 7'(х, у) = 0). Отсутствие асимптотической теории для общих жестких систем, аналогичной теории сингулярно-возмущенных систем, затрудняет разработку н оценку численных методов.
Теория объясняет, какой должна быть траектория н чего мы вправе требовать от численного метода, Опишем возможный вариант такой асимптотической теории. Основным ее объектом является многообразие Г, определяемое уравнениями (23) (Дх), Ф~(х)) = О, ! = 1, 2, „., ). Здесь Ф,"(х) — собственные векторы матрицы ~'„(х), соответствующие точкам жесткого спектра. Хотя это определение не конструктивно (так как оно оперирует с трудно вычисляемыми объектами), тем не менее в теоретическом анализе его использовать удалось.
Было показано, что в малой окрестности Г движение происходит так же, как и в окрестности поверхности у(х, у) = 0 для системы (12): все траектории очень быстро (со скоростью 0(Ь)) входят в О(г. г)-окрестность Г н движутся в ней со скоростью О(1). Удалось построить «обратную» замену переменных, сводящих общую жесткую систему к сингулярно-возмущенной, разложить х на быструю н медленную компоненты и выписать уравнения их эволюции. Имея достаточно ясное представление о том, как устроена траектория, оказалось возможным дать достаточно аккуратное обоснование некоторых разностных схем н получить оценку погрешности приближенного решения в терминах двух малых параметров: т и ггз й г71 жесткие систхмы оду 1//.т.
Это характерное обстоятельство: теория численного интегрирования жестких систем не может рассматривать предельного перехода при т- О, Поэтому теорема «аппроксимация + устойчивость (А, А(а), 1, н т.п.) «сходимостьь здесь не действует. Точнее, теорема справедлива, но бесполезна, так как нужно обосновать метод не в пределе прн т -ь О, а совсем в другой области значений т, когда малы оба упомянутых выше параметра. Теория, оперирующая только с аппроксимацией и А-устойчивостью, принципиально не полна.
Это иллюстрирует следующий пример явной А-устойчивой схемы, аппроксимирующей линейную жесткую систему х = Ахп (хк ы — х„)/т = р(х„, т) Ах„, (24) где р(х„, т) = (енх >' — 1)/(г(х„)т), а г(х) есть отношение Рэлея в точке гп г(х) = (Ах, х)/(х, х). Легко проверить аппроксимацию: прн т — О, очевидно, р = 1 + 0(т). Так же несложно проверить А-устойчивостьс если х — скаляр, Ах = — Хх, Ве й с О, то г(х) = Х, р = (е"' — 1)/(Хт) н нз (24) получаем х„«, = имх„. Можно ли на этом основании утверждать, что хотя бы для линейных жестких систем построен явный А-устойчивый алгоритм численного интегрирования? Видимо, нет.
Проанализируем вычисления по схеме (24). Если разложение точки х„в сумму по собственным векторам А содержит существенную жесгкую компоненту (при интегрировании в области пограничного слоя), то г(х„) ж — Ь н при Ьтм ! величина р ш 1//.т, т.е. формула (24) превращается в интегрирование с малым шагом т" = р(х, т)т 1//.. Следует только иметь в виду, что и время нужно интегрировать по формуле, аналогичной (24): г„,, = г„+ рт. Пройдя слой, траектория х попадает в область, где в разложении х„вклад жесткой компоненты пренебрежимо мал. В этом случае «(х„) = О(! ), р = 0(1) и делается шаг действительно с большим т. Но это сразу же приводит к росту в х„, жесткой компоненты, фактический шаг рт снова падает, и т.д. Трудно оценить эффективность такого адаптирующеюся алгоритма.
Как показали исследования В. И. Лебедева, явные схемы имеют некоторые возможности и при интегрировании жестких систем. Однако его теория основана на достаточно сложных построениях последовательности чередующихся малых и больших шагов. Она имеет самую непосредственную связь с устойчивыми последовательностями параметров в методе чебышевского ускорения итераций Гсм. Я 14). В-теория численного интегрирования.
Изложим основные понятия и результаты развиваемой в последние годы специальной теории численных методов для жестких систем (2). Класс изучаемых 8 — 1333 ЧИСЛЕННЫЕ МЕТОДЫ МАТЕМАТИЧЕСКОЙ ФИЗИКИ 1ч. и 226 систем выделяется важной количественной характеристикой правой части, называемой одиосторовнеа константой Лившица. Предполагается, что У удовлетворяет условию (25) (У(у) — У(х), у — х) ц 1зу — хяз, 1 О(1). Величина 1 считается имеющей порядок О(1), в то же время классическая константа Липшица или, что почти то же самое, Й1„Й может быть сколь угодно большой величиной: 1.ы 1.
Целью В-теории является получение таких оценок точности численного решения, которые не зависят от больших констант 1., а сформулированы в терминах только односторонней константы Липшица 1. Разумеется, эти оценки должны зависеть от гладкости искомого точного решения системы (2). В дальнейшем будем обозначать точное решение Х(1). Эта функция предполагается гладкой в том смысле, что Х(1) = О(!), Х(1) = О(! ), Х(Г) = О(! ),, (26) Другими словами, столько производных точного решения, сколько нужно при проведении тех или иных оценок, считаются величинами порядка О(1). Следовательно, речь идет об интегрировании системы вне слоев. Это предположение согласуется с предшествующим анализом.
Мы имеем дело с гладкой траекторией, со всех сторон окруженной существенно негладкими траекториями. При этом окружающие траектории содержат кратковременные участки, на которых их производные очень большие (О(Ь)), и являются гладкими вне этих тонких слоев. Одностороннее условие Липшнца ззрантнрует важное свойство множества траекторий системы. Пусть Х(1), У(1) — две такие траектории.
Оказывается, они не могут сильно расходиться с течением времени. В самом деле, оценим „'У, П У(1) — Х(1)!1з = $ (У вЂ” Х, У вЂ” Х) = = 2(у — Х, у — Х) = 2у(у) —,1(Х), у-Х) ц 21 Й у — ХЙ'. Отсюда !! у(1) — Х(1) 11 !! у(0) — Х(0) П «Т'. (27) Если бы мы использовали классическую константу Липшица, мы имели бы в экспоненте показатель 1.г, разрешающий существенно более сильно «разбегание» траекторий. При 1< 0 системы (2), (25) хорошо известны в теории дифференциальных уравнений под именем «диссипативных». В этом случае все траектории с ростом 1 сближаются, т.е.
обладают свойством «аттрактивности». 227 й Е71 ЖЕСТКИЕ СИСТЕМЫ ОДУ Общая теория В-сходимости. Анализ точности численного интегрирования жестких систем оперирует со следующими основными объектами: 1. Разностная схема, записанная, для определенности, в виде х„, = к„+ т Ф(е, к„, х„+,). (28) Здесь х„ — приближенное решение в точке Е„ = лт. Схема неявная; Ф, конечно, тем нли иным способом выражается через г', 2. Ограничение на сетку точного решения Х„ = Х(г„).
3. Погрешность аппроксимации (невязка а««,), которая получается при подстановке (Х„) в разностное уравнение: 4. Погрешность согласования (эта характеристика активно используется в западной литературе по численным методам) определяется следующим образом. Решим разностное уравнение (28), взяв в качестве х„ точное значение Х„. Получим некоторую величину е,.