Численные решения задачи на собственные значения для обыкновенных дифференциальных уравнений второго порядка с особыми точками

1.       ВВЕДЕНИЕ

Мы рассматриваем следующее радиальное уравнение Шредингера для кулоновских состояний

                                                                                               (1)

где r > 0 - радиальная переменная, l - неотрицательное целое (квантовое число орбитального углового момента), а - целое (а = 1 соответствует кулоновскому притяжению в атоме водорода), Rl (r) - кулоновская радиальная волновая функция, E -энергия. В общем случае, расстояния и энергия измеряются в единицах (эффективный боровский радиус) и  (эффективный ридберг), соответственно, где e - заряд электрона, m - эффективная масса электрона, или приведённая масса - в случае относительного движения двух притягивающихся частиц (экситон), к - статическая диэлектрическая проницаемость.

Уравнение (1) имеет две особые точки: регулярную особенность при r = 0 и иррегулярную особенность при r =

Уравнение (1) описывает кулоновские состояния атома водорода и водородоподобного атома в различных системах, например, состояния мелких донорных примесей в объёмных прямозонных полупроводниках или в полупроводниковых квантовых точках. В настоящей работе рассматривается задача на собственные значения и связанные с ней проблемы для уравнения (1). Энергии и волновые функции нескольких нижайших дискретных состояний вычисляются как в случае, когда известны точные аналитические решения задачи (см., например, [1,2]), что позволяет оценить точность расчётов, так и в задачах, решение которых возможны только с применением численных методов. Рассматривается задача о связанных состояниях водородоподобного атома в квантовой точке, о волновой функции основного состояния водородоподобных атомов в модели центрально-ячеечного потенциала нулевого радиуса [3,4], о связанных состояниях водородоподобного атома в модели кулоновского потенциала с жёсткой сердцевиной. Для решения этих проблем используются явные численно­аналитические методы, основанные на методе рекуррентных последовательностей [5, 6, 7].

Подстановка сводит (1) к системе двух уравнений первого порядка:

                                                                                                                               (2)

где 

Ясно (по определению, см. [8]), что точка r = есть регулярная особая точка для уравнения (1), а r = от - иррегулярная особенность (замена переменной z = 1/r). Методы построения решений уравнения (3), решения сингулярной задачи на собственные значения и задачи рассеяния для (3) развиты в [6,7] в общем случае, т.е. когда A(r) - это произвольная голоморфная в окрестности r = 0 n х n матричная функция, а (r) - n-мерный столбец. В данной работе рассматриваются решения уравнения, представленного в виде ОДУ второго порядка (1) (см. [5]).

2.          ФУНДАМЕНТАЛЬНАЯ СИСТЕМА РЕШЕНИЙ

2.1.      Окрестность регулярной особенности r = 0

Регулярное решение (1) имеет вид степенного ряда

                                                                                                                             (3)

Где коэффициенты a_n удовлетворяют следующему рекуррентному соотношению

                                                                                                  (4)

Ясно, что коэффициент a_l есть произвольная константа C, которая определяется нормировкой. Выражение для производной решения получается непосредственным дифференцированием (3).

Второе, иррегулярное решение (1) имеет вид алгебраической комбинации степенных рядов и логарифмической функции:

                                                                                                           (5)

Здесь  (r) - это решение (3), (4) с коэффициентом a_l, определяемым равенством

Коэффициенты b_n удовлетворяют следующему рекуррентному соотношению

                                                                               (6)

                                                                                                             (7)

В свою очередь, подстановка

                                                                                        (8)

сводит (7) к следующему уравнению для функции

                                                                                            (9)

A(r) - 2×2 матричная функция с диагональным матричным коэффициентом A₀:

                                                                          (10)

Для нахождения асимптотического разложения решений (8) при at r →∞ сделаем подстановку

z(r) = P (r)ω(r),                                                                                                                                       (11)

где матричная функция P (r) есть асимптотический ряд по обратным степеням r

причём полагаем, что P₀ - единичная матрица, те. . Подчиним матрицу P(r) такому требованию, чтобы функция ω(r) была решением уравнения

                                                                                                                             (12)

где B(r) - диагональная матричная функция: которая также есть асимптотический ряд по Матрицы P_n и B_n удовлетворяют следующим рекуррентным соотношениям

                                                       (13)

Поскольку при нашем выборе матрицы B₀ выполняется равенство

имеем:                                                                                                                      (14)

где

                                                                                     (15)

