Федоренко Р.П. Введение в вычислительную физику (1185915), страница 53
Текст из файла (страница 53)
Устойчивость же определяется уравнением в вариациях, которое для (22) имеет вид (Ь р — малое возмущение траектории <р) Ь р = з1п 2Р (а — 1/р) Ьр. (25) Устойчивость траектории зависит от знака и модуля а з1п 2у. В принципе возможна ситуация, когда е1п 2о = 1 и траектория сильно неустойчива: решение (25) ведет себя, как е1ч~', Однако из уравнения (22) видно, что р быстро уходит от такого значения и «неустойчивый» участок на траектории не может быть длительным. Ограничимся этими простыми соображениями, которые, видимо, можно превратить в достаточно строгий анализ.
Периодическая прогонка. Опишем полезный в приложениях алгоритм решения специальной системы уравнений высокого порядка, возникающей при решении краевой задачи для уравнения Штурма — Лиувилля (10.1) с периодическими краевыми условиями х(0) = х(Т), х(0) = х(Т). Разностное уравнение (10.2) можно считать определенным при всех. значениях л = О, 1, ..., Ж, если реализовать условие периодичности, отождествив выходящие за пределы сетки значения с сеточными: хи+, — — хо, х, = хи. Итак, мы приходим к системе уравнений, аналогичной (10.3): аохи ~охо + сох~ /о а„хл, — Ь хи + сихо = /„„ где л = 1, 2, „Ф вЂ” 1. Матрица системы отличается от знакомой нам трекдиагональной наличием двух ненулевых элементов: на последнем месте первой строки и на первом месте последней.
Для экономного (требующего 0(Ф) операций) решения такой системы А. А. Абрамо- $!81 ЖЕСТКИЕ ЛИНЕЬНЫЕ КРАЕВЫЕ ЗАДАЧИ вым построен алгоритм, обобщающий классическую прогонку. Он часто используется в практической вычислительной работе. Решение ищем в форме «прогоночного соотношения» х„,=Р„х„+Д„+Л„х . Очевидно, первое уравнение системы можно записать в этой форме, и мы получаем явные выражения для стартовых значений прогоночных коэффициентов: Р~ = соУЬо «2» = УВУЬо Я~ = ао~ЬВ Теперь стандартная процедура позволяет получить рекуррентную формулу.
Пусть Р„, Д„„Я„известны. Исключая из уравнения с индексом и значение х„ ,, имеем уравнение а„(Р„х + Д„ + Я„х, ) — Ь„х„ + с,х„ „ = У„, связывающее х„, х„+,, х„+и Этому уравнению можно придать стандартную прогоночную форму, разрешив его относительно х„. В результате мы получаем искомые соотношения: А=Ь вЂ” а Р, с ч е' с„ Р Я с'Р! А' ч+! А ' ~ч >! А Зту операцию можно продолжить вплоть до значения л = Ф вЂ” 1. Прогоночное соотношение х„, = Р х + (~и + ~их после подстановки в дс-е уравнение системы даст соотношение, связывающее хи с хо.
пРидадим емУ фоРмУ х„= аихо+ фл и бУдем искать решение в виде х„= п„хо+ ~„. Новые прогоночные коэффициенты а„, ~„(и = Ф, Л/ — 1, ..., 1) находим по рекуррентным формулам справа-налево, имея нх стартовые значения аи, р, . Для этого из прогоночного соотношения .т„, = Р„х„+ Ц„+ 11„хи, считая, что значения п„и Д„известны, ис- КЛЮЧИМ Х, ХАК Х„~ = Ро(ОИХО+ Р„) + 0„+ Д„(аИХО+ Ь). Приводя подобные члены, получаем рекуррентные соотношения и„, =Р„а„+Я„пи, Р„, = 0„+Р„(1„+ВЯЛ, Последнее такое соотношение имеет внд .Хо = аохо+ ро т.е. хо = роУ(1 "о).
Остальные значения х„найдем по формуле х„= а„хо + р„. Ч' пгиелиженньш методы вычислительной оизики [ч. и Пятиточечная прогонка. Опишем алгоритм решения системы уравнений с пятидиагональной матрнцей. Такие системы возникают при численном решении разностных уравнений, аппроксимирующих краевую задачу для уравнения четвертого порядка: 44х агх — +р — +д(х)=/, ОнгнТ, а4 айаг х(0) =Ао, х(0) =Во, х(Т) = А„х(Т) =В,. Ограничимся этой простейшей задачей. Вводя сетку, аппрокснмируем уравнение разностным: —, (х, — 4х„, + бх, — 4х,.„+ х„+г) + Ь4 + — р„(х„, — 2х„+ х„„) + д„х„= У„, Ьг л=2,З,...,ж — 2, Ь=Т(М. Краевые условия аппраксимируем самым простым способом: хо Ао х, — хо лВо хл = А„хл — хи-1 = ЬВс Придадим системе уравнений стандартную пятидиагональную форму: о о ~о~~ + охг = /о — Ь|хо + с~х~ — а1хг + с~хо = /1 а„х„г — Ь„х„, + с„х„— о(„х,.„, + е„х„г = ~„, ал,хл з — Ьл,хл г+ сл 1хл ~ — Ыл 1хл — — Ул алхл — Ьлхл, + с„хи = ~„, (н = 2, 3, ..., Ф вЂ” 2).
Формулы для коэффициентов системы оче- видны. Прогоночное соотношение имеет вид После несложных преобразований первые два уравнения (левые краевые условия) дают стартовые значения прогоночных коэффициентов (Р, Я, Д), и (Р, Тс, Д) . Стандартный вывод рекурреитных соотношений для прогоночных коэффициентов проводится в предположении, что в процессе прямой прогонки (слева-направо) коэффициенты (Р, Я, Д)„, и (Р, Я, Д)„г (и все предшествующие) уже найдены. С их помощью из стандартного и-го уравнения можно ис-, ключить х„г и х„, и получить связь между х„, х„+ы х„г, которая разрешается относительно х„. осгеднннне ьыстгых вг»щеннй 26! Несложные преобразования дают рекуррентные формулы: А = с' — Ь'Р, н «л' л У„'+ Ь.'д„ н+! А ' (««+! А а„+Ь„'я„ л+! А Эта операция продолжается стандартно до значения и = !!! — 2, т.е.
последнее прогоночное соотношение имеет вид хи г= Рн !хн !+Ян !хн+Дн Вместе с двумя последними уравнениями (правыми краевыми условиями) оно дает нам три линейных уравнения с тремя неизвестными хи г, хи „хл. Решив эту систему, процессом «обратной прогонки» мы вычислим все х„последовательно справа-налево. Предоставим читателю в качестве полезного упражнения внести необходимые изменения в том случае, когда краевые условия заданы с использованием вторых и третьих производных. Несколько больших изменений требует алгоритм в том случае, когда на одном конце задано одно краевое условие, на другом — три.
5 19. Осреднение быстрых вращений Рассмотрим важный в приложениях метод интегрирования специального класса обыкновенных дифференциальных уравнений. Приложения его столь разнообразны, что имеет смысл начать с абстрактной постановки задачи. Пусть имеется система уравнений г =,г(г), г(0) = 47 ! > 0 (1) б) г(0 чо) чо1 '!' чо в) г(7(«о) но) !го ~но (2) описывающая некоторое физическое явление «в главном» (факторамн, мало влияющими на эволюцию системы, пренебрегаем).
Известно общее решение — функция г(г, о ) (точнее, вектор-функция, но размерность г в дальнейшее явно входить не будет). И наконец (это существенное предположение, выделяющее узкий, но важный класс задач), пусть все траектории (1) периодичны с периодом Т(ц ), своим на каждой траектории. Итак, нам известна функция г(1, Чо), удовлетворяющая соотношениям а) г!(т, Со) = Лг(г, Оо)), 1т' г, О„ пгиелижеииые методы вычислительной»и»яки 1ч. н 262 Если бы начальные данные были заданы в момент га, общим решением (в силу независимости / ат г) была бы функция е(/ — га, са).
Систему (1) будем называть «невозмущенной», ее решение — «не- возмущенной траекторией». Пусть более полное описание явления, учитывающее влияние малых сил, приводит к системе, именуемой в дальнейшем «возмущенной»: (3) х = /(х) + еР(х), х(0) = дш е « 1. Нас интересует это более адекватное действительности описание явления. Здесь е — малый параметр, функции /, Р, Т будем считать «гладкими», т.е. они сами и их используемые в выкладках производные суть величины 0(1) (без этой оговорки предположение е «1 не имело бы смысла).
Система (3) не имеет явного решения и возникает вопрос; нельзя ли узнать что-либо о траектории (3); испалъзуя ее близосп к «интегрируемой» системе (1)? Если нас интересует ограниченный отрезок времени (например, три-пять периодов), ответ очевиден и ничего интересного в задаче нет. Из общего курса дифференциальных уравнений известна теорема о непрерывности решения задачи Коши по правой части, т.е. х(г, %>) = е(г, Се) + 0(е). (Для этого периодичность г не нужна.) Но что произойдет за «болъшой» интервал времени? Как «накопятся» последствия малого возмущения е Р(х) за время порядка 1/е? Здесь очевидного ответа нет.
Изложенная ниже достаточно сложная теория позволяет производить соответствующие расчеты. Речь идет о теории малых возмущений на больших временах. В задаче имеется малый параметр е и большой параметрг — время процесса О(1/е). Именно это последнее обстоятельство определяет нетривиалъный характер проблемы, решение которой удается продвинуть за счет использования важного свойства невозмущенной системы (1)— периодичности всех ее траекторий. Что касается «согласаванности» параметров е и г ю 0( 1/е), то она связана не с существом задачи, а просто с тем, что удается построить аппарат решения задачи (3), работакиций эффективно именно на временах О(1/е). В частных случаях удается распространить его действие иа времена 0(1/ез), а иногда и на весь интервал 10, м ]. В некоторых ситуациях удается построить метод, работающий на временах О(1/~е ), и это тоже представляет интерес.
Содержательные примеры. Рассмотрим пример, исторически положивший начало развитию и применению метода осреднения, Движение планет Солнечной системы достаточно точно описывается системой уравнений вида «х,. -3?'- — — /,(х~) + ~ ег? Рц(ха х ), 1= 1, 2, ..., /. (4) l 263 асгзлннннн выстгых нгхщзний 9 19) Здесь 1 — число планет, ! — номер планеты, х, — шестимерный вектор, описывающий состояние планеты-точки в фазовом пространстве, У; — сила притяжения Солнца, действующая на 1-ю планету, зы РП(хн х2) — сила взаимного тяготения 1-й и )-й планет, е,, — соответствуюший малый параметр.
Невозмущенная система х,=/,.(з,.), ~ 1, 2, ..., 1, (5) имеет известное решение — движение по эллипсам. Каждая планета имеет свой период Т, и, строго говоря, то, что будет излагаться ниже, неприменимо к данному примеру. Хорошую и эффективную теорию удается построить для одночастотной задачи, когда все компоненты невозмущеиной траектории имеют общий период. Обобщение этой теории на многочастотный случай (а именно таким является Солнечная система) связано с принципиальными и до настоящего времени еще не преодоленными трудностями. Тем не менее именно для расчета движения планет впервые без строгого обоснования («эврисгически») стали использоваться методы осреднения, которые берут начало в трудах классиков небесной механики, в частности Гаусса. Второй пример задачи (3) — расчет движения искусственного спутника Земли.
В этом случае 2 — сила притяжения Земли, «Р— малые силы, связанные с нестрогой сферичностью Земли, с сопротивлением крайне разреженной на высоте орбиты спутника атмосферы, с притяжением Луны и т.п. Наконец, третий пример — дрейф электрона в «скрещенных» магнитном и слабом электрическом полях. Может показаться, что для современного специалиста, вооруженного мощными ЭВМ, нижеследующее особого значения не имеет. В конце концов это обычная задача Коши, с которой «все ясно», существуют хорошие стандартные программы и можно «пробить» задачу мощью современных компьютеров. Однако речь идет об интегрировании задачи Коши на очень большом интервале времени, и здесь остро стоит вопрсх об оценке накопления вычислительных погрешностей.