Введение в вычислительную физику - Федоренко Р.П.
ISBN 5-7417-0002-0
Скачать (прямая ссылка):
«Тригонометрическая прогонка» основана на введении функции ip (г), связанной с R(r) соотношением
R(r) sin (p(r) — R'(r) cos ip(r) = 0. (4)
Поясним его происхождение, а заодно выясним «геометрический» смысл величины ip(r) (он будет полезен). На рис. 19 (качественно) изображена фазова* плоскость (R, R') и траектория {R(r), R1 (/•)}.
184
ПРИБЛИЖЕННЫЕ МЕТОДЫ ВЫЧИСЛИТЕЛЬНОЙ ФИЗИКИ
[Ч. II
Показан случай, когда траектория R(r) имеет шесть нулей (точек R(r) =0). Это означает, что показана пятая собственная функция
(первая, например, имеет только два нуля: г = 0 и г = г*). Представим точку (R, R') в полярных координатах:
R(r) = A(r) cos <p(r), R'(r) = A(r) sin <р(г). (5)
Соотношение (4) есть очевидное следствие (5). Угол <р(г) — это угол точки (R, R'). Краевому условию R(O) = 0 соответствует угол
<р(0) = л/2. Каждый нуль R(r) соответствует изменению ^(г) на —л. Отметим еще, что искомая функция R(r) определена с точностью до множителя: здесь выбрана нормировка R1(O) > 0, при которой угол <р(г) убывает; при нормировке /?'(0) < 0 угол ip(r) возрастает.
Теперь можно сформулировать краевое условие для <р(г*): Ч>(^*) = л/2 — &л, где к — номер собственного значения. Это очень удобный в практических расчетах факт: он позволяет явно задавать номер собственного значения. Получим дифференциальное уравнение для ip так, как это обычно делается в алгоритмах прогонки. Продифференцируем соотношение (4) по г:
R' sin <р + R<p' cos <р — R" cos <р + R'tp' sin <р = 0.
Подставляя в это выражение R" = (V(r) — X)R, имеем
i?'(l + ip') sin <р + R(tр + Я — K(r)) cos tp = 0. (6)
Соотношения (4) и (6) образуют систему линейных однородных уравнений для R и R'. Приравнивая нулю определитель, получаем уравнение для ip:
<р'+ sin2 <р 4-.(Я — K(r)) cos2 <р — О, Ip(O) = л/2. (7) ,
Численным интегрированием этой задачи Коши определена функ- / ция <р(г*, X). Ее график напоминает «лестницу» с почти вертикальными участками в окрестности значений л/2 — кп. Это надо учитывать при решении уравнения <р(г*, Я) = л/2 — кн.
Пусть найдено число Xk, при котором решение (7) удовлетворяет обоим краевым условиям. Нужно найти саму собственную функцию R(r), для чего достаточно определить амплитуду А(г). Выведем дифференциальное уравнение для А(г). Дифференцируем первое из соотношений (5): R' = A' cos <р — Atp' sin (р. Согласно второму из со-
СПЕКТРАЛЬНАЯ ЗАДАЧА ШТУРМА-ЛИУВИЛЛЯ
185
отношений (5) R' = A sin ф. С учетом этого получаем A' = А(I + <р') sin ip/cos tp.
Используя (7), вычисляем
I + <р' = I — sin2 ip + (K(V) — Я) cos2 ip = (I + V(r) — Я) cos2 (p.
В результате получаем
А' = A cos <р sin ip (I + V(г) — Я).
Значение /1(0) произвольно, например /1(0) = I. В любом случае решение ^(r)cos <р(г) нужно нормировать.
В заключение выясним качественную структуру траектории ip(r). При г, близких к нулю, величина Я — V(r) > 0 и R(r) осциллирует, <p'(r) < О и ip(r) монотонно убывает. На правой части интервала [0, г*] функция R(r) « C-n^, поэтому R'(r)/R(r) = = ctg <p(r) « const. Итак, (р(г) сначала монотонно убывает, затем скорость убывания замедляется, и график ip(r) становится почти горизонтальным.
Столь простая структура решения <р(г) может создать впечатление, что уравнение для <р можно интегрировать очень крупным шагом. В § 5 указывалось, что при выборе шага численного интегрирования следует ориентироваться на следующее простое соображение: каждый характерный участок решения нужно «покрыть» хотя бы пятью-десятью счетными точками. И кажется, что при интегрировании уравнения для ip достаточно взять' пять-десять точек на участке осцилляций и столько же точек на экспоненциальном участке, хотя, например, на участке осцилляций функция R(r) имеет, допустим,' пять полуволн. К сожалению, это не так: численный расчет «простой» траектории <р(г) требует такого же числа точек (такого же шага интегрирования), какого требует расчет «сложной» траектории R(r), если, конечно, используются стандартные методы, например методы Рунге—Кутты.
В самом деле, на участке осцилляций решение имеет (грубо, ориентировочно) вид <р(г) « л/2 — Cr. Постоянная С тем больше, чем больше номер собственного значения. Шаг численного интегрирования должен быть таким, чтобы за это время правая часть уравнения изменилась пренебрежимо мало. Ho функции sin2 <р и cos2 ц> за шаг h мало меняются лишь при условии Ch« 1, т.е. чем больше С (чем больше полуволн имеет решение R(r) на участке осцилляций), тем меньше должен быть шаг h.
И наконец, несколько слов о погрешности, связанной с переносом граничного условия /?(<») = О в точку г*. Мы рассмотрим его опять-таки на модельном для экспоненциальной части траектории
186
ПРИБЛИЖЕННЫЕ МЕТОДЫ ВЫЧИСЛИТЕЛЬНОЙ ФИЗИКИ
14. 11
уравнении R" — |Я|Л Точное решение должно иметь вид R(r) = С ехр (—г/|Я| ).
Приближенное решение R( г), обращающееся в нуль при г = г*, очевидно, есть
R=C [exp (—п/1 XI ) — exp (—2r*V| А| ) exp (/V | А-1 )].
На участке [0, г*] вклад в R(r) от «лишнего» слагаемого очень мал (порядка exp (—rV| Я| )).
