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