Научная статья на тему 'О вычислении атомных интегралов с экспоненциально коррелированными функциями'

О вычислении атомных интегралов с экспоненциально коррелированными функциями Текст научной статьи по специальности «Математика»

CC BY
158
33
i Надоели баннеры? Вы всегда можете отключить рекламу.
Ключевые слова
AТОМНЫЕ ИНТЕГРАЛЫ / ТРЕХЧАСТИЧНЫЕ КУЛОНОВСКИЕ СИСТЕМЫ / ПРОИЗВОДНЫЕ ФУНКЦИИ ФАДДЕЕВОЙ / УСТОЙЧИВЫЙ АЛГОРИТМ РАСЧЕТА / ATOMIC INTEGRALS / THREE-PARTICLE COULOMB SYSTEMS / DERIVATIVES OF FADDEEVA FUNCTION / STABLE ALGORITHM FOR CALCULATION

Аннотация научной статьи по математике, автор научной работы — Нечаев В. В., Зиганшина О. Д., Сучкова Н. К.

Рассмотрен новый тип корреляционных атомных интегралов, возникающих в вариационных расчетах энергии трехчастичных кулоновских систем. Подынтегральная функция в них наряду с линейным членом по межчастичному расстоянию под экспонентой дополнительно содержит квадратичный член. Показано, что эти интегралы аналитически выражаются через функцию Фаддеевой чисто мнимого аргумента и ее производные. Разработан устойчивый и быстрый алгоритм для вычисления производных функции Фаддеевой до двадцатого порядка. Даны тестовые значения исследованных специальных функцийA new type of correlation atomic integrals occurring in variation energy calculations of three-particle

i Надоели баннеры? Вы всегда можете отключить рекламу.
iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

Coulomb systems is studied. A integrand in them along with an interparticle distance linear term under an exponent additionally contains a quadratic term. It is demonstrated that these integrals are analytically expressed through Faddeeva function of a pure imaginary argument and its derivatives. A stable and fast algorithm for calculation of Faddeeva function derivatives to the twentieth order is developed. The test values of the studied special functions are provided

Текст научной работы на тему «О вычислении атомных интегралов с экспоненциально коррелированными функциями»

УДК 539.182/.184, 519.677

О ВЫЧИСЛЕНИИ АТОМНЫХ ИНТЕГРАЛОВ С ЭКСПОНЕНЦИАЛЬНО КОРРЕЛИРОВАННЫМИ ФУНКЦИЯМИ

В. В. Нечаев, О. Д. Зиганшина, Н. К. Сучкова

Саратовский государственный технический университет им. Ю. А. Гагарина Email: VL-Nechaev@yandex.ru

Рассмотрен новый тип корреляционных атомных интегралов, возникающих в вариационных расчетах энергии трехчастичных кулоновских систем. Подынтегральная функция в них наряду с линейным членом по межчастичному расстоянию под экспонентой дополнительно содержит квадратичный член. Показано, что эти интегралы аналитически выражаются через функцию Фад-деевой чисто мнимого аргумента и ее производные. Разработан устойчивый и быстрый алгоритм для вычисления производных функции Фаддеевой до двадцатого порядка. Даны тестовые значения исследованных специальных функций.

Ключевые слова: aтомные интегралы, трехчастичные куло-новские системы, производные функции Фаддеевой, устойчивый алгоритм расчета.

Calculation of Atomic Integrals with Exponentialy Correlated Functions V. V. Nechaev, O. D. Ziganshina, N. K. Suchkova

A new type of correlation atomic integrals occurring in variation energy calculations of three-particle Coulomb systems is studied. A integrand in them along with an interparticle distance linear term under an exponent additionally contains a quadratic term. It is demonstrated that these integrals are analytically expressed through Faddeeva function of a pure imaginary argument and its derivatives. A stable and fast algorithm for calculation of Faddeeva function derivatives to the twentieth order is developed. The test values of the studied special functions are provided. Key words: atomic integrals, three-particle Coulomb systems, derivatives of Faddeeva function, stable algorithm for calculation.

to n

При проведении высокоточных вариационных расчетов энергии трехчастичных кулоновских систем возникает необходимость в быстром вычислении следующих шестикратных интегралов [1]:

I(і,у,к) = |г/г/г* ехр(-от^ - /Згп -Ьг1 - сг2)іУ1іУ2, (1)

где, і, ], к - целые числа; г1, г2 - модули радиус-векторов точек 1 и 2 относительно некоторой точки трехмерного пространства; г12 - расстояние между точками 1 и 2; а, в, Ь, с - постоянные действительные коэффициенты; С¥1, с1У2- элементы объема для векторов г1, г2 соответственно.