и произвол в выборе диагональных элементов матриц . Для определённости мы их полагаем равными нулю:

                                                                                                                                 (16)

При получаем из (13):

                                                                                                          (17)

Интегрирование уравнения (12) с диагональной матрицей B даёт:

                                                                                   (18)
где C_i - произвольные константы, что приводит к следующему выражению для асимптотического разложения решений уравнения (9):

                                                                                                                    (19)

Следует отметить, что равенства (18), (19) дают выражение для асимптотического разложения решений (9) при r → ∞ и в случае произвольного числа уравнений в системе (9). В рассматриваемом случае n=2, E < 0 (состояния дискретного спектра) нас будут интересовать только сходящиеся при r → ∞ решения (7), (1). Соответствующие асимптотические разложения решения и его

производной имеют вид

                                                                               (20)

где ω₂ определяется (18) с произвольной константой C₂ и (B₀)₂ = - λ.

3.         СШИВКА РЕШЕНИЙ. ЗАДАЧА НА СОБСТВЕННЫЕ ЗНАЧЕНИЯ

3.1.      Стандартный случай

Стандартная задача на собственные значения для уравнения (1) на [0, ∞) состоит в следующем: нужно найти такие значения параметра E: E < 0 (энергии), при которых решения R_i (r) ограничены при r = 0 и R_i (r) → 0 при r → ∞, нужно также найти решения, соответствующие этим собственным значениям E, и нормированные условием:

                                                                                                                                    (21)

Для решения задачи необходимо в некоторой промежуточной точке r = ȓ сшить регулярное при r = 0 решение (3) и его производную (lhs - «левое» решение), соответственно, с решением, сходящимся при r → ∞, и его производной (см. (20)) (rhs - «правое» решение). Если обозначить через Ф₀( E) столбец, составленный из значений, принимаемых lhs (3) и его производной при a_i = 1 (4) в точке r = ȓ при данном E, а - через (E) - столбец из rhs и его производной при С₂ = 1, и определить 2х2 матрицу

то для выполнения сшивки должно выполняться равенство:

                                                                                                         (22)

Таким образом, процедура вычисления некоторого собственного значения E = E₀ сводится к численному решению уравнения

det A(E) = 0                                                                                                                                        (23)

После определения собственного значения Е₀, решая (22) находим связь между компонентами вектора х(Е₀), а затем, используя (21), - значения констант (С₂Е₀) и ai(E₀), чем полностью определяем собственную функцию и её производную на [0, ∞).

Следует отметить, что сшивка решений всегда может быть проведена с заданной точностью с помощью правильного выбора точки сшивки ȓ. Формально это можно выразить следующим образом: пусть X = Х(ȓ) обозначает некоторый математический объект, который нужно вычислить (например, энергия, собственная функция), тогда

При практических вычислениях, для оценки ошибок удобно использовать вариацию точки сшивки

Для оценки точности вычислений найдены энергии и волновые функции нескольких нижайших связанных состояний атома водорода (а =1) в (1). В этом случае известно точное аналитическое решение задачи [1, 2]. Энергии уровней (собственные значения) равны... Собственные функции (радиальные волновые функции) двух нижайших состояний с l = 0 имеют вид: (здесь первый индекс - главное квантовое число п, второй - l [1, 2]). Вычисленные значения энергий и собственных функций близки к точным с высокой точностью, в частности, для указанных состояний верными являются 14 значащих цифр в собственных значениях и 12 - в собственных функциях (подробное сравнение вычисленных и точных величин см. в разделе 3.2).

3.2.          Задача на собственные значения на конечном отрезке

В данном разделе рассматривается следующая сингулярная задача на собственные значения для уравнения (1) на конечном отрезке: нужно найти такие значения параметра Е, при которых ; нужно также найти решения (1), соответствующие этим собственным значениям Е, и нормированные условием:

                                                                                                                              (24)

Отметим, что к решению такой задачи на собственные значения сводится (в простой модели) задача о состояниях водородоподобного донора, расположенного в центре полупроводниковой сферической квантовой точки (КТ) (полупроводниковый шар малого (~ 10 нм) радиуса в диэлектрической матрице, см., например, [10]).

Задачу на собственные значения в этом случае следует решать численно - аналитическое решение можно найти только при тех значениях r₀, когда узел собственной функции задачи на бесконечном полуинтервале совпадает с радиусом КТ [10]. В данном случае, процедура вычисления собственных значений сводится к численному решению уравнения

                                                                                                                               (25)

