Диссертация (1150857), страница 6
Текст из файла (страница 6)
Для сохранения нормы волновой функции с достаточной точностью, а значит и для стабильности вычислений необходимо выполнение условия |max | Δ ≪ 1, где max – максимальноесобственное число матрицы H. При этом |max | возрастает с уменьшениемшага пространственной сетки. В результате, для стабильности метода, основанного на разложении (1.89), требуется достаточно маленький интервал Δи, соответственно, довольно большое число шагов по времени.— 43 —Одним из самых распространенных численных методов пропагации волновой функции, сохраняющих её норму, является метод Кранка-Николсон [63].Данный метод основан на использовании приближенного пропагатора]︂−1 [︂]︂[︂Δ −1Δ −11− S H .UCN ( + Δ, ) ≡ 1 + S H22(1.92)Разложив матрицу UCN ( + Δ, ) по степеням Δ, можно обнаружить, чтовсе члены такого разложения до второго порядка включительно совпадают ссоответствующими членами ряда в выражении (1.89).
Таким образом, методКранка-Николсон имеет точность вплоть до третьего порядка по шагу повремени Δ:( + Δ) = UCN () + (Δ3 ).(1.93)Главным преимуществом метода Кранка-Николсон по сравнению с подходом,основанным на разложении матричной экспоненты, является сохранение нормы вектора () пропагатором UCN и, как следствие, высокая стабильность.Выражение (1.93) приводит к следующей системе линейных алгебраическихуравнений (СЛАУ):[︂]︂[︂]︂ΔΔS + H ( + Δ) = S − H () + (Δ3 ).22(1.94)Для получения вектора ( + Δ) необходимо найти решение данной системы, при этом все матрицы в выражении (1.94) являются разреженными, чтооблегчает эту задачу.Методы решения СЛАУ можно разделить на прямые, основанные на разложении матрицы коэффициентов, и итерационные, которые строят последовательность приближенных решений.
Если матрица большая и разреженная,как в системе (1.94), то целесообразно использовать итерационные процедуры [62]. Среди итерационных методов наиболее эффективными для большинства задач являются различные варианты метода бисопряжённых градиентов— 44 —(BiCG, от англ. “biconjugate gradient”) и метода обобщённых минимальныхневязок (GMRES, от англ. “generalized minimal residual”) [70]. Их сходимостьзависит от обусловленности матрицы коэффициентов СЛАУ [62, 70]. Решениедифференциальных уравнений в частных производных с помощью локализованного конечного базиса, как правило, приводит к плохо обусловленнымматрицам, и система (1.94) не является исключением. Поэтому для сходимости итерационного метода необходимо предобуславливание СЛАУ:P−1 A = P−1 .(1.95)Здесь равенство выполняется с точностью до (Δ3 ),ΔA = S + H,2[︂]︂Δ = S − H (),2 = ( + Δ)(1.96)и матрица P−1 – предобуславливатель. Для лучшей сходимости надо подобрать P таким образом, чтобы матрица P−1 A была как можно ближе к единичной.
Если выбрать в качестве P матрицу A при = in :P=S+ΔH(in ),2(1.97)то P и A в произвольный момент времени будут отличаться только за счётзависящих от времени элементов матрицы потенциала. ПредобуславливательP−1 можно получить в явном виде обращением матрицы P, либо можно применить треугольное разложение P = LU, где L – нижняя треугольная матрица, U – верхняя треугольная матрица, при этом действие предобуславливателя на вектор сведётся к решению соответствующих треугольных СЛАУ.Однако, матрицы L и U, так же как и P−1 будут плотными, что крайне нежелательно. Поэтому целесообразно использовать неполное треугольное разложение (ILU, от англ. “incomplete LU factorization”):P ≃ L̃Ũ,(1.98)— 45 —которое является некоторой аппроксимацией точного разложения с разреженными матрицами L̃ и Ũ.
Подобную факторизацию можно осуществитьне единственным образом и, регулируя степень разреженности L̃ и Ũ и, соответственно, точность ILU-разложения, можно улучшить быстродействие процедуры решения СЛАУ [70]. Проведя необходимую ILU-факторизацию матрицы (1.97) один раз, далее её можно использовать для предобуславливаниясистемы (1.94) на каждом шаге по времени. Однако, при изменении шага повремени ILU-разложение придётся пересчитать.В представленной работе при расчётах процесса перезарядки, как правило, использовался метод, основанный на разложении матричной экспоненты (1.89), при этом учитывались первые четыре члена разложения.
Однако, было обнаружено, что для вычисления вероятностей рождения электронпозитронных пар (см. глава 2) данный метод слишком нестабилен, и потомуиспользовался метод Кранка-Николсон. Для решения СЛАУ при этом применялся BiCGS (от англ. “biconjugate gradient squared”) алгоритм, описанныйв работе [71], вместе с ILU-предобуславливателем.1.4Результаты расчётовВ данном разделе представлены результаты расчётов вероятности перезарядки в столкновениях голого ядра (снаряд) с водородоподобным ионом (мишень). Вычисления проводились с помощью метода, изложенного в разделе 1.3, как в приближении центрального столкновения (одноканальный расчёт), так и с учётом вращения межъядерной оси для разного числа каналов с различной проекцией углового момента электрона на межъядернуюось.Размеры рабочей области в случае столкновения U92+ −U91+ (1) состав-— 46 —ляли = 14000 фм по оси и = 3500 фм по оси . В остальных случаях и масштабировались пропорционально 1/ , где – заряд ядра.
Базисные сплайны Эрмита строились на равномерной сетке с = 61 и = 16.Таким образом, на каждый канал с фиксированной проекцией полного углового момента электрона приходилось по 16 × ( − 1) × = 15360 базисныхфункций. Было обнаружено, что такого количества достаточно для получения сходящихся результатов.Мишень была зафиксирована на оси в точке, сдвинутой на /4 от краярабочей области. Начальное состояние мишени находилось с помощью итерационной процедуры, описанной в параграфе 1.3.4.
В начальный момент времени in снаряд располагался на оси в точке с координатой (in ) = −1.2×(положение мишени считается началом координат), таким образом, он начинал движение за пределами рабочей области. Двигаясь вдоль оси , снарядпосле сближения с мишенью менял направление движения и в конечный момент времени out достигал точки (out ) = − /2. Волновая функция электрона развивалась во времени методом, основанным на разложении матричной экспоненты (см.
параграф 1.3.5), вплоть до момента out . Затем электронная плотность интегрировалась по половине пространства, содержащейснаряд, полученное значение считалось вероятностью перезарядки.Шаг по времени Δ был равномерным, число шагов варьировалось от 104до 105 . При расчётах для малых значений прицельного параметра шаг по времени приходилось уменьшать для преодоления неустойчивости метода пропагации.На рисунках 1.7, 1.8 и 1.9 представлены полученные результаты для вероятности перезарядки как функции прицельного параметра в столкновениях Xe54+ −Xe53+ (1), Pb82+ −Pb81+ (1) и U92+ −U91+ (1), соответственно. Вычисления проводились для энергии столкновения = 6 МэВ/а.е.м— 47 —(Pb, U) и 5.9 МэВ/а.е.м.
(Xe), что примерно соответствует энергии кулоновского барьера. Для всех трёх систем приведены значения, рассчитанные как в рамках приближения центрального столкновения ( = 1), таки с учётом вращения межъядерной оси. В последнем случае для столкновений Xe54+ −Xe53+ (1) и Pb82+ −Pb81+ (1) представлены данные, полученные при = 4, для столкновения U92+ −U91+ (1) – данные, полученныепри = 6. Также для Xe54+ −Xe53+ (1) и U92+ −U91+ (1) приведены соответствующие значения из работы [33], полученные методом ДФШ. В случае столкновения Xe54+ −Xe53+ (1) (см. рисунок 1.7) кривая для = 1качественно повторяет поведение двух других, но сдвинута относительноних в сторону увеличения значений прицельного параметра. Для столкновений Pb82+ −Pb81+ (1) и U92+ −U91+ (1) при больших значениях все кривыеведут себя сходным образом, но при малых значениях данные для приближения центрального столкновения отличаются от результатов многоканального расчёта кардинально.
Эта разница особенно ярко выражена для столкновений ионов урана. Отсюда можно сделать вывод о неприменимости приближения центрального удара для расчёта перезарядки в столкновениях тяжёлых ионов при малых значениях прицельного параметра. При этом результаты вычислений с учётом вращения межъядерной оси для Xe54+ −Xe53+ (1)и U92+ −U91+ (1) хорошо согласуются с соответствующими значениями из работы [33].С целью продемонстрировать работу метода при различных значенияхэнергии столкновения , также были проведены расчёты () в столкновении Xe54+ −Xe53+ (1) при = 10 и 20 МэВ/а.е.м.
Полученные результаты вместе с соответствующими значениями, рассчитанными для =5.9 МэВ/а.е.м, представлены на рисунке 1.10. Как можно увидеть из рисунка,количество осцилляций кривой () уменьшается с ростом . Данный факт— 48 —согласуется с адиабатической теорией перезарядки [46–48], в соответствиис которой () определяется осцилляциями электронной плотности междуядрами, и период этих осцилляций зависит только от межъядерного расстояния.Для того чтобы продемонстрировать сходимость расчётов по числу каналов с различным значением проекции углового момента электрона намежъядерную ось, на рисунке 1.11 показаны результаты вычислений () c = 2, 4, 6. Из графиков видно, что качественное отличие от значений, полученных в приближении центрального удара, проявляется уже при = 2.При этом разница между результатами для = 4 и = 6 практическиотсутствует.
Отсюда можно сделать вывод, что сходимость по числу каналовдостигнута уже при = 4.1.5ВыводыВ данной главе диссертации представлен новый метод расчёта вероятностейперезарядки в столкновениях тяжёлых многозарядных ионов. Подход основан на решении нестационарного двухцентрового уравнения Дирака в конечном базисе. Для построения базисных функций используются кубическиебазисные сплайны Эрмита, определенные на двумерной равномерной сеткев системе координат, вращающейся вместе с межъядерной осью.
Для учётавращения системы координат рассматривается несколько каналов с различной проекцией углового момента на межъядерную ось.С помощью разработанного метода были проведены расчёты вероятности перезарядки в низкоэнергетических столкновениях Xe54+ −Xe53+ (1),Pb82+ −Pb81+ (1) и U92+ −U91+ (1) для различных значений прицельного параметра. Вычисления проводились как в рамках приближения центрального— 49 —столкновения, так и с учётом вращения межъядерной оси с разным числомканалов с различной проекцией углового момента электрона на межъядерную ось. Из сравнения полученных результатов можно сделать вывод о том,что учёт вращения межъядерной оси необходим для корректного описанияпроцесса перезарядки в низкоэнергетических столкновениях тяжёлых ионов.При этом четырёх каналов c различными значениями уже достаточно дляполучения сходящихся результатов.
Данный факт делает разработанный метод более предпочтительным для расчёта различных процессов в подобныхстолкновениях, чем подходы, основанные на использовании трехмерной сетки.Результаты,полученныедлястолкновенийXe54+ −Xe53+ (1)иU92+ −U91+ (1) хорошо согласуются с соответствующими значениями изработы [33], рассчитанными с помощью метода ДФШ. При этом для столкновения U92+ −U91+ (1) удалось подтвердить нетривиальную особенностьзависимости вероятности перезарядки от прицельного параметра при малыхзначениях последнего.Основные результаты, представленные в данной главе, были опубликованы в работах [72] и [73].— 50 —1This work, 4 channelsThis work, 1 channelTupitsyn et al.0.8P(b)0.60.40.20050010001500200025003000350040004500500055006000b (fm)Рис. 1.7: Вероятность перезарядки в зависимости от прицельного параметра в столкновении Xe54+ −Xe53+ (1) при энергии = 5.9 МэВ/а.е.м.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