В исходном выражении (1) можно провести интегрирование по угловым переменным с использованием преобразования, предложенного в [2]. Перейдем к сферическим координатам для каждого трехмерного интеграла и свяжем полярную ось системы координат для внутреннего интеграла с радиус-вектором г1 для внешнего интеграла. Учитывая, что

Гп = г1 + г/ - 2г1г2 соя(в),

Г Г

іг12 = —— і (с оя (в))

Г12

получаем

где

I(i, j, k) = | d(px I ddx sin(e) I dp2R(i, j, k) = 8n2R(i, j, k),

о о о

“ “ r1 +r2

I( i, j, k) = | dr2 | dr2 | dr22r/+1 ri+1r1k2+1 exp (-ar2 - /3^2 - br2 - cr2) .

(2)

0 0 — -г2

Интеграл (2) сходится при условии

((а > 0) а (Ь + с > 0)) V ((а

Докажем это утверждение для 1(—1,—1,— 1). Воспользуемся периметрическими координатами [3]. Связь г1, г2 , г12 с новыми координатами в общем случае задается линейным преобразованием:

— = а2р2 + а12р12,

Г2 = а12р12 + а1р1 , г12 = а1р1 + а2р2,

= о) л (b + с > о) л (b + в> о) л (с + в> о)) .

(З)

где а1, а2, а12 - постоянные, отличные от нуля коэффициенты. Требование отличия коэффициентов от нуля необходимо для биективности преобразования. Положим а1 = а2 = а12 = 1. Якобиан перехода в этом случае имеет вид

0 а2 аг 3 = ёй а2 0 а12 = 2а1а2а12 = 2.

а а 0

©1) Нечаев В. В., Зиганшина О. Д., Сучкова Н. К., 2012

Основное достоинство периметрических координат заключается в том, что они позволяют перейти от зависимых пределов в интеграле (2) к независимым. В таких координатах условие треугольника

\Г2 - Гп\ < Г1 < Г2 + Г12 ,

\Г12 - Г\ < Г2 < Г12 + Л,

| г - т2 \ < тп < г + т2 выполняется автоматически для любых сочетаний периметрических координат у91 , Р2 , Р12 , т. е.

они независимы. Так же как и межчастичные расстояния, периметрические координаты положительны.

Введем вспомогательный производящий интеграл:

I(-1,- 1 ,-1) = -в—12 - ЬГ1 - ) іV^ .

■’ Г1 Г2 Г12

Проведем интегрирование по угловым переменным как в уравнении (2) и преобразуем его к периметрической системе координат:

га —]-г —

I(-1,-1,-1) = 8ж21 ехр(-аг^2 - /Зг12)йг12| ехр(-Ьг1)іг1 | ехр(-сг2)іг2 =

0 0 — -—21

: 16л2 Л ірір2 ехр(-а(р1 +р2)2 - (в + Ь)р1 - (в + с)р2)^ір12 ехр(-(Ь + с)р12).

(4)

Из этого выражения следует необходимое условие сходимости Ь + с > 0.

Теперь рассмотрим два случая:

1) а = 0.

Интеграл |^йрхйрг ехр(-(в + Ъ)рх - (Р + с)р2)

0 0

представляется в виде двух сомножителей

Iйрх ехр(-(в + с)рх) и |йр2 ехр(-(в + Ь)р2), откуда

00

вытекают еще два условия сходимости - в + Ь > 0 и в + с > 0.

2) а ф 0. Очевидно, интеграл (4) расходится при а < 0. При а < 0 интеграл

11йрхйр2 ехр{-а(рх +р2)2 - {Р + Ъ)рх - (в + с)р2)

0 0

всегда меньше, чем интеграл

IIёрхйр2 ехр(-арх -ар2 -(в + Ъ)рх - (в + с)р2) ■

0 0

Данный интеграл факторизуется на интегралы

|йрх ехр(-арх - (в + Ъ)Рх),

|йр2 схр(-ар1- (р + с)р2),

0

которые, при сделанном предположении а < 0, сходятся независимо от величин в + с и в + Ь. Снизу величина интеграла (4) ограничена нулем, что доказывает сходимость 1(—1,—1,— 1).