где - это регулярные при r = 0 решения (3). Для определённости, при вычислениях полагаем a_l = 1 в (4). После определения некоторого собственного значения Е₀, равенства (3) и (4) дают выражения для ненормированного решения. Для получения нормированной собственной функции нужно вычислить интеграл в левой части (24) от полученного решения (назовём его Q), тогда нормированное решение даётся выражениями (3), (4) с . Интегрирование в левой части (24) с функцией (3) даёт

                                                                                                           (26)

В Таблице 1 представлены результаты расчёта зависимости величин трёх наименьших собственных значений E = E_nl (индекс n соответствует таковому для задачи на бесконечном полуинтервале 0 ≤ r < ∞) от радиуса КТ г₀. Из таблицы, в частности, видно, что индекс n («главное квантовое число») при конечных г₀ уже не нумерует состояния: вырождение состояний с n = 2, l = 0,1 снимается.

Таблица 1. Зависимость энергий трёх нижайших кулоновских состояний от г₀

Как отмечено выше, аналитическое решение задачи для конечного отрезка можно найти тогда, когда узел собственной функции задачи на бесконечном полуинтервале совпадает с радиусом КТ. Из Таблицы 1 видно, что при г₀ = 2 вычисленное наименьшее собственное значение (энергия основного состояния) есть E₁₀ = -0,25 (с высокой точностью - по крайней мере 14 верных значащих цифр). Этому значению соответствует, в частности, решение (1) (ему соответствует также решение (1) , которое и имеет узел при r = 2, т.е. это и есть соответствующая собственная функция, нормированная условием (24), если положить . Точность метода можно оценить, в частности, сравнивая значение величины Q, вычисленное с помощью выражения (26): Q= 0.105306034687423 и с помощью точной формулы Q = 2 - 14e^-2 =0.105306034687422. Сравнение точных и найденных с помощью описанного метода значений собственной функции показывает совпадение 14 значащих цифр в этих величинах при всех г: 0 ≤ r ≤ 2.

3.3.     Волновая функция основного состояния в модели центрально­ячеечного потенциала нулевого радиуса

Модель центрально-ячеечного потенциала нулевого радиуса описывает влияние короткодействующего потенциала, создаваемого примесным атомом в полупроводнике (помимо кулоновского потенциала), на основное состояние в кулоновском потенциале [3]. Математически это означает следующее: нужно найти нормируемое условием (21) решение уравнения (1), соответствующее некоторому заданному значению E = Е0 < 0, не равному собственному значению (1). Ясно, что такое решение есть линейная комбинация регулярного (3) и иррегулярного (5) решений

                                                                                                                      (27)

при l = 0 (поскольку (r) ведёт себя в нуле как ), экспоненциально убывающая на бесконечности. Задача состоит в отыскании константы M при заданном Е0. Для этого нужно провести сшивку решений, т.е. найти нетривиальное решение уравнения (22), в котором матрица A и вектор х есть функции M, а столбец Ф0 = Ф0(M) составлен из значений, принимаемых (27) и его производной в точке r =ȓ при данном M. Процедура вычисления М при заданном Е0 сводится к численному решению уравнения det A(M) = 0.

После определения собственного значения M при данном Е0, решая (22), находим связь между компонентами вектораx (M), а затем, используя (21), - значения констант C2(M) и b-1(M), чем полностью определяем волновую функцию и её производную.

На рис. 1 представлена вычисленная зависимость коэффициента M в линейной комбинации регулярного и иррегулярного решений (1) от энергии Е.

Рис.1. Зависимость коэффициента M в (27) от параметра E (энергии) в (1)

3.4.      Задача на собственные значения в случае потенциала с жёсткой сердцевиной и кулоновским «хвостом»

Модели, использующие потенциал с жёсткой сердцевиной (hardcore potential), рассматриваются в различных приложениях, а методы решения задач квантовой механики при его наличии представляют собой самостоятельный интерес [11]. Модель потенциала с жёсткой сердцевиной и кулоновским хвостом может быть полезна для описания слабосвязанных состояний многоэкситонных комплексов, связанных на примеси в полупроводнике.

Рассматривается следующая задача на собственные значения для уравнения (1): нужно найти такие значения параметра Е: Е < 0 (энергии), при которых решение (1) Rl(r) = 0 при г ≤ гнс (гнс - радиус действия потенциала с жёсткой сердцевиной) и  Rl (r) → 0 при r → ∞; нужно также найти решения, соответствующие этим собственным значениям Е, и нормированные условием:

                                                                                                                        (28)