Диссертация (Единый аналитический и вычислительный подход к решению квантовой задачи трёх тел), страница 11

Вычисление матричных элементовоператоровПрименим построенные базисные функции (, , ) для решения урав­нения (2.1) методом Галёркина. При численной дискретизации, бесконечнаяобласть его задания [0, ∞] × [0, ∞] × [0, ] аппроксимируется конечной Ω =()()[0, max ] × [0, max ] × [0, ]. Область Ω разбивается на конечное число паралле­лепипедов так, что⋂︀∙ и , ̸= , пересекаются только по границам : ( ) ( ) = ∅.∙ Объединение всех параллелепипедов совпадает с областью Ω:65⋃︀=1 = Ω.Поскольку в уравнении (2.1) имеются внедиагональные блоки, связывающиеразличные компоненты Ψ ′ волновой функции, удобно использовать одно и тоже разбиение для всех компонент.

Базисные функции при этом могут как совпа­дать для разных компонент, так и зависеть от номера компоненты, ′ (, , ).Дискретное решение [ΨFEM ] уравнения (2.1) раскладывается по построенномуглобальному базису функций ′ (, , ) как[ΨFEM ]′ =∑︁ ′ ′ ( , , ).(2.19)=1В дальнейшем индекс номера компоненты у базисных функций будет опускать­ся.В силу свойств разбиения и аддитивности интеграла, матричные элементыпроизвольного оператора O(, , ) записываются в виде:ZT == ( , , )O( , , ) ( , , ) =Ω Z∑︁ ( , , )O( , , ) ( , , ) .(2.20)=1 Проекции функций (, , ) на элемент представляются как произведенияфункций по каждой из координат (2.18). Таким образом, матричные элемен­ты операторов, представимых в виде O( , , ) = O() ( )O() ( )O() ( ),можно записать как произведение одномерных интегралов:T = Z∑︁=1Z()()ℎ, ( )O() ( )ℎ, ( )2 ( )2 ( ) ×()()()()()ℎ, ( )O() ( )ℎ, ( )2 ( )2 ( ) ××()Zℎ, ( )O() ( )ℎ, ( ) sin .×(2.21)()Заметим, что все операторы, отвечающие кинетической энергии в уравнении (2.1),являются факторизуемыми.

Единственный оператор, не обладающий этим свой­66ством – оператор умножения на потенциал. Вычисление матричных элементовэтого оператора требует специального рассмотрения.2.3.4. Спектральное разложение по угловой переменнойВычислительные эксперименты по исследованию связанных состояний раз­личных квантовых систем показывают, что сходимость по угловой переменнойприближенных решений к точным является недостаточно быстрой. Гибкость ко­нечноэлементного подхода, однако, позволяет создать комбинированный метод,совмещающий МКЭ по пространственным переменным , , и спектральныйметод [126] по переменной .Разложим компоненты волновой функции Ψ , как функции угла , в рядпо сферическим гармоникам с нулевым вторым аргументом:[ΨFEM ] ( , , )=+ −1∑︁Φ ( , ) ( , 0) =+ −1∑︁=Φ ( , ) (cos ).=(2.22)В последнем уравнении использован факт пропорциональности таких сфериче­ских гармоник присоединённым полиномам Лежандра (cos ):√︃ ( , 0) = (cos ), =2 + 1 ( − )!.4 ( + )!(2.23)При выборе такого разложения для компонент, граничное условие (1.26) выпол­няется автоматически.

Нужно отметить, что количество элементов разложениядля каждой из компонент одинаково, но разлагаются эти компоненты по раз­ным наборам функций, т.е. базисные функции зависят от номера компоненты.При подстановке разложения (2.22) в уравнение (1.22), получается следу­67ющая бесконечная система двумерных дифференциальных уравнений:)︂ (︂)︂(︂)︂]︂( + 1) − 2211 21 21−−+ − ( + 1) 2 + 2Φ +222 √︀√︀+ (, )+ (, ) +1− (, )− (, ) −1+ 1 + 0Φ−Φ +1+122+ −1∑︁(, )Φ = Φ , = 0[1] .