Интегралы I (I, у', к) с I > -1, у > -1, к > - 1 определяются дифференцированием производящего интеграла 1(— 1,— 1,— 1) по параметрам в, Ь, с:

у+х / - У+х

I (і, у, к) = 1-

дЬ

_д_

дс

Интегралы (1) хорошо известны для случая а = = 0 [2-6]. Далее мы покажем, что в более общем случае а ф 0 они аналитически выражаются через функцию Фаддеевой [7] чисто мнимого аргумента и ее производные, а также рассмотрим устойчивый алгоритм их расчета.

Найдем аналитическое выражение для 1( 1, 1, 1):

^ ^ ч , .-2

I(-1,-1, -1) = 8л21 ехр(-агХ2 - в—Х2 )йгХ2| ехр(-Ьг1 )іг1 | ехр(-сг2 )йг2

0 0 — -—21

= %п21 ехр(-аг22 - Рг12 )р(гп )іг12 =

0

■((( У) - Я0( х ))=-

(6)

8пг

4п

а(Ь + с)(х - у)

где х = (Ь + в) / 2 Vа, у = (с + в) / 2^~а, і = в /^а ,

ау[а(х + у - і) I (х - у)

£0( У) - 80( х)

ф(х) = | ехр(-Ьг) I ехр( - су )іу = 2(Єхр( Ьх2Є2Хр( сХ)^ , ?0 (х) = єхр(х" ) } єхр( - г2 ) і2 .

0 Ь-х1 ( с Ь ^ х

Последняя функция возникает при интегрировании экспоненты от квадратичного полинома:

|ехр(-г2 - 2xz)dz = g0(х).

(7)

ё0(х)- это вспомогательная функция, которая = носит название масштабированной дополнительной функции ошибок. Некоторые свойства g0(x) можно найти в [8], где рассмотрена функция комплексного аргумента w(z) =

= ехр(- г2) ег^с (— гг) , так называемая функция Фаддеевой. Для чисто мнимого аргумента w(z) с точностью до множителя совпадает с ё0( х). Для случая комплексного аргумента эффективный алгоритм ее вычисления был разработан в работе [9] и недавно улучшен [10].

При дифференцировании интеграла

1(—1,—1,— 1) по параметрам в, Ь, с (5) возникают производные от функции ё0( х) высокого порядка

ё п(х ) = (ё0(х ))( „ ^

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

эффективное и точное вычисление которых представляет самостоятельную проблему. Специфика прикладных задач [1] такова, что необходимо рассчитывать десятки тысяч интегралов (1) для различных значений параметров а в, Ь, с с заданной точностью, поэтому использование имеющихся систем аналитических вычислений (Мар1е [11], Ма&ешайса [12] и др.) ограничено тестовыми расчетами для отдельных частных случаев. В связи с этим нами была поставлена задача изучения свойств ё „(х) и разработки быстрого, устойчивого алгоритма их вычисления.

Данные функции знакопостоянны (четные -положительны, нечетные - отрицательны):

g2)(х) > 0, g2х) < о, ] > о. (8)

Это свойство следует из формулы, получаемой дифференцированием (7):

gn (х) = (-2)п | zn ехр(-z2 - 2xz^,

о

так как подынтегральное выражение всегда неотрицательно. Функции ё „( х) подчиняются рекуррентным соотношениям:

gl( х) = 2 ^о( х) -1, (9)

gn (х) = 2^п-1(х) + 2(п - 1)gn-2(х) (10)

и обладают следующими свойствами на бесконечности:

Иш gn (х)х- ехр(-х2) = 2п4Л,

х——-^ !

Иш gn (х)хп+1 = (-1)пП-. (11)

х—+^ 2

При х = 0 можно получить простые соотношения g2„(0) = 2п-1(2п- 1)\\ '-1П , g2n+l(0) = -2п(2п)!! .

Функцию g„ (х) можно также представить в виде

gn (х) = gо(х)Кп (х) + Кп (х) , (12)

где К„(х), Ь„(х) - полиномы, как и g„( х), подчиняющиеся рекуррентной формуле (9), причем К0(х) = = 1, К1(х) = 2х ; ь0(х) = 0, Ь,(х) = —1. К„(х) выражаются через полиномы Эрмита комплексного аргумента К„(х) = (—1)„ Н„ (гх). Первое слагаемое правой части отношения

gо( х) = -Кп ( х)/Кп ( х) + gn ( х)/Кп ( х) (13)

является оценкой сверху (нечетные „) и снизу (четные „) для g 0(х), а при больших „, х оно же может использоваться для подсчета g 0( х). Для доказательства этого рассмотрим систему рекурсий и ограничений

хп = апхп-2 + Ьхп-1, (14)

Уп = апУп-2 - Ьу„-1, (15)

Ь, хп, Уп > 0 , (16)

ап > ап-1 >... > шах{0,1 - Ь}, (17)

т. е. Ь, х„ , неотрицательны,у„ , а„ положительны, {ап } строго монотонно возрастает.

Лемма 1. Последовательность {хп } при

х0 < х1 Ь > 1 - а2 х0/ х1 Ь > 1 - а2 х0/х1, (18) удовлетворяющая (14), (16), (17), строго монотонно возрастает, начиная со второго члена.

Доказательство. Из условий (18) и соотношения (14) следует х1 - х0 > 0 , х2 -х1 = а2х0 + (Ь- 1)х1 > 0 .Строгая монотонность для следующих членов доказывается по индукции

хп - хп-1 = апхп-2 - ап-1 хп-3 +

+ Ь(хп-1 - хп-2 ) > ап-1(хп-2 - хп-3 ) +

+ Ь(х П-1 - хп-2 ) > 0.

Здесь использовано условие (17) для замены а„

на а„—1.

Лемма 2. Последовательности {х27 /у27}, {х27+1/у2}+1}, определяемые системой рекурсий (14), (15) и ограничений (16), (17), являются неубывающими.

Доказательство. Система неравенств

- +

Ьхп-1

хп = апхп-2 + Ьхп-1 >

Уп апУп-2 - ЬУп-1 Уп-2 апУп-1 Уп-2

доказывает лемму.

Лемма 3. Если последовательность частичных сумм

*п =Х 1/ а 7 (19)

7

неограничена сверху, то последовательности {хп }, {хп / Уп }, определяемые системой рекурсий

(14), (15) и ограничений (16), (17), (18), неограни-чены сверху.

Доказательство. Из уравнения (15) и условий (16), (17) следует, что

Уп < У

[я /2]

.П‘

2

где т = шоё(п,2), [п/2] - целая часть числа п/2.

Последовательность {хя} строго монотонно возрастает по условию (18) и

[п/2]

Хп > (ап + Ь)Хп-2 , хп > х0 П (а27+т + Ь) '

(20)

7=1

Тогда

х

[п/2]

> кП

Уп 7=1

/

1 +

V

2 ]+т

[п /2]

> Ьк X

7=1

(21)

2 7+1

где к = тт{хиХо}/Ут .

Последовательность частичных сумм в (21) является неограниченной при я ^ ж, как и последовательность (19). Последовательность {хп} неограниченна в силу выполнения неравенств (17) и (20).

Следствие 1. Прямая рекурсия (9), (10) для функции gn (х) неустойчива при х > 0, а знак ошибки одинаков для всех gn (х).

Доказательство. Обозначим через (х) сумму точной функции gn(X) и

ошибки ЕП (х ): ^п (х) = gn (х) + £п (х) . В соот-

ветствии с рекурсией для первого члена (9) ошибка £1 (х) = 2х£0 (х) имеет тот же знак, что и £0(х), согласно (10) этот знак сохраняется для всех еп (X). Поведение данной функции подчиняется уравнению (14). Последовательность {(-1)” gn(х)} удовлетворяет (15). Таким образом, ошибка еп (х) не ограничена сверху (лемма 3) и строго монотонно возрастает (лемма 1) при 2х > 1 -2е0(х)/2х£1(х) = 1 -1/х. Поскольку действительных корней уравнение 2 х2 = х -1 не имеет, ошибка растет во всей рассматриваемой области (х > 0) . Относительная ошибка \еп (х)/ gn (х) возрастает с увеличением п (лемма 2) и также неограниченна (лемма 3).

Следствие 2. Последовательность дробно-полиномиальных функций

{- Ь” (х) / К” (х)} (22)

сходится при х > 0 к g0(х). При этом последовательность четных членов сходится к пределу сверху, а нечетных - снизу.

Доказательство. Перепишем (13) в виде

- Ьп (х) / Кп (х) = g0(х)(1 - gn(х)/Ьп (х))-1.

Знакопостоянная последовательность {- Ьп (х)} удовлетворяет (14). В качестве последовательности, присутствующей в (15), снова

будем рассматривать {(- 1)п gn (х)} .Тогда из лемм 2 и 3 следует

НЬп(х)/gn(х)\ = жи Ит(-Ьп(х)/Кп(х) = gо(х).

пп

С учетом неравенств (8) получаем

- Ь2 7 ( х)/ К2 7 ( х) < gо( х) <- Ь2 7+1( х)/ К2 7+1( х) ,

\Ш1(-Ь27 (х)/К27 (х) = gо(x) + 0 ,

^т-Ь27+1(х) / К27+1(х) = gо(х) - 0 .

В силу леммы 2 последовательности

{Ь27(х) /g27(х)} , {Ь27+1(х) /g27+1(х)}

монотонны, а значит, и последовательности

{-Ь27(х) /К27(х)} , {-Ь27+1 (х) /К27+1(х)}

монотонно приближаются к пределу g 0( х) с разных сторон.

В соответствии со следствием 1 знак ошибки аппроксимации при использовании прямой рекурсии для приближенного вычисления g 0(х) сохраняется для всех gn(х). Это значит, что если g'0 (х) является верхней (нижней) оценкой для g 0(х), то g'n (х), вычисленная по прямой рекурсии (9), (10), также будет оценкой сверху (снизу) для gn (х). Оценка ошибки аппроксимации g 0( х) (22) может быть получена как разность текущего дробно-полиномиального разложения и разложения на единицу (или любое нечетное число) более высокого

п°рядка £п = ^ + Ьп / Кп | < |- Ьп+1 / Кп+1 + Ьп / Кп .

Здесь для краткости обозначено Ьп = Ьп (х),

Кп = Кп (х), gn = gn (х), £п = £п (х).

Рассмотрим выражение под знаком модуля

- Ьп+1 / Кп+1 + Ьп / Кп = (ЬпКп+1 - Ьп+Кп )/(Кп+1Кп). Для его упрощения преобразуем числитель отношения с учетом выполнения рекуррентного соотношения (9) для полиномов Ьп , Кп:

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Мп = ЬпКп+1 - КК =

= Ьп (Кп+1 - 2 хКп) - Ьп+1 (Кп - 2 хЬп) =

= 2п(ЬпКя-1 - Ьп-К) = -2пМп-1.

То есть Мп = (-2)пп! М0, М0 = Ь0 К1 - Ь1К0 = 1 или МП = (-2)пп!. Тогда выражение для ошибки примет вид

£п = ^0 + Ьп /Кп\< 2ппМ(Кп+К) < 2пп!/К\ . Функция gn также может быть представлена в виде дробно-полиномиального разложения. Из (12) и (13) следует представление gn в виде предела последовательности дробно-полиномиальных выражений

gn = Ит ^пт / Кт , (23)

где

^пт КтЬп КпЬт . (24)

Члены последовательности (КтЬп - КпЬт)/Кт являются верхними (при четном т ) и нижними

оценками (при нечетном т) для gn в соответствии со следствиями 1, 2 к лемме 3.

Заметим что ЬпКп+1 - Ья+К = (-2)пп!,

ЬяКя+2 - Ья+2Кп = Ьп (Кп+2 - 2(П + 1)Кп ) -

- Кп (Ьп+2 - 2(П + 1)Ьп ) =

= 2х(ЬяКя+1 - КпЬп+1) = 2х(-2)яп!.

При п Ф т величина Ь К - Ь К является полип т т п

номом степени |п - т\ -1. Использование данного полинома позволяет избежать ошибок вычитания, возникающих при непосредственном использовании прямой рекурсии. Для представления этого полинома можно использовать многочленную рекурсию.

Обозначим 2пт = ЬпКт - ЬтКп =

= Ьп (Кт - (2 хГПКп) - К я (Ьт - (2 хГПЬя) .

Из разности

т

Кт -(2х)т-яКя = X(К -2хК-х)(2х)т- =

= X 2(1 -1)(2 х )т-К-

непосредственно следует рекурсия

2ят = X 20' -1)(2 х)т-‘ 2п

(25)

2ПП = 0, 2яя+1 = (-2)пп! 2йХ. (26)

Коэффициенты полинома Znm являются целыми числами и легко могут быть рассчитаны без потери точности при применении к ним рекурсии (25), (26). В результате, Znm / Кт при т > п задает верхние и нижние границы для gn.

Для полиномов Znm может быть определена и другая рекурсия. Для её получения определим минимальное число р, не превышающее п и той же четности, что и т:

р = п + той(т -п,2). (27)

Преобразуем:

т т

Кт -ЕтрКр = £ Вт‘(К - 2(1 - 1)К-2) = £2хБт‘К-1,

'=р+2 '=р+2

где Вт,' = 2т)/2(т-1)\\/(1 -1)!!. (28)

Тогда с учетом равенства

2пт = Ьп (Кт - ВтРКр ) - Кп (Ьт - ВтРЬр ) + В^

получаем вторую рекуррентную формулу:

т

2 т = Втр2пр + X 2 хВт2п-1. (29)

,=р+2

При расчете, а также при исследовании свойств gn(х) полезно использовать отношение последовательных функций gn:

Гп (х) = gn (х)/gn-1( х) < 0 , П > 1 .

Из формулы (10) для этой функции при п > 1 можно получить выражения для прямой и обратной рекурсий:

Гя (х) = 2х +

Гп-1( х) =

2(п -1):

Гп-1( х)

2(п -1)

Гп(х) - 2 х

(30)

(31)

Соотношению (30) соответствует конечная непрерывная дробь

2(п -1) 2 (п - 2)

2

гя (х) = 2 х +^——...-------- -----. (32)

2 х + 2 х + 2 х +1/ g0(x)

Уравнению (31) соответствует бесконечная

непрерывная дробь

гп(х) = -2п- 2(п+) 2(я+2)... = -К(2к/2х). (33) 2 х + 2 х + 2 х + к=я

При этом прямая рекурсия (30) и обратная (31) должны быть устойчивы для х < 0 и х > 0 соответственно.

Из формулы (30) и гп (х) < 0 следует п -1

2х <\гп(х)\ < 2х +:

(х < 0, п > 2). (34)

Аналогично, из формулы (31) получаем

пп

-----П---< гп (х) < я , (х > 0, п > 1). (35)

х + (п +1)/2 х ь 1 х

Из этих свойств, в частности, находим

к(х)\ < ^1(х)

п!

(х > 0, п > 1),

х

к(х)\ < ^1(х)(2х'Г2 , (х < 0, п > 2).

Обобщением свойств (34), (35) являются следующие неравенства:

21 21+1

- К(2(к + п)/2х) < гп(х) < - К (2(к + )/2х) к=0 к=0

при х > 0, п > 1 и

21+1 21

2х + К(2(п-к)/2х) < гп(х) < 2х + К(2(п-к)/2х)

к=1 к=1

при х > 0, п > 21 + 2, для левой части неравенства и п > 21 + 1 для правой (I = 0,1,...).

Приведем еще одно соотношение, вытекающее из (10):

gn (х) / gn-2 (х) = 2(п - 1) + 2хгп-1 (х) ■ (36)

1. При х > 0 с учетом гп-1(х) < 0 и

gn (х)/ gn-2 (х) > 0 получим

0 < gn (х)/ gn-2( х) ^ 2(п -1).

2. При х < 0 с учетом гп-1(х) < 0 и (34) из (36) следует

0 < gn (x)/gn-2(x) < 4(п - ^ (п > 2) . (37)

Очевидно, последнее неравенство справедливо и для х > 0.

При разработке программы вычисления gn( х) мы ограничились п < 20. Численные расчеты

х

функций gn(х) в широком диапазоне х показывают, что прямая рекурсия (10) неустойчива (см. следствие 1). В этом случае можно воспользоваться обратной рекурсией, которая устойчива, но не является уточняющей. Она требует знания точных значений функций g 19( х), g 20( х). В нашей работе мы использовали для их вычисления на интервале 0 < х < 6 интерполяцию кубическими сплайнами и Паде-аппроксимацию. Для х > 6 мы использовали разложение в непрерывную дробь [13] отношения двух последовательных функций (33), обратную рекурсию для гп (х) (31) и рациональную Чебышевскую аппроксимацию для g0(х) [14]. При х < 0 применялась прямая рекурсия для функций гп (х) в сочетании со сплайн-интерполяцией для g 0( х) и g 1( х). Было найдено, что для больших отрицательных х рациональнее пользоваться выражением функций gn (х) через их значения в положительной области

gя (-х) = ехр(х2) | dz ехр(-z2) =

= ехр(х2) I | dz ехр(-z2) - |dz ехр(-z2)

V -ж х

= у/лехр(х2) - g0(х),

( dя Л

gn(-х) = (-1)П 'Пехр(х2) -gя (х)

(38)

Соотношение симметрии (38) позволяет получить более точные значения по сравнению с асимптотической формулой (11) и дает возможность избежать различного представления g-функции на отрезках полуоси (-да, 0].

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Разработанный алгоритм обеспечивал 9 верных знаков для требуемых вспомогательных функций при всех х. В действительности точность расчета зависит только от точности используемых аппроксимаций стартовых функций в уравнениях (9), (10). Заметим, что для важного случая Ь = с (х = у) в выражениях (5), (6) для интегралов (1) имеется неопределенность типа 0/0. Применяя правило Лопиталя для раскрытия неопределенности, можно найти дополнительные формулы в данном специальном случае. Однако более общим и эффективным оказывается алгоритм, основанный на разложении интеграла (6) в ряд Тейлора вблизи особой точки Ь = с.

Вычислим производную п-го порядка от выражения в квадратных скобках формулы (6) по переменной г = х - у :

(П) П

=1с

£0( У) - 8о( х)

7=0

[8о(У) - 8о(х)]

(П-7)

= I с

7=0

■ (-1)7 1'!

J [(-^Г^я-7 (У) - £п-7 (х)]

Г 1Л я-7 Ч2У

Используя правило Лопиталя, определим производную в точке г = 0:

£о(У)- £о(х)

(П)

.у(--7!

П (7 + 1)!

_ (-\)П+1 -1

2П+1 £п+1

(- 1Г'£п+1^-£«(|

ХС7 (-Ц = (-\)П+1 - 1 _ ^ П 7 +1 2П+1(п + 1ГП+

2

Откуда получим

£о( У)- £о( х) х - У

= у (-1) Гх - у

^0(27 +1)! I 2

2 7

8.2 7+1

х + У 2

Тогда разложение интеграла (6) в ряд Тейлора можно записать в виде

/(-1-1-1) =

4 л~

(-1)

ос а(х -0^(2;' +1)! ^ 2

-Е-

<§2 у+1

(39)

где ^ = х + у и ^ - I = (Ь + е)/(2^~а) > 0 согласно условию (3).

Этот ряд знакопостоянен по свойству (8) и уже не имеет особенности. Члены функционального ряда (39) неотрицательны. Для отношения соседних членов с использованием (37) можно получить следующее неравенство:

Нт g2 7+1( */2) (г/2)2

^ж g27-1(* /2) 27(27 +1)

.< „т-7*! = 0< 1.

'2.1(2 7 +1)

Следовательно, по признаку Даламбера ряд (39) сходится. Если при его вычислении ограничиться к -1 членом, то сумму остаточного ряда можно оценить с помощью системы неравенств. Перепишем (37) в виде (х)| < 4(х2 + п -1)|^„_2 (х)|. Выберем целое четное число / > х", тогда

\ёп (х) < 2"-,п (I + п -1)!! /(/ + т -1)!! \ёт (х)| , где т = то&(я,2). Учитывая (8) и формулу Стирлинга для оценки факториала, получаем

К = у-г^К/гУ л2 7 <^а(/2) у 22 7 у+2 /2

к у (27+1)! 1,2) (I + т -1)!!у (27 +1)!! 12

2 7

П

г

г=0

—х

<

2j er2

< g2j+1K /І) T 2 J (і + iJ)1 '2 f rJ 2 < - gl (s/2) T

< (і + m-1)!! У J! f 2J (і + m-1)!!у ^ 2J

< -giK/2') (er(l + -І-

(і + m -1)!! ( 2k

\ J У

j-l /2

I

j=k

er

2j

Для k > er / 2 получаем окончательную

оценку

rJ^!2- і er2 il + 1

(ег2/2к )к-1 /2 (I + т -1)!! ^ 2к)) 1 - (ег2/2к) '

При дифференцировании (39) по параметрам в, Ь, с (см. 5) получается разложение для соответствующего интеграла (1).

В таблице представлены результаты вычисления вспомогательных функций gn(x) по программе, реализованной на языке Fortran.

Тестирование всех интегралов было проведено при помощи системы аналитических вычислений Маріє [11].

Значения функции g 0( х ) и ее производных gn (х ) при различных значениях аргумента

n \ x -0.5 0.5 3.0 10.0

0 1.730234434E+0 5.456413608E-1 1.586356399E-1 4.975365939E-2

1 -2.730234434E+0 -4.543586392E-1 -4.818616082E-2 -4.926812176E-3

2 6.190703301E+0 6.369240823E-1 2.815431483E-2 9.710752718E-4

3 -1.711 164104E+1 -1.180510475E+0 -2.381875428E-2 -2.857432653E-4

4 5.425586084E+1 2.641034019E+0 2.601336334E-2 1.115863254E-4

5 -1.911489891E+2 -6.803049778E+0 -3.446985417E-2 -5.421961426E-5

б 7.337075976E+2 1.960729041E+1 5.331450836E-2 3.147096890E-5

10 3.011285623E+5 3.044513508E+3 9.302643509E+1 1.329000489E-5

15 -1.517837837E+9 -5.745251213E+6 -1.782099170E+2 -3.489317344E-5

20 1.657950330E+13 2.725817937E+10 1.167747294E+5 4.277646952E-4

Заключение

Таким образом, в настоящей работе получено аналитическое выражение двухчастичного шестикратного интеграла, встречающегося в высокоточных вариационных расчетах трехчастичных кулоновских систем. Разработан устойчивый и быстрый алгоритм для вычисления производных от функции Фаддеевой чисто мнимого аргумента до двадцатого порядка.

Благодарности

Считаем своим приятным долгом выразить благодарность кандидату физико-математических наук Д. А. Шершакову за многочисленные обсуждения, кандидату физико-математических наук В. А. Дубровскому за содействие при выполнении работы.

Список литературы

1. Shershakov D. A., Nechaev V. V., Berezin V. I. Exponential basis functions with quadratic dependence on interelectron distance for variational calculations of two-electron atoms // J. Phys. B. 2000. Vol. 33, № 1. P. 123-130.

2. Calais J.-L., Lowdin P. O. A simple method of treating atomic integrals containing function of r // J. Mol. Spectr. 19б2. Vol. 8, № 3. P. 203-211.

3. Pekeris C. L. Ground state of two-electron atoms // Phys. Rev. 1958. Vol. 112, № 5. P. 1б49-1б58.

4. SackR. A., Roothan C. C. J., Kolos W. Recursive evalu-tion of some atomic integrals // J. Math. Phys. 19б7. Vol. 8, № 5. P. 1093-1094.

5. Эфрос В. Д. Задача трех тел. Обобщенное экспоненциальное разложение, произвольные состояния в коррелированном базисе и энергия связи мезомо-лекул // Журн. эксперим. и теорет. физ. 198б. Т. 90, № 1. С. 10-24.

6. Ley-Koo E., Bunge C. F., Jauregui R. Evalution of relativistic atomic integrals using perimetric coordinates // Intern. J. Quant. Chem. 1997. Vol. б3, № 1. P. 93-97.

7. ФаддееваВ. Н., ТерентьевН.М. Таблицы значений функции w(z) = ez (1 + 2i / yfnj et dt) от комплексного аргумента. М., 1954. 2б8 с.

8. Абрамовиц М., Стиган И. Справочник по специальным функциям. М., 1979. 832 с.

9. Gautschi W. Efficient computation of the complex error function // SIAM J. Numer. Anal. 1970. Vol. 7, № 1. P. 187-198.

Ю. Н. Зайко. Динамика дуальной фазы

10. Poppe G. P. M., Wijers C. M. More efficient computation of the complex error function // ACM Trans. Math. Soft. 1990. Vol. 16, № 1. P. 38-46.

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

11. Maplesoft, «Maple», Version 15, Waterloo Maple Inc. (2012). URL: http://www.maplesoft.com (дата обращения: 01.07.2012).

12. Wolfram Research, Inc., «Mathematica», Version 8.0,

Champaign, IL (2012). URL: http://www.wolfram.com (дата обращения: 01.07.2012).

13. Джоунс У., Трон В. Непрерывные дроби / пер. с англ. М., 1985. 414 с.

14. Cody W. J. Rational Chebyshev approximations for the error function // Math. Comput. 1969. Vol. 23, № 107. P. 631-637.

УДК 537.8

ДИНАМИКА ДУАЛЬНОЙ ФАЗЫ

Ю. Н. Зайко

Поволжский институт управления им. П. А. Столыпина филиал РАНХ и ГС при Президенте РФ, Саратов E-mail: zyrnick@rambler.ru

В статье рассмотрена динамика дуальной фазы, связанной с дуальным характером дионов - гипотетических частиц, обладающих одновременно электрическим и магнитным зарядом. Получены уравнения, описывающие динамику дуальной фазы из условия, чтобы уравнения Максвелла сохраняли свой вид под действием дуальных преобразований, соответствующий эффективному электрическому характеру заряда диона. Найдены частные решения этих уравнений в вакууме и в поле электрического диполя. Показано, что найденные ранее решения уравнений Максвелла в виде сферических волн с нулевым орбитальным моментом соответствуют монопольному излучению диона. Показан голдстоуновский характер дуальной фазы.

Ключевые слова: дион, магнитный заряд, голдстоуновская мода.

Dynamics of Dual Phase Yu. N. Zayko

This article presents dynamics of dual phase, which is connected with dual character of dions - hypothetical particles possessing both electrical and magnetic charges. For description of dual phase the set of equations is received from condition that Maxwell equations conserve their electrical character under dual transformations, what corresponds to effective electrical charge of dion. Some special solutions of these equations for example in vacuum, and in field of electrical dipole are found. It was shown, that solutions of Maxwell equations for spherical waves with zero orbital moment momentum correspond to monopole radiation of dion. Goldstone character of dual phase was shown.

Key words: dion, magnetic charge, Goldstone mode.

Введение

Уравнения Максвелла допускают формулировку, симметричную относительно электрических и магнитных зарядов [1]:

rotH----------= eJ, divE = e p ,

dt

- dH - -

rotE Л--------=-g.J, divH = gp.

dt

(1)

Здесь Е и Н - напряженности электрического и магнитного полей, р и 3 - плотность числа и тока частиц соответственно (в системе Хевисайда-Лоренца). Предполагается, что частицы, называемые дионами и служащие источниками полей, несут одновременно электрический (е) и магнитный (^ ) заряды. Как известно, с помощью дуального преобразования полей

Ё = cos в ■ Ё + sin в ■ Й,

H = -sin 9-Ё + cos в-Й,

О - е а - 8

(2)

'=4

Є2 + g2.

Ч 4

где в - параметр дуального преобразования (дуальная фаза), можно исключить магнитный заряд из уравнений (1), т.е. привести их к виду [1]:

. й д E r rot H-------------= qJ,

дt

д H j

rot E + —— = - gJ,

д t

div E = qp.

div H = 0,

(3)

описывающему электродинамику частиц с эффективным электрическим зарядом q1. Таким образом, вводимый подобным способом магнитный заряд оказывается ненаблюдаемым [1] 2. При этом параметр в считается постоянным, хотя и неопределенным в силу неопределенности е и g. С этой точки зрения как заряды (электрический и магнитный), так и поля получают свое наименование в результате соглашения, и мы всегда можем изменить его и называть электрическое

1 Аналогично уравнения движения дуально заряженных частиц приводятся к уравнениям движения частиц, несущих эффективный электрический заряд q.

2 Аналогично может быть «исключен» и электрический заряд [1].

© Зайко Ю. Н., 2012

i Надоели баннеры? Вы всегда можете отключить рекламу.