Федоренко Р.П. Введение в вычислительную физику (1185915), страница 37
Текст из файла (страница 37)
Решение ищем в виде ~Р= У, (6, Р)Я(г)/г, где У, есть известная сферическая функция (1, т — целые числа). Обозначая 2 = 26 Е р(г) — 2Е у(г) + Н~+ '), г =„т т !зг пуиБлиженные методы Вычислите!!ьной Физики 1ч. и имеющее важное следствие: оно определяет знак точек спектра Х; нз него следует Х «О, В самом деле, при достаточно большом г потенциал Ф" становится очень малым и решения уравнения (2) качественно совпадают с решениями уравнения Я" + 1Я = О, т.е.
Я(г) С!е'~ "+ Сге "У ". Если Х > О, это — функция типа апъ'Х г, что, очевидно, несовместимо с нормировкой (3). (Однако эти решения ограничены, поэтому все Х > 0 образуют сплошной спектр, которым мы не интересуемся.) Если Х < О, нормировка достигается прн С, = О. Интервал (О и г < ее) можно (грубо) представить в виде двух частей. При г, близких к нулю, знак множителя $'(и) — Х определяется потенциалом у'(г) и, так как У(г) < О, может быть отрицатель- и ным. В этом случае решения урави НВ ! Л-ехР Пения Я" = (гг — Х)Я имеют коле! бательный характер. При больших г знак множителя 1г(г) — г! определяется величиной ( — г.) > О, т.й Р!У>-г<О Р!У>-1»О функция Я(г) имеет экспоненциРнс. !З альный характер типа е " (рнс.
18), причем, как мы увидим в дальнейшем, врассчетах придется иметь дело с достаточно большими значениями л~ ( Ц . Чтобы правильно оценить целесообразность метода, который будет изложен ниже, начнем с анализа трудностей, встречающихся при попытке решения задачи стандартными средствами. Используем метод «пристрелки». Значение Я(0) =О. Зададим Я'(О) произвольно, например Я'(0) = 1, и решим (приближенно) задачу Коши для уравнения (2), считая г! тем илн иным образом заданным. Таким образом, мы имеем функцию Я(г, Х). Искомые собственные значения — это те дискретные величины Х, при которых функция Я(г, Х) прн г-» ее имеет асимптотику Се "~!', т,е.
не содержит второй растущей компоненты общего решения. На практике просто назначают достаточно большое значение г' н ставят условие Я(г', Х) = О, Это есть уравнение для собственных значений, которое решается, например, методом Ньютона или просто подбором: при А! получаем Я(г', г!!) > О, прн Хз имеем Я(г', гг) < О, возьмем г!э = 0.5(Х!+ Лг), и т,д. Однако численная реализация этой процедуры наталкивается на серьезные затруднения. Дело в том, что на экспоненциальном участке решение задачи Коши неустойчиво.
Малое отклонение численного решения от точного приведет к появлению в решении ком- спвктгАльидя «АдА«А шттгмА-лиувилля З г51 1зз поненты С,е"' " (пусть даже с малым коэффициентом С,), и при больших г эта компонента станет основной частью, решения. Кроме того, погрешности численного интегрирования в этом случае также дают вклад в приближенные решения порядка Йге'" " (Ь вЂ” шаг численного интегрирования, р — порядок точности используемого метода). Трудности подобного рода преодолеваются методом прогонки.
В этом случае решение ищется в виде «прогоночного соотношения» Я(г) = а(г)А'(г). Для искомой функции а(г) легко выводится (см. 9 9) уравнение а = 1 — а (г'(г) — Х). Левое краевое условие очевьщно: а(О) = О, Условием на правом конце является условие выхода на асимптотнку Я » е «, которое можно аппроксимировать условием а(г') = — 1И вЂ” Х. Реализация такого подхода связана с двумя затруднениями. На участке, где 1г(г) — Х < 0 н решение имеет колебательный характер, обязательно есть точки г, в которых Я' = 0 и, следовательно, а = «».
На участке, где Я(г) экспоненциально убывает, тоже возникают осложнения. Разберемся в существе дела, пренебрегая величиной Р(г) по сравнению с Х. Уравнение а' = 1 — ~ Ц аз имеет два положения равновесия: а = «= 1/~/Щ . Несложный анализ поля направлений показывает, что ветвь а = 1/ьгГЦ является асимптотнчески устойчивой, а ветвь а = — 1/~/ ГЦ, наоборот, — неустойчивой. Нужно попасть именно на эту вторую ветвь. Перейдем к описанию алгоритма, справляющегося с этими трудностями, — к алгоритму тригонометрической прогонки, Перед изложением существа дела сделаем несколько замечаний о характере вычислительной проблемы.
Для приложений необходимо несколько первых собственных значений (нумерация собственных значений делается в порядке возрастания )Х ~). Такие расчеты надо делать для нескольких 1, т.е. задачу нужно решать многократно. Следовательно, алгоритм ее приближенного решения должен быть достаточно эффективным, экономичным. Оператор (2) самосопряжен, поэтому все Х„действительны. Напомним еще осцилляционную теорему: естественная нумерация собственных значений связана с числом нулей функции Я(г), «Тригонометрическая прогонка» основана на введении функции р(г), связанной с Я(г) соотношением К(г) з1п р(г) — К(г) соз р(г) = О.
Поясним его происхождение, а заодно выясним «геометрический» смысл величины р(г) (он будет полезен). На рис. 19 (качественио) изображена фазовая плоскость (Я, /1') и траектория (Я(г), А'(г)). 1ч. и пеиедиженные методы еынислительной»изики 1З4 Рис. 19 Численным интегрированием этой задачи Коши определена функция Р(г', Х). Ее график напоминает «лестницу» с почти вертикальными участками в окрестности значений л/2 — /сн. Это надо учитывать при решении уравнения ~р(г', Х) = и/2 — хл, Пусть найдено число Х, при котором решение (7) удовлетворяет обоим краевым условиям. Нужно найти саму собственную функцию .к(г), для чего достаточно определить амплитуду А(г). Выведем дифференциальное уравнение для А(г). Дифференцируем первое иа соотношений (5): Я' = А' соз Р— А Р' з1п р.
Согласно второму иэ со- Показан случай, когда траектория л(г) имеет шесть нулей (точек 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:» ! а ~ часто можно воспользоваться подходящими асимптотнческнми методами, т.е.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