. . , = , + 1, . . . .(2.24)+[︂(︂=Для упрощения записи, здесь опущен индекс . Матричные элементы потенци­ала (, ) определяются какZ(, ) = 2 sin (, 0) (, , ) (, 0).(2.25)02.3.5. Вычисление матричных элементов потенциалаМатричные элементы потенциалов (, ) могут быть просто вычисле­ны для потенциалов, зависящих только от расстояний между частицами. Дляэтого разложим полный потенциал (, , ) в ряд по полиномам Лежандра откосинуса угла : (, , ) =∞∑︁ (, ) (cos ), где (, ) ==02 + 12Z1 cos (, , ) (cos ).−1(2.26)Тогда вычисление матричных элементов (, ) сводится к вычислению инте­гралов(, ) = 2∞∑︁Z1 (, )=0 cos (, 0) (cos ) (, 0).(2.27)−1Для вычисления этих интегралов, воспользуемся дважды представлением дляпроизведения сферических гармоник [79]:√︃1 1 (, 0)2 2 (, 0) =∑︁(21 + 1)(22 + 1) 01 +21 02 0 ,1 +2 (, 0).1 1 2 24(2 + 1)(2.28)68Сумма в данном представлении конечна в силу треугольного свойства коэффи­циентов Клебша-Гордана.Нужно отметить, что представление (2.28) приводит к обобщению хорошоизвестной формулы для произведения полиномов Лежандра [98] на присоеди­нённые полиномы Лежандра:21() =()21∑︁1 +21 2(), 1 2 где√︃1 122 =(21 + 1)(22 + 1) 1 +2 01 +2.1 02 0 1 1 2 24(2 + 1) 1 1 2 2Используя соотношение (2.28) для (cos ) (, 0),1 ∑︁ (cos ) (, 0) =0√︃(2 + 1)(2 + 1) 0 00 0 (, 0),4(2 + 1)и свойство ортогональности, найдём для интеграла следующее значение:Z11 cos (, 0) (cos ) (, 0) =2√︂2 + 1 0 .2 + 1 00 0(2.29)−1Таким образом, матричный элемент (2.25) может быть записан в виде:(, ) =+∑︁√︂ (, )=|−|2 + 1 0 .2 + 1 00 0(2.30)Важно, что для фиксированных значений индексов , и сумма в полученномвыражении содержит конечное число слагаемых, следовательно, в вычисленияне вносятся дополнительные погрешности.

Если в разложении (2.22) использу­ется слагаемых, то для вычисления всех матричных элементов по форму­ле (2.25) требуется вычисление ( + 1)/2 интегралов по угловой перемен­ной для фиксированных , для каждой компоненты . При использованиипредставления (2.30), общее требуемое количество угловых интегралов равно2( + ) − 1, причём эти интегралы не зависят от номера компоненты. В ре­альных применения, ≫ 1 и ≫ , в этом пределе количество вычислений69угловых интегралов сокращается примерно в /4 раз, причём их можно вы­числять один раз, а не для каждой компоненты волновой функции отдельно.Для произвольного потенциала, разложение (2.26) должно быть полученос помощью численного интегрирования.

Однако, для некоторых специальныхтипов потенциалов может быть найдено аналитическое представление коэффи­циентов (, ). Например, такое представление может быть получено для важ­ного случая парного кулоновского взаимодействия, функциональная часть ко­торого равна (, , ) =11= √︀.| x − y|( )2 + ( )2 − 2 cos Коэффициенты в этом случае равны1 (, ) = (),| | где =, () =⎧⎨ ,⎩ −(+1) ,|| ≤ 1, > 1.(2.31)Такое представление существенно увеличивает как скорость, так и точностьвычисления матричных элементов кулоновского потенциала.2.4. Оценки погрешности численного методаВ рамках МКЭ доступны оценки погрешностей как для решений уравне­ний, так и для их производных [109].

Рассмотрим семейство регулярных аф­финных конформных элементов с максимальным диаметром элемента ℎ, вклю­чающее все произвольные полиномы степени не выше . Предположим также,что каждая компонента точного решения Ψ уравнения (2.1) обладает достаточ­ной гладкостью, а именно принадлежит пространству Соболева H+1 (Ω). Тогдапри выполнении некоторых дополнительных, естественных в рассматриваемомслучае предположений [109], существует такая не зависящая от ℎ постоянная , что||Ψ − ΨFEM ||, Ω ≤ ℎ+1− |Ψ|+1, Ω ,70 = 0, 1.(2.32)Здесь ΨFEM – дискретное решение, || · ||, Ω – норма в векторном по индексукомпонент пространстве H (Ω), а | · |+1, Ω – полунорма порядка + 1 в такомже пространстве H+1 (Ω). Оценка (2.32) является априорной, то есть позво­ляет найти оценку погрешности без знания дискретного решения ΨFEM . Приэтом она неконструктивна в том смысле, что вычислить (и даже оценить) еёправую часть не представляется возможным.

Однако на основе формул такоготипа можно построить апостериорные оценки погрешностей [127], которые ужеявляются конструктивными.2.4.1. Экстраполяционные формулыАприорная оценка (2.32) может быть также использована для построенияэкстраполяционных формул. С помощью таких формул значения некоторыхматричных элементов решения, найденные для нескольких фиксированных па­раметров численной схемы, могут быть пересчитаны для других (обычно беско­нечно больших или малых) значений параметров, предлагая значения, прибли­жающие сошедшиеся по этим параметрам значения.Простая экстраполяционная формула была предложена в работе [128] длявычисления уровней атомных систем:ext = −1 +(−1 − −2 )( − −1 ),2−1 − −2 − (2.33)где – энергия, вычисленная с полиномами степени, не большей, чем .

Форму­ла (2.33) даёт точную экстраполяцию по , если собственное значение уравненияведёт себя как = ext + const ℎ2(2.34)при больших значениях . Формула (2.33) была использована при некоторыхвычислениях уровней энергии антипротонного гелия в 3.2 и резонансов атомагелия в 4.2.Сравнивая выражения (2.32) и (2.34), легко заметить, что множитель пе­ред ℎ+1− в (2.32) зависит от , так что экстраполяция при → ∞ (2.33)71может давать достаточно большую погрешность.

Более аккуратная оценка по­лучается, если искать экстраполяцию в виде функции, зависящей от количествастепеней свободы, т.е. количества глобальных базисных функций , в МКЭ.Собственные значения можно искать в виде [111]: = ˜ext + ()−2 ,(2.35)где () – количество базисных функций, отвечающее решению задачи с по­линомами степени , а – неизвестный параметр. Вычислив энергии длятрёх значений = − 2, − 1 и , затем можно исключить константы и ˜ext нелинейное уравнениеиз уравнений (2.35) и найти для энергии ˜ext − =˜ext − −1(︃˜ext − −1˜ext − −2(−1)/ ()))︃ loglog(((−2)/(−1)).(2.36)Это уравнение сводится к уравнению (2.33) в случае, когда показатель степенив правой части равен 1, т.е.

( − 1)/ () = ( − 2)/ ( − 1). Уравне­ние (2.36) применялось при вычислениях тримера гелия в разделе 3.3.1, гдезначения показателя степени оказалось примерно равным 0.8. Таким образом,экстраполяционные формулы (2.33) и (2.36) могут приводить к различающим­ся результатам. Обсуждение этого различия можно найти в 3.3.1.

Важно такжеотметить, что приведённые экстраполяционные формулы могут быть использо­ваны не только для вычисления поправок к энергиям, но и, например, к сред­неквадратичным радиусам и другим матричным элементам.2.4.2. Оценки погрешности и адаптивный подходПрименение численных методов для решения физических задач требует,как правило, проведения серии расчётов с различными параметрами численнойсхемы. Анализируя семейство полученных приближенных решений, можно какоценить их точность, так и использовать полученную информацию для выбо­ра оптимальных параметров при расчёте нового элемента семейства.

В рамках72МКЭ, применяются два основных итеративных подхода к уточнению решений:уменьшение (или в общем случае – изменение) размера конечных элементов(ℎ-уточнение) и увеличение степени полиномов на КЭ (-уточнение). Эти под­ходы независимы друг от друга, так что они могут применяться совместно врамках так называемого ℎ-уточнения [111]. В этом подходе, предварительноеℎ-уточнение обеспечивает равномерное распределение ошибок по КЭ всей три­ангуляции, а требуемая окончательная точность достигается за счёт финаль­ного -уточнения.