Fedorenko-RP-Vvedenie-v-vychislitelnuyu-fiziku (810773), страница 36
Текст из файла (страница 36)
Согласно второму иэ со- Показан случай, когда траектория л(г) имеет шесть нулей (точек 14(г) = 0). Это означает, что показана пятая собственная функция (первая, например, имеет только два нуля: г = 0 и г= г'), Представим точку (Я, Я') в полярных координатах: Я(г) = А(г) соз р(г), К(г) = А(г) зш Р(г). (5) Соотношение (4) есть очевидное следствие (5). Угол Р(г) — это угол точки (Я, Я').
Краевому условию /1(0) = 0 соответствует угол Р(0) = н/2. Каждый нуль Я(г) соответстл вУет изменению Р(г) на — х. Отметим 1В(г),а'1г1) еще, что искомая функция Я(г) определена с точностью до множвтеля: здесь выбрана нормировка /т'(0) ) О, при которой угол Р(г) убывает; при нормировке Я'(0) < 0 угол ~р(г) возрастает. Я Теперь можно сформулировать краевое условие для р(г"): Р(г') = л/2 — йп, где х — номер собственного значения. Это очень удобный в практических расчетах факт: он позволяет явно задавать номер собственного значения. Получим дифференциальное уравнение для Р так, как это обычно делается в алгоритмах прогонки. Продифференцируем соотношение (4) по г: Я' згп р + М Р' соз Р— л" соз Р + Я' Р' зш ~р = О. Подставляя в это выражение Я" = (Г(г) — Х)Я, имеем Я'(1+ р') зш р+й(Р+9ь — Р(г)) соз Р=О.
(6) Соотношения (4) и (б) образуют систему линейных однородных уравнений для Я и Я . Приравнивая нулю определитель, получаем уравнение для Р: Р' + з)пг ~ +.(Х вЂ” Р(.)) ,г,р = О, р(О) = и/г. (7) ]85 $15] СПЕКТРАЛЬНАЯ ЗАДАЧА ШТУРМА-ЛИУВИЛЛЛ отношений (5) /г' = А з]п р. С учетом этого получаем А'=А(1+ р') зш р/соз р. Используя (7), вычисляем 1+ р'= 1 — з]пз р+ (]'(г) — ]ь) созз р= (1+ ~(г) — Х) созз р, В результате получаем А' А соз ~р з]п ~р (1 + ]У(у) — А). Значение А(0) произвольно, например А(0) =1, В любом случае решение А(г)соз р(г) нужно нормировать.
В заключение выясним качественную структуру траектории р(г). При г, близких к нулю, величина Х вЂ” ]У(г) > 0 и /1(г) осциллирует, р'(г) < О и р(г) монотонно убывает. На правой части интервала [О, г'[ функция к(г) ж е /~, поэтому /«(г)/Ус(г) = = с]к р(г) жсопзк Итак, р(г) сначала монотонно убывает, затем скорость убывания замедляется, и график р(г) становится почти горизонтальным.
Столь простая структура решения р(г) может создать впечатление, что уравнение для р можно интегрировать очень крупным шагом. В й5 указывалось, что при выборе шага численного интегрирования следует ориентироваться на следующее простое соображение: каждый характерный участок решения нужно «покрыть» хотя бы пятью-десятью счетными точками. И кажется, что нри интегрировании уравнения для р достаточно взять' пять- десять точек на участке осцилляций и столько же точек на экспоненциальном участке, хотя, например, на участке осцнлляций функция Ю(г) имеет, допустим,'пять полуволн.
К сожалению, это не так: численный расчет «простой» траектории р(г) требует такого же числа точек (такого же шага интегрирования), какого требует расчет «сложной» траектории /з(г), если, конечно, используются стандартные методы, например методы Рунге — Кутты. В самом деле, на участке осцилляций решение имеет (грубо, ориентировочно) вид р(г) ж я/2 — Сг. Постоянная С тем больше, чем больше номер собственного значения.
Шаг численного интегрирования должен быть таким, чтобы за это время правая часть уравнения изменилась пренебрежимо мало. Но функции з]пз р и созз р за шаг /з мало меняются лишь при условии С/з~], т.е. чем больше С (чем больше полуволн имеет решение л(г) на участке осцилляций), тем меньше должен быть шаг 6.
И наконец, несколько слов о погрешности, связанной с переносом граничного условия /1(м) = 0 в точку г'. Мы рассмотрим его опять-таки на модельном для экспоненциальной части траектории 1зб нгивлиженные методы вычислительной»ивнев 1ч.
и уравнении Я" = ~ ЦЯ. Точное решение должно иметь вид 1«(г) = С ехр ( — гт'1Ц ). Приближенное решение Я(г), обращающееся в нуль при г= г*, очевидно, есть Л = С (ехр ( — г»'1 ХГ) — ехр (-2г т'ГЦ ) ехр (л/~ХГ)). Ка участке 10, г'1 вклад в М(г) ог «лишнего» слагаемого очень мал (порядка ехр ( — г"~/Щ)). Вычисление точек комплексного спектра. Спектральная задача существенно осложняется в случае несамосопряженного дифференциального оператора, когда собственные значения могут оказаться комплексными. Рассмотрим несложную задачу — —" + а — "" + Ьх = Хх„О и г а 1, ,~~в нг (б) с краевыми условиями х(0) = х(1) =О.
Пусть коэффициенты а, Ь вЂ” какие-то несложные функции или даже постоянные комплексные числа (илн х — вектор, а, Ь вЂ” матрицы). Исследование спектра (вообще говоря, комилексного) несколько облегчается тем, что обычно задачей численного анализа является определение точек спектра, расположенных в некоторой окрестности нуля. Размеры этой окрестности, конечно, зависят от модуля а.
В области ~ Х1:» ! а ~ часто можно воспользоваться подходящими асимптотнческнми методами, т.е. изучать спектр задачи (8) как слабое возму- щенке хороню изученного спектра задачи при а = О. Эти аналитические методы достаточно эффективны при больших ~ Ц, но теряют точность и иногда просто непригодны в области малых ~ Х 1 . Однако именно точки спектра с малыми ~ Х ~ нередко представляют наибольшую прикладную ценность, и численные методы в этой области удачно дополншот аналитические исследования. Таким образом, речь идет о приближенном вычислении относительно небольшого числа точек спектра.
И здесь в зависимости ат ситуации может быть использован метод «пристрелкн», когда интегрируется задача Коши с начальными данными х(0, Х) = О, х(0, Х) = 1 и этим алгоритмом численного интегрирования определяется функция комплексного переменного Ф(Х) щ х(1, Х). В другой ситуации может оказаться целесообразным применение метода прогонки, когда, например, решение ищется в форме х(г) = а(г) х(г), а для «прогоночного коэффициента» а(г) обычным способом получается уравнение (типа уравнения Риккзти), содержащее параметр Х. Интегрируя зту задачу (конечно, численно), вычисляем а(г, Х) н определяем Ф(Х) щ а(1, А). Если на правом конце за- 1 ий спвкт»Альнля злдАч» штл'мА-лнтвнлля дано более общее краевое условие (например, х(1) + рх(1) = 0), оно вместе с прогоночным соотношением х(1) = а(1) х(1) образует систему линейных уравнений относительно х(1) и х(1), а Ф(А) определяется как детерминант этой системы.
Так нли иначе мы получаем функцию комплексного переменного, значения которой вычисляются процедурой интегрирования задачи Коши (в комплексных числах). Точки спектра исходной задачи суть нули этой функции. При определенных условиях, которые хорошо изучены в теории обыкновенных дифференциальных уравнений, Ф(й) — аналитическая функция (это следствие простой формы зависимости уравнения от А), не имеющая полюсов прн ограниченных А.
Поэтому для подсчета числа ее нулей в некоторой области С следует вычислить известный интеграл по контуру нли, проще говоря, вращение векторного поля (Ве Ф(А), 1ш Ф(Х)) при обходе контура дС. Конечно, это — громоздкая операция: надо «покрыть» контур сеткой точек (Х,.) и в каждой точке Х, вычислить Ф(Х,.), т.е. проинтегрировать задачу Коши. Несколько облегчает работу то, что сетка (»,) не должна быть особенно густой, Точнее, дело обстоит так. Расчет начинается с достаточно широкой области С, имеющей (для простоты и определенности) форму прямоугольника.
Вычислив вращение вдоль контура (оно будет равно 2лл, где л — число точек спектра в С), делят его пополам (пополам делится .та сторона прямоугольника, которая на данном этапе процесса локализации корней Ф(Л) длиннее). Вычислив вращение вдоль контура одного из полученных меньших прямоугольников, определяют число точек спектра в двух частях исходной области, н т.д. Когда прямоугольник велик (и его контур достаточно длинен), шаг сетки на контуре может быть взят достаточно большим.
Например, в некоторых расчетах, проводившихся по этой схеме, было принято считать «нормальной» ситуацию, в которой прн переходе от», к Х,+, значения ага Ф(А) изменялись в пределах интервала [я/6, и/З~. Если изменение было меньшим, шаг увеличивался, если ббльшим, — происходил возврат в точку Хп шаг Л изменения Х уменьшался (напрнмер, Л:= Л/2) н делался переход в точку Х,. », = Х,. + Л на контуре.
Таким образом, сетка Я,.) не задавалась заранее, а «генерировалась» простым алгоритмом с адаптацией (с регулированием шага в зависимости от градиента функции агк Ф(Х)). Конечно, такая тактика сопряжена с некоторым риском: прн вычислении Ф(Х,) н Ф(Х,»>) приращение аргументов определено с точностью до 2кп (А — любое целое), причем число А вычислитель назначает сам. Поясним сказанное подробнее. Вычислив комплексное число х + /у, обычно обращаются к подпрограмме, входящей в стандартную ьвз ПРИБЛИ.'КЕИНЫЕ МЕТОДЫ ВЫЧИСЛИТЕЛЬНОЙ ФИЗИКИ !ч. и библиотеку и имеющей на языке РОВТВАХ нмя АТАХ2(х, у). Результатом является главное значение агс1й(у/х); к нему из каких-то дополнительных соображений нужно добавить 2АЛ.
В рассматриваемом случае, используя предположение о том, что шаг А = Х~+, — Х~ «достаточно мал», число й выбирают таким, чтобы изменение агйФ(Х,,) по сравнению с агяФ(Х,.) было минимальным. Здесь, конечно, есть риск ошибиться на 2йи. Можно получить достаточно надежное подтверждение правильности этого решения (нли обнаружить его ошибочность), «пройдя» участок ( Х,, Аз+,) с существенно меньшим шагом, но это слишком «дорого» н не делается без достаточных к тому оснований. Вероятность ошибочности такого способа вычисления вращения векторного поля существенным образом зависит от расстояния контура до какой-то точки спектра.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
