Физика плазмы, 2021, T. 47, № 6, стр. 483-498

О восстановлении двумерной функции распределения ионов в пробкотроне по измерениям спектров коллективного томсоновского рассеяния

Е. Д. Господчиков ab, Т. А. Хусаинов a*, А. Г. Шалашов a

a Институт прикладной физики РАН
Нижний Новгород, Россия

b Нижегородский государственный университет им. Н.И. Лобачевского
Нижний Новгород, Россия

* E-mail: hta@ipfran.ru

Поступила в редакцию 07.12.2020
После доработки 20.01.2021
Принята к публикации 25.01.2021

Полный текст (PDF)

Аннотация

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

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

1. ВВЕДЕНИЕ

Рассеяние мощного зондирующего миллиметрового излучения на тепловых флуктуациях магнитоактивной плазмы (так называемое коллективное томсоновское рассеяние) позволяет получать информацию о распределении ионов плазмы по скоростям с хорошим пространственным и временным разрешением [1]. Возможности анализа спектров рассеянного сигнала для диагностики функции распределения быстрых ионов в тороидальных магнитных ловушках были продемонстрированы экспериментально на токамаках TFTR [2], JET [3], TEXTOR [47], ASDEX-U [6, 8], на стеллараторе W7AS для диагностики температуры тепловых ионов и нижнегибридной неустойчивости плазмы [911], на стеллараторе LHD для диагностики распределения как быстрых, так и тепловых ионов [12] и, совсем недавно, на новейшем стеллараторе W7X для диагностики температуры тепловых ионов [13]. Наряду с оптическими методами, спектроскопией нейтронов и гамма-квантов, коллективное рассеяние миллиметрового излучения является одним из основных способов диагностики функции распределения быстрых ионов в токамаках [14]; в частности, этот метод рассматривается как основной способ диагностики термоядерных альфа-частиц в токамаке-реакторе ИТЭР [15].

Последний прогресс в методах удержания термоядерной плазмы в существующих прямых магнитных ловушках [1620] и планирование новых физических задач для прямых магнитных ловушек следующего поколения [20, 21] привели к тому, что сейчас активно совершенствуются методы диагностики плазмы в прямых магнитных ловушках. Важной возможностью, позволяющей проводить диагностику функции распределения энергичных ионов, является реализация в прямой магнитной ловушке опробованного на тороидальных магнитных ловушках метода регистрации спектров коллективного томсоновского рассеяния [22].

Следуя работе [22], отметим основные особенности коллективного томсоновского рассеяния в плазме крупномасштабной прямой магнитной ловушки на примере крупнейшей на данный момент подобной установки ГДЛ (газодинамическая ловушка), функционирующей в ИЯФ им. Г.И. Будкера СО РАН. Плазма в этой установке делится на две фракции: фоновую плазму, удерживаемую в газодинамическом режиме, и высокоэнергичные ионы, появляющиеся в результате бомбардировки фоновой плазмы пучками нейтральных атомов с энергией порядка 20–30 кэВ [23]. Энергичные ионы, во-первых, составляют существенную часть всех ионов плазмы (20% в центральном сечении и до 100% в точках разворота); во-вторых, они удерживаются в пробкотроне в бесстолкновительном адиабатическом режиме. Частота соударений быстрых ионов с электронами и ионами фоновой плазмы много меньше частоты баунс-осцилляций между точками разворота в магнитном поле. В отличие от тороидальных систем модуль магнитного поля в открытой ловушке меняется сильно (например, в установке ГДЛ пробочное отношение превышает 30), поэтому функция распределения энергичных ионов будет существенно разной в разных поперечных сечениях ловушки [22]. Это, с одной стороны, накладывает дополнительные ограничения на пространственное разрешение метода диагностики, а с другой стороны, может быть использовано для восстановления функции распределения частиц.

Напомним, что диагностика плазмы по коллективному рассеянию заключается в регистрации спектральной плотности мощности рассеянного излучения $d{{P}_{s}}{\text{/}}d\omega $, которая в свою очередь пропорциональна спектральной плотности флуктуаций плотности плазмы ${{S}_{f}}({\mathbf{k}},\omega )$, где $\omega = {{\omega }_{{\text{s}}}} - {{\omega }_{{\text{i}}}}$ и ${\mathbf{k}} = {{{\mathbf{k}}}_{{\text{s}}}} - {{{\mathbf{k}}}_{{\text{i}}}}$ – разность частот и волновых векторов зондирующего и рассеянного излучения. В эксперименте обычно направление k фиксировано геометрией рассеяния, а ω меняется (попадает в полосу спектроанализатора). Поэтому мы имеем одномерный спектр рассеяния, который непосредственно позволяет восстановить одномерную функцию распределения ионов по проекциям скорости на направление k

${{S}_{f}} \propto F(\omega {\text{/}}k),\quad F(u) = \int {f({\mathbf{v}})\delta (u - {\mathbf{vk}}{\text{/}}k){{d}^{3}}{\mathbf{v}}} ,$
где $f({\mathbf{v}})$ – локальная функция распределения ионов или электронов по трехмерному пространству скоростей, рассчитанная в рассеивающем объеме. Для предложенного в [22] эксперимента на установке ГДЛ угол рассеяния и полоса анализа по ω выбираются так, что спектральная плотность флуктуаций определяется одномерной функцией распределения быстрых ионов.

В дальнейшем будем предполагать, что для вращающейся в магнитном поле частицы все направления скорости поперек магнитного поля равновероятны, поэтому ансамбль таких частиц можно характеризовать двумерной функцией распределения $f({{\upsilon }_{{||}}},{{\upsilon }_{ \bot }})$ по поперечной и продольной по отношению к магнитному полю компонентам вектора скорости. Воспользовавшись этой естественной симметрией, мы можем переписать выражение для одномерной функции распределения в виде

(1)
$\begin{gathered} F\left( {u,\theta } \right) = \int\limits_0^{2\pi } {\int\limits_{ - \infty }^\infty {\int\limits_0^\infty {f\left( {{{\upsilon }_{{||}}},{{\upsilon }_{ \bot }}} \right)} } } \times \\ \, \times \delta \left( {u - {{\upsilon }_{{||}}}\cos \theta - {{\upsilon }_{ \bot }}\sin \theta \cos \varphi } \right){{\upsilon }_{ \bot }}d{{\upsilon }_{ \bot }}d{{\upsilon }_{{||}}}d\varphi , \\ \end{gathered} $
где θ – угол между вектором k и направлением магнитного поля. Измерение спектра рассеянного сигнала при заданном угле рассеяния не позволяет восстановить полную функцию распределения $f\left( {{{\upsilon }_{{||}}},{{\upsilon }_{ \bot }}} \right)$, однако такое восстановление в принципе возможно по двумерной функции $F(u,\theta )$, то есть если мы проведем измерения при некотором наборе углов рассеяния. Такой подход к восстановлению функции распределения близок к классической задаче томографии, которая в данном случае может быть решена стандартными математическими методами решения некорректных обратных задач. Метод был предложен и теоретически проанализирован [24], продемонстрирован на синтетических данных [25] и успешно реализован при обработке реальных экспериментальных данных [14].

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

2. ИНТЕГРАЛЬНОЕ УРАВНЕНИЕ

Рассмотрим формулу (1) для произвольной точки вдоль силовой линии магнитного поля:

$\begin{gathered} F\left( {u,R,\theta } \right) = \int\limits_0^{2\pi } {\int\limits_{ - \infty }^\infty {\int\limits_0^\infty {{{f}_{R}}\left( {\upsilon _{{||}}^{'},\upsilon _{ \bot }^{'}} \right) \times } } } \\ \, \times \delta \left( {u - \upsilon _{{||}}^{'}\cos \theta - \upsilon _{ \bot }^{'}\sin \theta \cos \varphi {\kern 1pt} '} \right)\upsilon _{ \bot }^{'}d\upsilon _{ \bot }^{'}d\upsilon _{{||}}^{'}d\varphi {\kern 1pt} ', \\ \end{gathered} $
где ${{f}_{R}}(\upsilon _{{||}}^{'},\upsilon _{ \bot }^{'})$ – функция распределения в сечении, где напряженность магнитного поля равна $B = R \cdot {{B}_{{\min }}}$, локальное пробочное отношение R играет роль координаты вдоль силовой линии магнитного поля. Здесь и далее мы будем использовать обозначения $\upsilon _{{||, \bot }}^{'}$ для скорости в заданном сечении R и ${{\upsilon }_{{{\text{||,}} \bot }}}$ – для скорости в центральном сечении $R = 1$, отвечающем минимальному значению магнитного поля. Для иона, осуществляющего бесстолкновительное движение в плавнонеоднородном магнитном поле, можно ввести два интеграла движения

$\begin{gathered} \upsilon _{ \bot }^{{'2}}{\text{/}}R = {\text{const, }} \\ \upsilon _{{{\text{||}}}}^{{'{\text{2}}}} + \upsilon _{ \bot }^{{'{\text{2}}}} + 2e{{Z}_{{\text{i}}}}{{m}_{{\text{i}}}}\Phi {\kern 1pt} ' = {\text{const}}{\text{.}} \\ \end{gathered} $

Первый из них отвечает сохранению магнитного момента, второй – сохранению энергии с учетом электростатического амбиполярного потенциала $\Phi {\kern 1pt} '$. Характерные значения перепада амбиполярного потенциала вдоль силовой линии определяются температурой мишенной плазмы [18, 26], поэтому при рассмотрении динамики энергичных ионов им можно пренебречь. В результате для быстрых ионов связь между скоростями ${{\upsilon }_{{{\text{||,}} \bot }}}$ и $\upsilon _{{{\text{||,}} \bot }}^{'}$ можно найти в явном виде

(2)
$\begin{gathered} \upsilon _{ \bot }^{'} = \sqrt R {{\upsilon }_{ \bot }}, \\ \upsilon _{{||}}^{'} = {{\upsilon }_{{||}}}\;\sqrt {1 - \left( {R - 1} \right)\upsilon _{ \bot }^{2}{\text{/}}\upsilon _{{||}}^{2}} . \\ \end{gathered} $

В силу теоремы Лиувилля фазовый объем сохраняется, поэтому функция распределения ионов в данном сечении R выражается через функцию распределения в центральном сечении $R = 1$ ловушки следующим образом:

${{f}_{R}}\left( {\upsilon _{{||}}^{'},\upsilon _{ \bot }^{'}} \right) = {{f}_{1}}\left( {{{\upsilon }_{{||}}},{{\upsilon }_{ \bot }}} \right),$
где скорости в минимуме магнитного поля выражены через скорости в сечении R в соответствии с преобразованиями (2).

Для простоты выберем продольную геометрию рассеяния $\theta = 0$ и возьмем интеграл по $\varphi {\kern 1pt} '$. Получим

$F\left( {u,R,0} \right) = 2\pi \int\limits_{ - \infty }^\infty {\int\limits_0^\infty {{{f}_{R}}\left( {\upsilon _{{||}}^{'},\upsilon _{ \bot }^{'}} \right)\delta \left( {u - \upsilon _{{||}}^{'}} \right)\upsilon _{ \bot }^{'}d\upsilon _{ \bot }^{'}d\upsilon _{{||}}^{'}} .} $

Далее перейдем в подынтегральном выражении к скоростям в минимуме магнитного поля

$\begin{gathered} F\left( {u,R} \right) = \int\limits_{\upsilon _{{||}}^{2} \geqslant \left( {R - 1} \right)\upsilon _{ \bot }^{2}} {{{f}_{1}}\left( {{{\upsilon }_{{||}}},{{\upsilon }_{ \bot }}} \right) \times } \\ \, \times \delta \left( {u - {{\upsilon }_{{||}}}\sqrt {1 - {{\left( {R - 1} \right)\upsilon _{ \bot }^{2}} \mathord{\left/ {\vphantom {{\left( {R - 1} \right)\upsilon _{ \bot }^{2}} {\upsilon _{{||}}^{2}}}} \right. \kern-0em} {\upsilon _{{||}}^{2}}}} } \right)\frac{{2\pi R{{\upsilon }_{ \bot }}d{{\upsilon }_{ \bot }}d{{\upsilon }_{{||}}}}}{{\sqrt {1 - {{\left( {R - 1} \right)\upsilon _{ \bot }^{2}} \mathord{\left/ {\vphantom {{\left( {R - 1} \right)\upsilon _{ \bot }^{2}} {\upsilon _{{||}}^{2}}}} \right. \kern-0em} {\upsilon _{{||}}^{2}}}} }}. \\ \end{gathered} $

Формально интеграл берется по области $\upsilon _{{||}}^{2} \geqslant \left( {R - 1} \right)\upsilon _{ \bot }^{2}$, в которую переходит полуплоскость $\upsilon _{{||}}^{'} \in \Re $, $\upsilon _{ \bot }^{'} \geqslant 0$. Однако за пределами этой области δ-функция всюду равна нулю, поэтому мы можем расширить интеграл на все пространство скоростей ${{\upsilon }_{{||}}},\;{{\upsilon }_{ \bot }}$. Также в силу ожидаемой в открытой ловушке симметричности функции распределения для захваченных ионов по ${{\upsilon }_{{||}}}$, мы можем рассматривать только область $u \geqslant 0$ для $F\left( {u,R} \right)$; при этом подынтегральная δ-функция будет равна нулю в области ${{\upsilon }_{{||}}} < 0$. С учетом этого мы можем переписать интеграл в более простом виде

$\begin{gathered} F\left( {u,R} \right) = \int\limits_0^\infty {\int\limits_0^\infty {{{f}_{1}}\left( {{{\upsilon }_{{||}}},{{\upsilon }_{ \bot }}} \right)\delta \left( {u - {{\upsilon }_{{||}}}\sqrt {1 - \left( {R - 1} \right)\upsilon _{ \bot }^{2}{\text{/}}\upsilon _{{||}}^{2}} } \right) \times } } \\ \, \times \frac{{2\pi R{{\upsilon }_{ \bot }}d{{\upsilon }_{ \bot }}d{{\upsilon }_{{||}}}}}{{\sqrt {1 - \left( {R - 1} \right)\upsilon _{ \bot }^{2}{\text{/}}\upsilon _{{||}}^{2}} }}. \\ \end{gathered} $

Возьмем интеграл по ${{\upsilon }_{{||}}}$, получим

(3)
$F\left( {u,R} \right) = 2\pi R\int\limits_0^\infty {{{f}_{1}}\left( {\sqrt {{{u}^{2}} + \left( {R - 1} \right)\upsilon _{ \bot }^{2}} ,{{\upsilon }_{ \bot }}} \right){{\upsilon }_{ \bot }}d{{\upsilon }_{ \bot }}} .$

Это основное соотношение, которое мы будем анализировать в предположении, что по результатам диагностики плазмы нам известна непрерывная функция $F\left( {u,R} \right)$. Конечно, при использовании данных реального эксперимента для построения такой функции придется воспользоваться интерполяцией по дискретному набору измерений.

Для компактности последующих выкладок введем новые обозначения

$f(x,y) = {{f}_{1}}\left( {{{\upsilon }_{0}}\sqrt x ,{{\upsilon }_{0}}\sqrt y } \right),\quad x = \upsilon _{{||}}^{2}{\text{/}}\upsilon _{{\text{0}}}^{{\text{2}}},\quad y = \upsilon _{ \bot }^{2}{\text{/}}\upsilon _{0}^{2},$
$\begin{gathered} G\left( {\xi ,w} \right) = {{(\pi R)}^{{ - 1}}}F\left( {{{\upsilon }_{0}}\sqrt \xi ,R} \right), \\ \xi = {{u}^{2}}{\text{/}}\upsilon _{0}^{2},\quad w = R - 1, \\ \end{gathered} $
где ${{\upsilon }_{0}}$ – любая характерная скорость, введенная только для нормировки. В итоге, мы приходим к интегральному уравнению

(4)
$G\left( {\xi ,w} \right) = \int\limits_0^\infty {f\left( {\xi + wy,y} \right)} dy.$

Математическая постановка задачи восстановления состоит в том, чтобы по заданной функции $G\left( {\xi ,w} \right)$ найти функцию $f\left( {x,\,y} \right)$. При этом будем считать, что переменные $\xi ,\;w,\;x,\;y$ изменяются в интервале $\left[ {0,\,\infty } \right)$. Функция $f\left( {x,\,y} \right)$ является действительной положительно определенной кусочно-гладкой функцией, интегралы от которой по любой конечной области ограничены. Отметим, что реальная функция распределения в силу релятивистских ограничений также должна обращаться в нуль начиная с определенных значений переменных $x,\;y$, однако для удобства вычисления аналитических интегралов мы не будем требовать выполнения этого условия, а вместо этого оставим за собой свободу вводить ограничения на закон спадания функции $f\left( {x,\,y} \right)$ при $x,\,y \to \infty $.

3. ПРЕОБРАЗОВАНИЕ ФУРЬЕ

Поскольку все аргументы в уравнении (4) заданы на положительной полуоси, для поиска решения естественно воспользоваться косинус-преобразованием Фурье. Определим

(5)
$\begin{gathered} {{F}_{c}}\left( {{v},y} \right) = \int\limits_0^\infty {f\left( {x,y} \right)\cos \left( {x{v}} \right)} dx, \\ f\left( {x,y} \right) = \frac{2}{\pi }\int\limits_0^\infty {{{F}_{c}}\left( {{v},y} \right)\cos \left( {x{v}} \right)} d{v}. \\ \end{gathered} $

Необходимым условием для этих равенств является сходимость интегралов $\int_0^\infty {f\left( {x,y} \right)} dx$ для всех значений $y \geqslant 0$. После преобразования (5) уравнение (4) перепишется в виде

(6)
$G\left( {\xi ,w} \right) = \frac{2}{\pi }\int\limits_0^\infty {\int\limits_0^\infty {{{F}_{c}}\left( {{v},y} \right)\cos \left( {{v}\xi + {v}wy} \right)} dyd{v}} .$

Мы видим, что косинус-преобразование Фурье сохраняет интегральный характер уравнения, однако приводит его к нормальному виду, допускающему стандартное решение. Это решение приведено в следующем разделе, но прежде чем переходить к нему, сделаем одно отступление методического характера.

Функции $G\left( {\xi ,\,w} \right)$ и $f\left( {x,\,y} \right)$ принципиально заданы нам только при положительных значениях своих параметров. Давайте на время забудем об этом и распространим все наши зависимости в область $\xi ,\,w,\,x \in \left( { - \infty ,\,\infty } \right)$, но по-прежнему будем считать $y \geqslant 0$. Тогда мы можем применить к уравнению (4) комплексное преобразование Фурье по переменным $\xi ,\;w,\;x$:

$\begin{gathered} \tilde {f}\left( {{{k}_{x}},y} \right) = \int\limits_{ - \infty }^\infty {f\left( {x,y} \right){{e}^{{i{{k}_{x}}x}}}} dx, \\ f\left( {x,y} \right) = \frac{1}{{2\pi }}\int\limits_{ - \infty }^\infty {\tilde {f}\left( {{{k}_{x}},y} \right){{e}^{{ - i{{k}_{x}}x}}}} d{{k}_{x}} \\ \end{gathered} $
$\begin{gathered} \tilde {G}\left( {{{k}_{\xi }},{{k}_{w}}} \right) = \int\limits_{ - \infty }^\infty {\int\limits_{ - \infty }^\infty {G\left( {\xi ,w} \right){{e}^{{i{{k}_{\xi }}\xi + i{{k}_{w}}w}}}d\xi dw} } , \\ G\left( {\xi ,w} \right) = \frac{1}{{4{{\pi }^{2}}}}\int\limits_{ - \infty - }^\infty {\int\limits_\infty ^\infty {\tilde {G}\left( {{{k}_{\xi }},{{k}_{w}}} \right){{e}^{{ - i{{k}_{\xi }}\xi - i{{k}_{w}}w}}}d{{k}_{\xi }}d{{k}_{w}}} } \\ \end{gathered} $

Уравнение (4) при этом примет вид

$\begin{gathered} \tilde {G}\left( {{{k}_{\xi }},{{k}_{w}}} \right) = \\ \, = \frac{1}{{2\pi }}\int\limits_{ - \infty }^\infty {\int\limits_{ - \infty }^\infty {\int\limits_{ - \infty }^\infty {\int\limits_0^\infty {\tilde {f}\left( {{{k}_{x}},y} \right){{e}^{{i{{k}_{\xi }}\xi + i{{k}_{w}}w - i{{k}_{x}}\left( {\xi + wy} \right)}}}dyd{{k}_{x}}d\xi dw} } } } = \\ \end{gathered} $
$\begin{gathered} \, = 2\pi \int\limits_{ - \infty }^\infty {\int\limits_0^\infty {\tilde {f}\left( {{{k}_{x}},y} \right)\delta \left( {{{k}_{\xi }} - {{k}_{x}}} \right)\delta \left( {{{k}_{w}} - {{k}_{x}}y} \right)dyd{{k}_{x}}} } = \\ \, = 2\pi \int\limits_0^\infty {\tilde {f}\left( {{{k}_{\xi }},y} \right)\delta \left( {{{k}_{w}} - {{k}_{\xi }}y} \right)dy} . \\ \end{gathered} $

Далее учитываем, что интеграл берется только для положительных значений y, и приходим к

$\tilde {G}\left( {{{k}_{\xi }},{{k}_{w}}} \right) = \left\{ \begin{gathered} {{2\pi \tilde {f}\left( {{{k}_{\xi }},{{{{k}_{w}}} \mathord{\left/ {\vphantom {{{{k}_{w}}} {{{k}_{\xi }}}}} \right. \kern-0em} {{{k}_{\xi }}}}} \right)} \mathord{\left/ {\vphantom {{2\pi \tilde {f}\left( {{{k}_{\xi }},{{{{k}_{w}}} \mathord{\left/ {\vphantom {{{{k}_{w}}} {{{k}_{\xi }}}}} \right. \kern-0em} {{{k}_{\xi }}}}} \right)} {\left| {{{k}_{\xi }}} \right|,\quad {{k}_{\xi }}{{k}_{w}} \geqslant 0}}} \right. \kern-0em} {\left| {{{k}_{\xi }}} \right|,\quad {{k}_{\xi }}{{k}_{w}} \geqslant 0}} \hfill \\ 0,\quad {{k}_{\xi }}{{k}_{w}} < 0 \hfill \\ \end{gathered} \right..$

При ${{k}_{\xi }}{{k}_{w}} \geqslant 0$ это соотношение дает нам простую алгебраическую связь между фурье-образами заданной функции G и искомой функции f. Однако при ${{k}_{\xi }}{{k}_{w}} < 0$ мы получили дополнительное ограничение на функцию $\tilde {G}$, которое необходимо учитывать при доопределении ее прообраза $G\left( {\xi ,w} \right)$ на область отрицательных аргументов. Чтобы понять суть этого ограничения, мы можем перенести неопределенность и ограничения функции G при $\xi ,w < 0$ в неопределенность и ограничения функции $\tilde {G}$ в области ${{k}_{\xi }}{{k}_{w}} \geqslant 0$. Воспользуемся действительностью $G\left( {\xi ,w} \right)$ и $f\left( {x,y} \right)$, то есть $\tilde {G}\left( {{{k}_{\xi }},{{k}_{w}}} \right) = \tilde {G}{\kern 1pt} *\left( { - {{k}_{\xi }}, - {{k}_{w}}} \right)$ и $\tilde {f}\left( {{{k}_{x}},y} \right) = $ $ = \tilde {f}{\kern 1pt} *\left( { - {{k}_{x}},y} \right)$, где звездочкой обозначена операция комплексного сопряжения. Без потери общности будем считать, что функция распределения доопределена в нефизическую область симметричным образом: $f\left( {x,y} \right) = f\left( { - x,y} \right)$. Тогда комплексное преобразование Фурье переходит в косинус-преобразование (6), $\tilde {f}\left( {{{k}_{x}},y} \right) = \tilde {f}\left( { - {{k}_{x}},y} \right)$ и $\tilde {G}\left( {{{k}_{\xi }},{{k}_{w}}} \right) = \tilde {G}\left( { - {{k}_{\xi }}, - {{k}_{w}}} \right)$. Используя эти условия, имеем

$\begin{gathered} G\left( {\xi ,w} \right) = \frac{1}{{4{{\pi }^{2}}}}\int\limits_0^\infty {\int\limits_0^\infty {\tilde {G}\left( {{{k}_{\xi }},{{k}_{w}}} \right) \times } } \\ \, \times \left( {{{e}^{{i{{k}_{\xi }}\xi + i{{k}_{w}}w}}} + {{e}^{{ - i{{k}_{\xi }}\xi - i{{k}_{w}}w}}}} \right)d{{k}_{\xi }}d{{k}_{w}} = \\ \, = \frac{1}{{2{{\pi }^{2}}}}\int\limits_0^\infty {\int\limits_0^\infty {\tilde {G}\left( {{{k}_{\xi }},{{k}_{w}}} \right)\cos \left( {{{k}_{\xi }}\xi + {{k}_{w}}w} \right)d{{k}_{\xi }}d{{k}_{w}}} } . \\ \end{gathered} $

Получившееся уравнение эквивалентно уравнению (6). Таким образом, если мы корректно определим функцию $G\left( {\xi ,\,w} \right)$ в области $\xi ,\,w < 0$, то с помощью комплексного преобразования Фурье мы можем элементарно решить уравнение (4). Но сама задача корректного доопределения $G\left( {\xi ,\,w} \right)$ все равно сводится к решению интегрального уравнения (6), порождаемого косинус-преобразованием Фурье.

4. РЕШЕНИЕ ЧЕРЕЗ ИНТЕГРАЛЬНОЕ ПРЕОБРАЗОВАНИЕ МЕЛЛИНА

Уравнение (6) сформулировано для неизвестной функции двух переменных. Несмотря на простой вид, нам не удалось найти это уравнение в математических справочниках. Тем не менее, его решение удается построить, опираясь на технику, развитую для решения одномерных линейных интегральных уравнений 1-го рода [27].

Введем интегральное преобразование Меллина [28]

$\hat {f}\left( s \right) = \int\limits_0^\infty {f\left( x \right){{x}^{{s - 1}}}dx} ,\quad f\left( x \right) = \frac{1}{{2\pi i}}\int\limits_{\sigma - i\infty }^{\sigma + i\infty } {\hat {f}\left( s \right){{x}^{{ - s}}}ds} .$

Интеграл во втором выражении берется в смысле главного значения, а параметр σ принадлежит интервалу $\sigma \in \left( {{{\sigma }_{1}},\,{{\sigma }_{2}}} \right)$, границы которого определяются сходимостью интегралов

$\int\limits_0^1 {f\left( x \right){{x}^{{{{\sigma }_{1}} - 1}}}dx} < \infty ,\quad \int\limits_1^\infty {f\left( x \right){{x}^{{{{\sigma }_{2}} - 1}}}dx} < \infty .$

Сходимость приведенных интегралов для некоторых ${{\sigma }_{1}},\;{{\sigma }_{2}}$ при ${{\sigma }_{1}} < \,{{\sigma }_{2}}$ является критерием применимости преобразования Меллина для функции $f\left( x \right)$.

Умножим левую и правую части уравнения (6) на ${{\xi }^{{{{s}_{\xi }} - 1}}}{{w}^{{{{s}_{w}} - 1}}}$ и проинтегрируем по ξ и w от нуля до бесконечности, получим

(7)
$\begin{gathered} \hat {G}\left( {{{s}_{\xi }},{{s}_{w}}} \right) = \frac{2}{\pi }\int\limits_0^\infty {\int\limits_0^\infty {{{F}_{c}}\left( {{v},y} \right)} } \times \\ \, \times \left[ {\int\limits_0^\infty {\int\limits_0^\infty {\cos \left( {{v}\xi + {v}wy} \right)} \,{{\xi }^{{{{s}_{\xi }} - 1}}}{{w}^{{{{s}_{w}} - 1}}}d\xi dw} } \right]dyd{v}. \\ \end{gathered} $

В левой части уравнения (7) стоит двумерное преобразование Меллина для функции $G\left( {\xi ,\,w} \right)$

(8)
$\hat {G}\left( {{{s}_{\xi }},{{s}_{w}}} \right) = \int\limits_0^\infty {\int\limits_0^\infty {G\left( {\xi ,w} \right){{\xi }^{{{{s}_{\xi }} - 1}}}{{w}^{{{{s}_{w}} - 1}}}d\xi dw} } ,$
обратное преобразование к которому также имеет вид

$G\left( {\xi ,w} \right) = - \frac{1}{{4{{\pi }^{2}}}}\int\limits_{{{\sigma }_{\xi }} - i\infty }^{{{\sigma }_{\xi }} + i\infty } {\int\limits_{{{\sigma }_{w}} - i\infty }^{{{\sigma }_{w}} + i\infty } {\hat {G}\left( {{{s}_{\xi }},{{s}_{w}}} \right){{\xi }^{{ - {{s}_{\xi }}}}}{{w}^{{ - {{s}_{w}}}}}d{{s}_{w}}d{{s}_{\xi }}} } .$

Диапазон изменения параметров ${{\sigma }_{\xi }} = \operatorname{Re} \left( {{{s}_{\xi }}} \right)$, ${{\sigma }_{w}} = \operatorname{Re} \left( {{{s}_{w}}} \right)$ мы определим после того, как закончим преобразование правой части уравнения (7). Выражение в квадратных скобках в (7) является преобразованием Меллина для ядра интегрального преобразования (6). Преобразуем это выражение, перейдя к новым переменным интегрирования

$\bar {\xi } = {v}\xi ,\quad \bar {w} = {v}yw.$

Теперь выражение в квадратных скобках можно посчитать, используя табличные преобразования Меллина для тригонометрических функций [28]:

$\left[ {...} \right] = \frac{1}{{{{{v}}^{{{{s}_{\xi }} + {{s}_{w}}}}}{{y}^{{{{s}_{w}}}}}}}\int\limits_0^\infty {\int\limits_0^\infty {\cos \left( {\bar {\xi } + \bar {w}} \right)} {{{\bar {\xi }}}^{{{{s}_{\xi }} - 1}}}{{{\bar {w}}}^{{{{s}_{w}} - 1}}}d\bar {\xi }d\bar {w}} = $
(9)
$\begin{gathered} \, = \frac{1}{{{{{v}}^{{{{s}_{\xi }} + {{s}_{w}}}}}{{y}^{{{{s}_{w}}}}}}}\int\limits_0^\infty {\int\limits_0^\infty {\left[ {\cos \left( {\bar {\xi }} \right)\cos \left( {\bar {w}} \right) - \sin \left( {\bar {\xi }} \right)\sin \left( {\bar {w}} \right)} \right]} \times } \\ \, \times {{{\bar {\xi }}}^{{{{s}_{\xi }} - 1}}}{{{\bar {w}}}^{{{{s}_{w}} - 1}}}d\bar {\xi }d\bar {w} = \frac{1}{{{{{v}}^{{{{s}_{\xi }} + {{s}_{w}}}}}{{y}^{{{{s}_{w}}}}}}}\Gamma \left( {{{s}_{\xi }}} \right)\Gamma \left( {{{s}_{w}}} \right) \times \\ \end{gathered} $
$\, \times \cos \left( {\frac{\pi }{2}\left( {{{s}_{\xi }} + {{s}_{w}}} \right)} \right),$
где $\Gamma (...)$ обозначает гамма-функцию, а выкладки корректны при $\operatorname{Re} \left( {{{s}_{\xi }}} \right),$ $\operatorname{Re} \left( {{{s}_{w}}} \right) \in \left( {0,1} \right)$. После проделанных преобразований выражение (7) принимает вид

$\begin{gathered} \hat {G}\left( {{{s}_{\xi }},{{s}_{w}}} \right) = \frac{2}{\pi }\Gamma \left( {{{s}_{\xi }}} \right)\Gamma \left( {{{s}_{w}}} \right)\cos \left( {\frac{\pi }{2}\left( {{{s}_{\xi }} + {{s}_{w}}} \right)} \right) \times \\ \, \times \int\limits_0^\infty {\int\limits_0^\infty {{{F}_{c}}\left( {{v},y} \right)} } {{{v}}^{{ - {{s}_{\xi }} - {{s}_{w}}}}}{{y}^{{ - {{s}_{w}}}}}dyd{v}. \\ \end{gathered} $

Оставшиеся интегралы представляют преобразование Меллина для фурье-образа функции распределения с аргументами

${{s}_{y}} = 1 - {{s}_{w}},\quad {{s}_{{v}}} = 1 - {{s}_{\xi }} - {{s}_{w}}.$

Таким образом, после преобразований Меллина изначальное интегральное уравнение представляет собой просто алгебраическое произведение образов с заменой переменных

(10)
$\begin{gathered} \hat {G}\left( {{{s}_{\xi }},\,{{s}_{w}}} \right) = \frac{2}{\pi }\Gamma \left( {{{s}_{\xi }}} \right)\Gamma \left( {{{s}_{w}}} \right)\cos \left( {\frac{\pi }{2}\left( {{{s}_{\xi }} + {{s}_{w}}} \right)} \right) \times \\ \, \times {{{\hat {F}}}_{c}}\left( {1 - {{s}_{\xi }} - {{s}_{w}},1 - {{s}_{w}}} \right), \\ \end{gathered} $
где

(11)
${{\hat {F}}_{c}}\left( {{{s}_{{v}}},{{s}_{y}}} \right) = \int\limits_0^\infty {\int\limits_0^\infty {{{F}_{c}}\left( {{v},y} \right){{{v}}^{{{{s}_{{v}}} - 1}}}{{y}^{{{{s}_{y}} - 1}}}d{v}dy} } ,$
(12)
${{F}_{c}}\left( {{v},y} \right) = - \frac{1}{{4{{\pi }^{2}}}}\int\limits_{{{\sigma }_{{v}}} - i\infty }^{{{\sigma }_{{v}}} + i\infty } {\int\limits_{{{\sigma }_{y}} - i\infty }^{{{\sigma }_{y}} + i\infty } {{{{\hat {F}}}_{c}}\left( {{{s}_{{v}}},{{s}_{y}}} \right){{{v}}^{{ - {{s}_{{v}}}}}}{{y}^{{ - {{s}_{y}}}}}d{{s}_{y}}d{{s}_{{v}}}} } .$

Определим необходимые ограничения на $\operatorname{Re} \left( {{{s}_{y}}} \right),$ $\operatorname{Re} \left( {{{s}_{{v}}}} \right)$. Будем отталкиваться от уже имеющихся ограничений, возникших для интегрального преобразования (9): $\operatorname{Re} \left( {{{s}_{\xi }}} \right),\;\operatorname{Re} \left( {{{s}_{w}}} \right) \in \left( {0,1} \right)$. Для преобразования Меллина (8) для функции $G\left( {\xi ,\,w} \right)$ такой диапазон параметров требует, чтобы положительно определенная функция $G\left( {\xi ,\,w} \right)$ не имела особенностей при $\xi = 0$ и $w = 0$, а также спадала на бесконечности не медленнее, чем ${1 \mathord{\left/ {\vphantom {1 {\left( {\xi w} \right)}}} \right. \kern-0em} {\left( {\xi w} \right)}}$. Условия для преобразования Меллина (11) для функции ${{F}_{c}}\left( {{v},y} \right)$ определить немного сложнее, поскольку она не является положительно определенной. Выразим нужный нам интеграл через $f\left( {x,\,y} \right)$

$\begin{gathered} \int\limits_0^\infty {{{F}_{c}}\left( {{v},y} \right){{{v}}^{{{{s}_{{v}}} - 1}}}d{v}} = \int\limits_0^\infty {\int\limits_0^\infty {f\left( {x,y} \right)\cos \left( {x{v}} \right){{{v}}^{{{{s}_{{v}}} - 1}}}d{v}dx} } = \\ \, = \Gamma \left( {{{s}_{{v}}}} \right)\cos \left( {{{\pi {{s}_{{v}}}} \mathord{\left/ {\vphantom {{\pi {{s}_{{v}}}} 2}} \right. \kern-0em} 2}} \right)\int\limits_0^\infty {f\left( {x,y} \right){{x}^{{ - {{s}_{{v}}}}}}dx} . \\ \end{gathered} $

Интеграл по $v$ берется явно и соответствует преобразованию Меллина для $\cos \left( {x{v}} \right)$ [28]. Сходимость интеграла по ν требует дополнительного условия $\operatorname{Re} \left( {{{s}_{{v}}}} \right) > 0$; если оно выполнено, то с учетом условий $\operatorname{Re} \left( {{{s}_{\xi }}} \right),\;\operatorname{Re} \left( {{{s}_{w}}} \right) \in \left( {0,1} \right)$ интеграл по x будет сходиться – при условии, что $f\left( {x,\,y} \right)$ не имеет особенностей при $x \to 0$ и спадает быстрее, чем ${1 \mathord{\left/ {\vphantom {1 x}} \right. \kern-0em} x}$ при $x \to \infty $. Рассматривая аналогично преобразование по y, получим, что для сходимости требуется, чтобы $f\left( {x,\,y} \right)$ не имела особенностей при $y \to 0$ и спадала быстрее, чем ${1 \mathord{\left/ {\vphantom {1 y}} \right. \kern-0em} y}$ при $y \to \infty $, при этом дополнительных ограничений на $\operatorname{Re} \left( {{{s}_{y}}} \right) \in \left( {0,\,1} \right)$ не возникает. Отметим, что накладываемые на функцию распределения условия выполняются автоматически, если наложить естественное физическое требование конечности числа частиц и суммарной кинетической энергии в системе. Далее можно показать, что если выполняются вышеперечисленные условия для f, то выполняются и аналогичные условия для функции G. В итоге, для применимости преобразований Меллина достаточно потребовать, чтобы константы в преобразованиях Меллина принадлежали следующим интервалам:

(13)
$\begin{gathered} \operatorname{Re} \left( {{{s}_{y}}} \right) \in \left( {0,1} \right),\quad \operatorname{Re} \left( {{{s}_{{v}}}} \right) \in \left( {0,1} \right), \\ \operatorname{Re} \left( {{{s}_{y}} - {{s}_{{v}}}} \right) \in \left( {0,1} \right). \\ \end{gathered} $

Решение уравнения (6) можно представить в виде следующей последовательности линейных преобразований:

$\begin{gathered} G(\xi ,w)\xrightarrow{{(8)}}\hat {G}({{s}_{\xi }},{{s}_{w}})\xrightarrow{{(10)}} \\ \, \to {{{\hat {F}}}_{c}}({{s}_{v}},{{s}_{y}})\xrightarrow{{(12)}}{{F}_{c}}(v,y)\xrightarrow{{(5)}}f(x,y), \\ \end{gathered} $
над стрелками указан номер формулы, описывающей преобразование. В явном виде эта цепочка выражается как

$\begin{gathered} f\left( {x,y} \right) = \int\limits_0^\infty {\int\limits_{{{\sigma }_{{v}}} - i\infty }^{{{\sigma }_{{v}}} + i\infty } {\int\limits_{{{\sigma }_{y}} - i\infty }^{{{\sigma }_{y}} + i\infty } {\int\limits_0^\infty {\int\limits_0^\infty {KG\left( {\xi ,w} \right)d\xi dwd{{s}_{y}}d{{s}_{{v}}}} } } } } d{v}, \\ K = - \frac{1}{{4{{\pi }^{2}}}}\frac{{{{\xi }^{{{{s}_{y}} - {{s}_{{v}}} - 1}}}{{w}^{{ - {{s}_{y}}}}}{{v}^{{ - {{s}_{{v}}}}}}{{y}^{{ - {{s}_{y}}}}}\cos \left( {x{v}} \right)}}{{\Gamma \left( {{{s}_{y}} - {{s}_{{v}}}} \right)\Gamma \left( {1 - {{s}_{y}}} \right)\sin \left( {\pi {{s}_{{v}}}{\text{/}}2} \right)}}. \\ \end{gathered} $

Смена порядка интегрирования позволяет получить решение в стандартной для линейного оператора форме в виде интегральной свертки

(14)
$f\left( {x,y} \right) = \int\limits_0^\infty {\int\limits_0^\infty {P\left( {\xi ,w,x,y} \right)G\left( {\xi ,w} \right)d\xi dw} } $
с ядром

$P\left( {\xi ,w,x,y} \right) = \int\limits_0^\infty {\int\limits_{{{\sigma }_{{v}}} - i\infty }^{{{\sigma }_{{v}}} + i\infty } {\int\limits_{{{\sigma }_{y}} - i\infty }^{{{\sigma }_{y}} + i\infty } {Kd{{s}_{y}}d{{s}_{{v}}}d{v}} } } .$

Последнее выражение можно упростить, поскольку, как уже отмечалось, интеграл по ${v}$ берется явно и соответствует преобразованию Меллина для $\cos \left( {x{v}} \right)$, получающиеся при этом тригонометрические функции сокращаются и остается выражение

(15)
$\begin{gathered} P\left( {\xi ,w,x,y} \right) = - \frac{1}{{4{{\pi }^{2}}}} \times \\ \, \times \int\limits_{{{\sigma }_{{v}}} - i\infty }^{{{\sigma }_{{v}}} + i\infty } {\int\limits_{{{\sigma }_{y}} - i\infty }^{{{\sigma }_{y}} + i\infty } {\frac{{{{\xi }^{{{{s}_{y}} - {{s}_{{v}}} - 1}}}{{x}^{{{{s}_{{v}}} - 1}}}{{{\left( {wy} \right)}}^{{ - {{s}_{y}}}}}\Gamma \left( {1 - {{s}_{{v}}}} \right)}}{{\Gamma \left( {{{s}_{y}} - {{s}_{{v}}}} \right)\Gamma \left( {1 - {{s}_{y}}} \right)}}d{{s}_{y}}d{{s}_{{v}}}} } . \\ \end{gathered} $

Заметим, что это выражение можно рассматривать как определение функции не четырех, а двух переменных, поскольку $P{{\xi }^{2}}$ зависит только от комбинаций $wy{\text{/}}\xi $ и $x{\text{/}}\xi $. В дальнейшем это свойство автомодельности нам не потребуется, однако оно позволяет табулировать ядро $P(\xi ,w,x,y)$ с помощью двумерной таблицы или интерполировать его с помощью функции двух переменных. На этом может основываться эффективный практический способ вычисления двумерной свертки (14).

Докажем, что получившийся обратный интегральный оператор действительно позволяет восстановить функцию $f\left( {x,\,y} \right)$. Применим оператор (14) с ядром (15) к функции $G\left( {\xi ,\,w} \right)$ наиболее общего вида (4):

$\int\limits_0^\infty {\int\limits_0^\infty {G\left( {\xi ,w} \right)P\left( {\xi ,w,x,y} \right)} } d\xi dw = $
(16)
$\begin{gathered} \, = - \frac{1}{{4{{\pi }^{2}}}}\int\limits_0^\infty {\int\limits_0^\infty {\left( {\int\limits_0^\infty {f\left( {\xi + wy{\kern 1pt} ',y{\kern 1pt} '} \right)dy{\kern 1pt} '} } \right) \times } } \\ \, \times \int\limits_{{{\sigma }_{{v}}} - i\infty }^{{{\sigma }_{{v}}} + i\infty } {\int\limits_{{{\sigma }_{y}} - i\infty }^{{{\sigma }_{y}} + i\infty } {\frac{{{{\xi }^{{{{s}_{y}} - {{s}_{{v}}} - 1}}}{{x}^{{{{s}_{{v}}} - 1}}}{{{\left( {wy} \right)}}^{{ - {{s}_{y}}}}}\Gamma \left( {1 - {{s}_{{v}}}} \right)}}{{\Gamma \left( {{{s}_{y}} - {{s}_{{v}}}} \right)\Gamma \left( {1 - {{s}_{y}}} \right)}} \times } } \\ \end{gathered} $
$\, \times d{{s}_{y}}d{{s}_{{v}}}d\xi dw.$

Далее делаем замену переменных $\xi = x{\kern 1pt} '\; - wy{\kern 1pt} '$ и меняем порядок интегрирования так, чтобы интегралы по $x{\kern 1pt} '$ и $y{\kern 1pt} '$ брались последними,

$\begin{gathered} \int\limits_0^\infty {\int\limits_0^\infty {\int\limits_0^\infty {{{\xi }^{{{{s}_{y}} - {{s}_{{v}}} - 1}}}{{w}^{{ - {{s}_{y}}}}}f\left( {\xi + wy{\kern 1pt} ',y{\kern 1pt} '} \right)dy{\kern 1pt} 'd\xi dw} } } = \\ \, = \int\limits_0^\infty {\int\limits_0^\infty {\int\limits_0^{x{\kern 1pt} '{\text{/}}y{\kern 1pt} '} {{{{\left( {x{\kern 1pt} '\; - wy{\kern 1pt} '} \right)}}^{{{{s}_{y}} - {{s}_{{v}}} - 1}}}{{w}^{{ - {{s}_{y}}}}}f\left( {x{\kern 1pt} ',y{\kern 1pt} '} \right)dwdx{\kern 1pt} 'dy{\kern 1pt} '} } } . \\ \end{gathered} $

В получившемся выражении интеграл по w является известным табличным интегралом [28]

(17)
$\begin{gathered} \int\limits_0^{x{\kern 1pt} '{\text{/}}y{\kern 1pt} '} {{{{\left( {x{\kern 1pt} '\; - wy{\kern 1pt} '} \right)}}^{{{{s}_{y}} - {{s}_{{v}}} - 1}}}{{w}^{{ - {{s}_{{v}}}}}}dw} = \\ \, = y{\kern 1pt} {{'}^{{{{s}_{y}} - {{s}_{{v}}} - 1}}}{{\left( {\frac{{x{\kern 1pt} '}}{{y{\kern 1pt} '}}} \right)}^{{ - {{s}_{{v}}}}}}\frac{{\Gamma \left( {1 - {{s}_{y}}} \right)\Gamma \left( {{{s}_{y}} - {{s}_{{v}}}} \right)}}{{\Gamma \left( {1 - {{s}_{{v}}}} \right)}}, \\ \end{gathered} $
сходящимся при $\operatorname{Re} \left( {{{s}_{y}} - {{s}_{v}}} \right) > 0$, $\operatorname{Re} \left( {1 - {{s}_{v}}} \right) > 0$, $x{\kern 1pt} '{\text{/}}y{\kern 1pt} ' > 0$. Подставляя (17) в исходное выражение (16) получим, что все гамма-функции сокращаются, а из оставшейся части набираются прямое и обратное преобразование Меллина [28], образующие дельта-функцию:

$\begin{gathered} \int\limits_0^\infty {\int\limits_0^\infty {G\left( {\xi ,w} \right)P\left( {\xi ,w,x,y} \right)} } d\xi dw = \\ \, = - \frac{1}{{4{{\pi }^{2}}}}\int\limits_0^\infty {\int\limits_0^\infty {\int\limits_{{{\sigma }_{{v}}} - i\infty }^{{{\sigma }_{{v}}} + i\infty } {\int\limits_{{{\sigma }_{y}} - i\infty }^{{{\sigma }_{y}} + i\infty } {\frac{1}{{x{\kern 1pt} '}}{{{\left( {\frac{x}{{x{\kern 1pt} '}}} \right)}}^{{{{s}_{{v}}} - 1}}}\frac{1}{{y{\kern 1pt} '}}{{{\left( {\frac{y}{{y{\kern 1pt} '}}} \right)}}^{{ - {{s}_{y}}}}} \times } } } } \\ \end{gathered} $
$\begin{gathered} \, \times f\left( {x{\kern 1pt} ',y{\kern 1pt} '} \right)d{{s}_{y}}d{{s}_{{v}}}dx{\kern 1pt} 'dy{\kern 1pt} ' = \int\limits_0^\infty {\int\limits_0^\infty {f\left( {x{\kern 1pt} ',y{\kern 1pt} '} \right)\delta \left( {x - x{\kern 1pt} '} \right)'} } \times \\ \, \times \delta \left( {y - y{\kern 1pt} '} \right)dx{\kern 1pt} 'dy{\text{'}}\;{\kern 1pt} = f\left( {x,y} \right). \\ \end{gathered} $

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

Нужно отметить, что математическое доказательство существования и единственности решений уравнения (6) для произвольного G достаточно трудоемко, поскольку это уравнение со знакопеременным ядром [27]. Однако того факта, что в качестве G мы берем только функции, являющиеся образом для некоторой достаточно быстро убывающей функции f, возникающей из физической задачи, достаточно для того, чтобы наше решение существовало и определялось соотношениями (14) и (15).

5. РАЗЛОЖЕНИЕ РЕШЕНИЯ ПО ОРТОГОНАЛЬНЫМ ПОЛИНОМАМ

Хотя нам и удалось предъявить аналитическое выражение для обратного оператора в квадратурах, тем не менее, такое решение является во многом формальным и неудобным для реального применения. Поэтому мы предложим здесь еще один способ обращения уравнения (4).

По данной функции $G\left( {\xi ,\,w} \right)$ мы можем восстановить дискретный набор функций ${{G}_{k}}\left( x \right)$, определенный как

${{G}_{0}}\left( x \right) = G\left( {x,0} \right),$
${{G}_{1}}\left( x \right) = \mathop {\lim }\limits_{w \to 0} \int\limits_\infty ^x {\frac{\partial }{{\partial w}}G({{\xi }_{1}},w)d{{\xi }_{1}}} ,$
${{G}_{k}}\left( x \right) = \mathop {\lim }\limits_{w \to 0} \int\limits_\infty ^x {\int\limits_\infty ^{{{\xi }_{1}}} {...} \int\limits_\infty ^{{{\xi }_{{k - 1}}}} {\frac{{{{\partial }^{k}}}}{{\partial {{w}^{k}}}}G\left( {{{\xi }_{k}},w} \right)d{{\xi }_{k}}...d{{\xi }_{2}}d{{\xi }_{1}}} } .$

Легко показать, что

(18)
${{G}_{k}}\left( x \right) = \int\limits_0^\infty {{{y}^{k}}f\left( {x,y} \right)dy} ,$
то есть эти функции описывают зависимости от переменной x моментов функции распределения по переменной y (напомним, что x и y есть, соответственно, кинетическая энергия продольного движения и циклотронного вращения частицы).

Из функции ${{G}_{k}}$ можно составить линейные комбинации таким образом, чтобы под интегралом сформировались какие-либо ортогональные на интервале $\left[ {0,\,\infty } \right)$ полиномы. Рассмотрим в качестве базового примера полиномы Лагерра [29]

${{L}_{n}}\left( y \right) = \frac{{{{e}^{y}}}}{{n!}}\frac{{{{d}^{n}}}}{{d{{y}^{n}}}}\left( {{{e}^{{ - y}}}{{y}^{n}}} \right) = \sum\limits_{k = 0}^n {\left( {\begin{array}{*{20}{c}} n \\ k \end{array}} \right)\frac{{{{{\left( { - 1} \right)}}^{k}}}}{{k!}}{{y}^{k}}} ,$
с условием ортогональности

$\int\limits_0^\infty {{{L}_{n}}\left( y \right){{L}_{k}}\left( y \right){{e}^{{ - y}}}dy} = {{\delta }_{{nk}}}.$

Здесь $\left( {\begin{array}{*{20}{c}} n \\ k \end{array}} \right)$ обозначает биномиальный коэффициент. Введем обозначение ${{L}_{n}}\left[ {{{G}_{k}}\left( x \right)} \right]$ – формальный полином Лагерра ${{L}_{n}}\left( y \right)$, в котором проведена замена ${{y}^{k}} \to {{G}_{k}}\left( x \right)$. С учетом (18) получаем набор интегральных уравнений для $n = 0,1,2...$

(19)
${{L}_{n}}\left[ {{{G}_{k}}\left( x \right)} \right] = \int\limits_0^\infty {{{L}_{n}}\left( y \right)f\left( {x,y} \right)dy} .$

Переменная x стала здесь просто параметром, а по переменной y мы можем разложить функцию распределения по полной системе полиномов Лагерра:

(20)
$f\left( {x,y} \right) = \sum\limits_{n = 0}^\infty {{{f}_{n}}\left( x \right){{L}_{n}}\left( y \right)\exp \left( { - y} \right)} ,$
где ${{f}_{n}}\left( x \right)$ – коэффициенты разложения, независимо определенные для каждого значения x. В силу ортогональности полиномов Лагерра уравнения (19) принимают вид

${{L}_{n}}\left[ {{{G}_{k}}\left( x \right)} \right] = {{f}_{n}}\left( x \right).$

Таким образом, мы построили решение уравнения (6) в виде

(21)
$f\left( {x,y} \right) = \sum\limits_{n = 0}^\infty {\sum\limits_{k = 0}^n {\left( {\begin{array}{*{20}{c}} n \\ k \end{array}} \right)\frac{{{{{\left( { - 1} \right)}}^{k}}}}{{k!}}{{G}_{k}}(x)} {{L}_{n}}\left( y \right)\exp \left( { - y} \right)} .$

Получившийся ряд имеет важное свойство, обеспечивающее его эффективность при решении физической задачи. Частичная сумма из первых N членов по n дает функцию ${{f}_{N}}$, обладающую тем же набором первых N моментов по переменной y, что и искомая функция, то есть, если в соотношения (18) подставить ${{f}_{N}}$, то они будут выполняться строго при $k = 0,1...,N - 1$ для всех x. Соответственно, остаточный член $f - {{f}_{N}}$ обладает нулевыми первыми N моментами. Например, учет только первых трех членов дает верные среднее значение и дисперсию распределений по y, зависящие от x как от параметра.

Обсудим сходимость получившегося ряда. Для этого достаточно рассмотреть аналог теоремы Парсеваля для разложения (20), следующий из условий ортогональности полиномов Лагерра,

$\begin{gathered} \int\limits_0^\infty {{{f}^{2}}\left( {\xi ,y} \right){{e}^{y}}dy} = \\ \, = \sum\limits_{m = 0}^\infty {\sum\limits_{n = 0}^\infty {{{f}_{m}}\left( \xi \right){{f}_{n}}\left( \xi \right)\int\limits_0^\infty {{{L}_{m}}\left( y \right)} } {{L}_{n}}\left( y \right){{e}^{{ - y}}}dy} = \\ \, = \sum\limits_{m = 0}^\infty {f_{m}^{2}\left( \xi \right)} . \\ \end{gathered} $

Сумма квадратов коэффициентов разложения ограничена, только если сходится интеграл в левой части, то есть если функция $f\left( {x,\,y} \right)$ спадает по y быстрее, чем $\exp \left( {{{ - y} \mathord{\left/ {\vphantom {{ - y} 2}} \right. \kern-0em} 2}} \right)$. Для оценки скорости сходимости рассмотрим функцию вида

(22)
$f = g\left( x \right)\exp \left( { - ay} \right).$

Коэффициенты ${{f}_{n}}\left( x \right)$ разложения (20) для этой функции

${{f}_{n}}\left( x \right) = \int\limits_0^\infty {L\left( y \right)\exp \left( { - ay} \right)dy} g\left( x \right) = {{a}^{{ - 1}}}{{\left( {1 - 1{\text{/}}a} \right)}^{n}}g(x).$

Видно, что при выполнении условия сходимости амплитуда коэффициентов экспоненциально уменьшается с ростом n. То есть для бесконечно гладкой функции коэффициенты ряда (21) уменьшаются в геометрической прогрессии Иначе говоря, ряд если сходится, то сходится быстро. В случае, когда есть разрыв значения функции или ее первой производной, скорость сходимости также можно оценить на конкретных примерах. Рассмотрим разложение по полиномам Лагерра для модельных функций ${{g}_{0}} = St\left( {1 - y} \right)$, имеющей разрыв значения, и ${{g}_{1}} = \left( {1 - y} \right)St\left( {1 - y} \right)$, имеющей разрыв первой производной; здесь

$St\left( y \right) = \left\{ {\begin{array}{*{20}{c}} {1,\quad y > 0} \\ {0,\quad y < 0} \end{array}} \right..$

На рис. 1 показаны частичные суммы первых 11 и 101 членов ряда для функции ${{g}_{0}}$ (слева) и ${{g}_{1}}$ (справа). Видно, что для функции с разрывом значения учет даже сотни слагаемых ряда не дает удовлетворительной сходимости, а для функции с разрывом производной удовлетворительная сходимость достигается уже для десяти слагаемых, а сотня слагаемых повторяет функцию практически идеально.

Рис. 1.

Частичные суммы первых 11 и 101 членов ряда (21) для функций ${{g}_{0}}$ (слева) и ${{g}_{1}}$ (справа). Сплошной черной линией показаны исходные функции, серыми штриховыми линиями показаны суммы первых 11 членов ряда, черными штриховыми – суммы первых 101 члена ряда.

Заметим, что расходимость ряда (21) для функций распределения со степенными хвостами вытекает из самой схемы построения нашего решения, которая требует сходимости набора интегралов (18). Для восстановления функций распределения со степенными хвостами $f\left( {x,y} \right)\sim 1{\text{/}}{{y}^{\alpha }}$ мы не можем использовать члены ряда с $k > \alpha - 1$. Тем не менее, даже в тех случаях, когда наш ряд формально не сходится из-за обращения в бесконечность высоких моментов исходной функции, мы можем использовать конечное число конечных членов ряда (21) для удовлетворительного приближенного решения исходной задачи. С физической точки зрения это не удивительно, поскольку у приближенной функции распределения все младшие моменты гарантированно совпадают с точными значениями для восстанавливаемой функции f. Процедура построения приближенного решения с использованием конечного числа членов ряда рассмотрена в следующем разделе.

6. ВОССТАНОВЛЕНИЕ ФУНКЦИИ ПО КОНЕЧНОМУ ЧИСЛУ ЧЛЕНОВ РЯДА

Для того чтобы донести идею, развиваемую ниже, начнем с простого частного случая. У нас есть определенный произвол при выборе системы ортогональных полиномов. В частности, мы можем воспользоваться вместо ${{L}_{n}}(y)$ полиномами ${{L}_{n}}(\beta y)$, где $\beta > 0$ не зависит от y. С учетом перемасштабированного соотношения ортогональности, нетрудно получить набор различных эквивалентных записей ряда (21) с дополнительным свободным параметром β:

(23)
$\begin{gathered} f\left( {x,y} \right) = \sum\limits_{n = 0}^\infty {\sum\limits_{k = 0}^n {\left( {\begin{array}{*{20}{c}} n \\ k \end{array}} \right)\frac{{{{{\left( { - 1} \right)}}^{k}}}}{{k!}}{{\beta }^{{k + 1}}}{{G}_{k}}(x)} \times } \\ \, \times {{L}_{n}}\left( {\beta y} \right)\exp \left( { - \beta y} \right). \\ \end{gathered} $

Этот ряд наследует основное свойство разложения (21) – частичная сумма из первых N членов дает приближенную функцию ${{f}_{N}}$, воспроизводящую точно первые N моментов по переменной y, то есть

$\int\limits_0^\infty {{{y}^{k}}{{f}_{N}}\left( {x,y} \right)dy} = \int\limits_0^\infty {{{y}^{k}}f\left( {x,y} \right)dy} ,\quad k = 0,1,...,N - 1;$
при этом остаточный член обладает нулевыми первыми N моментами.

Видно, что подбирая β можно обеспечить сходимость полученного ряда для любой экспоненциально спадающей функции f, включая те, для которых исходный ряд (21) расходится. Более того, нетрудно видеть, что все функции вида (22) будут точно описываться одним членом ряда (23), если $\beta = a$. В реальной задаче мы априори не знаем функцию f, более того, асимптотика этой функции при $y \to \infty $ скорее всего не является экспоненциальной. Тем не менее, свобода выбора параметра $\beta $ и тут позволяет оптимизировать сходимость бесконечного ряда (23) для f и уменьшить ошибку приближенного решения ${{f}_{N}}$, связанную с отбрасыванием остаточного члена. Такое отбрасывание в реальной задаче в первую очередь связано с естественной ограниченностью экспериментальных данных, которые мы не можем бесконечно дифференцировать, надеясь извлечь из этого какую-то информацию. При этом возникает задача оптимизации восстановления функции распределения в ситуации, когда доступно только конечное число членов ряда ${{G}_{n}}$. Для ее решения можно предложить следующий подход. Учтем первые N членов разложения, а выбором β добьемся, чтобы $(N + 1)$-й член разложения был равен нулю. Это предполагает решение алгебраического уравнения

(24)
$\sum\limits_{k = 0}^{N + 1} {\left( {\begin{array}{*{20}{c}} {N + 1} \\ k \end{array}} \right)\frac{{{{{\left( { - 1} \right)}}^{k}}}}{{k!}}{{G}_{k}}(x){{\beta }^{k}}} = 0$
для величины β при каждом значении x. Полученная приближенная функция будет состоять из N членов разложения и точно воспроизводить $N + 1$ первых моментов функции распределения по переменной y. Если вспомнить смысл переменных x и y, то уравнение (24) можно интерпретировать с физической точки зрения как уравнение для эффективной поперечной температуры, наиболее точно описывающей распределение частиц с заданной продольной энергией.

Приведенный метод можно усовершенствовать, если в качестве системы ортогональных полиномов взять обобщенные полиномы Лагерра, определенные как

$L_{n}^{{(\alpha )}}\left( y \right) = \frac{{{{y}^{{ - \alpha }}}{{e}^{y}}}}{{n!}}\frac{{{{d}^{n}}}}{{d{{y}^{n}}}}\left( {{{e}^{{ - y}}}{{y}^{{n + \alpha }}}} \right) = \sum\limits_{k = 0}^n {\left( {\begin{array}{*{20}{c}} {n + \alpha } \\ {n - k} \end{array}} \right)\frac{{{{{\left( { - 1} \right)}}^{k}}}}{{k!}}{{y}^{k}}} ,$
где $\left( {\begin{array}{*{20}{c}} {n + \alpha } \\ {n - k} \end{array}} \right)$ обозначает обобщенные биномиальные коэффициенты [30]. Возьмем в качестве базиса разложения набор полиномов $L_{n}^{{(\alpha )}}\left( {\beta y} \right)$, удовлетворяющих соотношению ортогональности

$\int\limits_0^\infty {L_{n}^{{(\alpha )}}\left( {\beta y} \right)L_{k}^{{(\alpha )}}\left( {\beta y} \right){{y}^{\alpha }}{{e}^{{ - \beta y}}}dy} = \frac{{\Gamma \left( {n + \alpha + 1} \right)}}{{{{\beta }^{{\alpha + 1}}}n!}}{{\delta }_{{nk}}}.$

Тогда разложение искомой функции примет вид

(25)
$\begin{gathered} f\left( {x,y} \right) = \sum\limits_{n = 0}^\infty {{{f}_{n}}\left( x \right)L_{n}^{{(\alpha )}}\left( {\beta y} \right){{y}^{\alpha }}{{e}^{{ - \beta y}}}} , \\ {{f}_{n}}(x) = \frac{{{{\beta }^{{\alpha + 1}}}n!}}{{\Gamma \left( {n + \alpha + 1} \right)}}L_{n}^{{\left( \alpha \right)}}\left[ {{{\beta }^{k}}{{G}_{k}}\left( x \right)} \right], \\ \end{gathered} $
где $L_{n}^{{\left( \alpha \right)}}\left[ {{{\beta }^{k}}{{G}_{k}}\left( x \right)} \right]$ мы понимаем как выражение, получаемое из полинома $L_{n}^{{(\alpha )}}\left( y \right)$ в результате формальной замены ${{y}^{k}} \to {{\beta }^{k}}{{G}_{k}}\left( x \right)$. Полученный ряд сходится в широком диапазоне параметров α и β для достаточно быстро спадающей на бесконечности функции $f(x,y)$. Это позволяет нам, варьируя эти параметры, оптимизировать приближение функции распределения конечным числом членов ряда, при этом мы можем варьировать эти параметры для каждого значения x независимо.

Как и раньше, частичная сумма ${{f}_{N}}$ первых N членов ряда (25) точно воспроизводит первые N моментов искомой функции при любых α и β. Потребуем теперь, чтобы $(N + 1)$-й и $(N + 2)$-й члены ряда строго равнялись нулю. Это приведет к алгебраической системе

(26)
$\left\{ {\begin{array}{*{20}{c}} {L_{{N + 1}}^{{\left( \alpha \right)}}\left[ {{{\beta }^{k}}{{G}_{k}}\left( x \right)} \right] = 0} \\ {L_{{N + 2}}^{{\left( \alpha \right)}}\left[ {{{\beta }^{k}}{{G}_{k}}\left( x \right)} \right] = 0} \end{array}} \right.,$
определяющей параметры $\alpha (x)$ и $\beta (x)$ для каждого значения x независимо. Соответствующая этому выбору параметров приближенная функция ${{f}_{N}}$ будет точно воспроизводить уже $N + 2$ момента функции распределения.

Разложение по обобщенным полиномам Лагерра позволяет учесть специфику распределения частиц в адиабатической магнитной ловушке, а именно, распределения с незаполненным конусом потерь. Напомним, что y соответствует кинетической энергии поперечного движения (вращения) частицы. Поэтому распределения с конусом потерь хорошо моделируются функциями вида ${{y}^{\alpha }}\exp \left( { - \beta y} \right)$, которые можно описать точно одним членом ряда (25). Можно ожидать, что при восстановлении функции распределения с обедненным распределением ионов в области потерь по реальным данным разложение (25) позволит обеспечить заданную точность меньшим числом членов ряда, чем разложение (23), опирающееся на “менее физические” базисные функции.

Рассмотрим на примерах, как работает восстановление функции распределения через разложение по полиномам Лагерра с одним варьирующимся масштабом β и через обобщенные полиномы Лагерра с двумя варьирующимися параметрами α и β. Процедура проверки все время одинаковая: по заданной аналитически функции $f(x,y)$ рассчитываем диагностический отклик $G(\xi ,w)$, вычисляем моменты ${{G}_{k}}(x)$, а затем применяем решения (23) или (25) с конечным числом членов ряда. Мы всегда рассматриваем три члена с учетом принудительно зануляемых, то есть восстановленная приближенная функция содержит максимум два слагаемых, если варьируется один параметр разложения, и всего одно слагаемое, если варьируются два параметра разложения. Даже в таких жестких условиях уравнения (24) или (26) гарантируют, что первые три момента, то есть число частиц, полная энергия циклотронного вращения и ее дисперсия восстанавливаются точно для каждой группы частиц, имеющих заданную продольную скорость.

На рис. 2 представлены результаты восстановления распределений вида $f\left( {x,\,y} \right) = {{f}_{x}}\left( x \right){{f}_{y}}\left( y \right)$. В данном случае множитель, отвечающий за зависимость от x, восстанавливается точно конечной суммой ряда, при этом параметры α и β одинаковы при всех x. Рисунок 2а относится к элементарному случаю анизотропного максвелловского распределения, для которого первый член разложения (25) дает точный ответ при $\alpha = 0$ и $\beta = 2$ (одинаковые графики слева и справа). Если же зафиксировать “неправильную” поперечную температуру $\beta = 1$ и оптимизировать параметр α, то первые два члена разложения обеспечивают относительную точность на уровне 5% во всей “тепловой области” ${\text{|}}x + y{\text{|}} < 3$ (график в центре, относительная точность определяется по максимальному относительному сдвигу линий уровня). Случай (б) отвечает максвелловскому распределению с конусом потерь, для которого первый член разложения (25) дает точный ответ при $\alpha = 3$ и $\beta = 2$ (график справа). Остальные два метода (слева – подгонка поперечной температуры β при $\alpha = 0$ и в центре – подгонка α при навязанной неправильной температуре $\beta = 1$) позволяют удовлетворительно восстановить функцию распределения всего по двум членам разложения. Случай (в) тоже отвечает максвелловскому распределению с конусом потерь, но более сложного вида. Теперь мы уже не можем точно восстановить функцию распределения одним членом ряда. Видно, что оптимизация одновременно по двум параметрам α и β, приводящая к простой формуле $f\left( {x,\,y} \right)$$1.41 {{y}^{{1.56}}}\exp \left( { - x - 1.72 y} \right) $ (график справа), дает лучшее приближение, чем двухчленные приближения, получаемые в результате оптимизации по каждому из параметров по отдельности (слева и в центре). Таким образом, на приведенных “учебных” примерах мы убедились, что разложение (25) с дополнительным свободным параметром позволяет существенно повысить точность восстановления по сравнению с однопараметрическим разложением (23).

Рис. 2.

Примеры восстановления по первым членам ряда (25) для модельных функций распределения: а) $f\left( {x,y} \right) = \exp \left( { - x - 2y} \right)$, б) $f\left( {x,y} \right) = {{y}^{3}}\exp \left( { - x - 2y} \right)$, в) $f\left( {x,y} \right) = \left( {0.5\,y + {{y}^{2}} + 0.3\,{{y}^{3}}} \right)\exp \left( { - x - 2y} \right)$ для трех способов оптимизации: по β при $\alpha = 0$ (слева), по α при β = 1 (в центре) и по обоим параметрам α и β (справа). Сплошные черные линии соответствуют линиям уровня исходной функции, штриховые – линиям уровня восстановленных функций. Получившиеся в результате оптимизации параметры разложения приведены на каждом графике сверху.

Теперь рассмотрим пример, моделирующий физическую задачу. На рис. 3 приведены линии уровня функции распределения, описывающей горячие ионы, рождающиеся в результате обдирки нейтрального пучка, введенного в плазму под углом 45° к магнитному полю, и последующей кулоновской релаксации, – случай, отвечающий условиям эксперимента на установке ГДЛ. Исходная функция распределения задается либо анизотропным максвелловским распределением, сдвинутым в точку, отвечающую скорости нейтральных атомов (график слева),

(27)
$\begin{gathered} {{f}_{1}}\left( {x,y} \right) = \\ \, = \exp \left( { - \frac{{{{{\left( {\sqrt x + \sqrt y - \sqrt {2\varepsilon } } \right)}}^{2}}}}{{2\delta \varepsilon }} - \frac{{{{{\left( {\sqrt x - \sqrt y } \right)}}^{2}}}}{{2\delta \varepsilon \delta {{\theta }^{2}}}}} \right), \\ \end{gathered} $
либо таким же распределением, но обрезанным со стороны высоких энергий (график справа),

(28)
${{f}_{2}}\left( {x,y} \right) = {{f}_{1}}\left( {x,y} \right) {\text{St}}\left( {{{\varepsilon }_{{\text{m}}}} - x - y} \right).$
Рис. 3.

Восстановление функции распределения пучка ионов по первому члену ряда (25) в случае оптимизации по обоим параметрам α и β. Левый график соответствует гладкому распределению (27), правый график – распределению с разрывом (28). Параметры распределений $\varepsilon = 30$, $\delta \varepsilon = 1$, $\delta \theta = 0.3$, ${{\varepsilon }_{m}} = 40$. Сплошные черные линии соответствуют линиям уровня исходной функции, штриховые – линиям уровня восстановленных функций.

Второе распределение моделирует особенности релаксации функции распределения за счет кулоновских столкновений с частицами фоновой плазмы, при которых доминирующими эффектами являются диффузия по направлениям скорости и потери энергии за счет трения. В отличие от рассмотренных выше, данные распределения уже не являются факторизованными, поэтому при восстановлении мы провели оптимизацию по α и β отдельно для каждого значения x. Результаты восстановления изображены на графиках штриховыми линиями уровня. Точность восстановления функции ${{f}_{1}}$ с использованием всего одного члена ряда (25) составляет около 1%. Видно, что для функции с разрывом ${{f}_{2}}$ восстановление несколько хуже, однако его точность все равно достаточно велика. Например, наиболее трудно восстанавливаемый этим методом параметр – положение разрыва, отвечающего максимальной энергии ${{\varepsilon }_{{\text{m}}}}$, – восстанавливается с точностью около 5% при определении максимальной энергии по уровню функции распределения 0.01. Отметим также, что получаемое из малого числа членов оптимизированного ряда приближение для функции с разрывом положительно определено, в отличие от результата, получаемого путем суммирования большого числа членов неоптимизированного ряда (рис. 1). Напомним, что метод автоматически обеспечивает верные значения для числа частиц, поперечной энергии и ее дисперсии при каждом значении продольной скорости частицы.

В целом, обобщая наш опыт применения указанного метода к различным функциям распределения, можно сформулировать следующее эмпирическое заключение. Если выбором α и β добиться того, что два последовательных члена ряда (25) равны нулю, то сумма всех последующих неучтенных членов ряда будет мала по сравнению с ${{f}_{N}}$, или, другими словами, сумма первых ненулевых членов будет хорошо повторять искомую функцию. Это утверждение верно и для отдельных функций с “патологическими” особенностями, например, для функций со степенными хвостами, для которых все старшие моменты начиная с некоторого расходятся.

7. ЗАКЛЮЧЕНИЕ

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

Задача восстановления сводится к решению интегрального уравнения, для которого предложено два метода решения. Оба предъявленных аналитических метода: через преобразование Меллина и через разложение по ортогональным полиномам, – допускают реализацию на практике, однако метод разложения по ортогональным полиномам выглядит несколько более предпочтительно с точки зрения физической обоснованности и удобства реализации. При использовании разложения по ортогональным полиномам каждый член ряда отвечает за восстановление соответствующего момента функции распределения. Естественно, что поиск решения в виде ряда имеет и свои недостатки, к которым можно отнести необходимость вычисления производных все более высокого порядка от экспериментальных данных для вычисления членов ряда. При этом возникает проблема оптимизации процедуры восстановления функции распределения при принципиальной доступности из экспериментальных данных только конечного числа членов бесконечного ряда. Возможности такой оптимизации рассмотрены в данной статье на примере случая, когда нам доступны только три члена ряда, что требует вычисления первой и второй производной от экспериментальных данных. Показана возможность очень хорошего восстановления модельных функций распределения за счет выбора масштаба и порядка обобщенных полиномов Лагерра.

Работа выполнена при поддержке Российского научного фонда, проект № 19-72-20139.

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

  1. Bindslev H. // Rev. Sci. Instrum. 1999. V. 70. P. 1093. https://doi.org/10.1063/1.1149521

  2. Machuzak J.S., Woskov P.P., Gilmore J., Bretz N.L., Park H.K., Bindslev H. // Rev. Sci. Instrum. 1997. V. 68. P. 458. https://doi.org/10.1063/1.1147605

  3. Bindslev H., Hoekzema J.A., Egedal J., Fessey J.A., Hughes T.P., Machuzak J.S. // Phys. Rev. Lett. 1999. V. 83. P. 3206. https://doi.org/10.1103/PhysRevLett.83.3206

  4. Porte L., Bindslev H., Hoekzema F., Machuzak J., Woskov P., Van Eester D. // Rev. Sci. Instrum. 2001. V. 72. P. 1148. https://doi.org/10.1063/1.1322617

  5. Meo F., Bindslev H., Korsholm S.B., Tsakadze E.L., Walker C.I., Woskov P., Vayakis G. // Rev. Sci. Instrum. 2004. V. 75. P. 3585. https://doi.org/10.1063/1.1788853

  6. Korsholm S.B., Bindslev H., Meo F., Leipold F., Michelsen P.K., Michelsen S., Nielsen S.K., Tsakadze E.L., Woskov P., Westerhof E., FOM ECRH team, Ooster-beek J.W., Hoekzema J., Leuterer F., Wagner D., ASDEX Upgrade ECRH team // Rev. Sci. Instrum. 2006. V. 77. 10E514. https://doi.org/10.1063/1.2217921

  7. Stejner M., Nielsen S.K., Korsholm S.B., Salewski M., Bindslev H., Furtula V., Leipold F., Meo F., Michelsen P.K., Moseev D., Bürger A., Kantor M., de Baar M. // Rev. Sci. Instrum. 2010. V. 81. 10D515. https://doi.org/10.1063/1.3475540

  8. Meo F., Bindslev H., Korsholm S.B., Furtula V., Leute-rer F., Leipold, F., Michelsen P. K., Nielsen S.K., Salewski M., Stober J., Wagner D., Woskov P. // Rev. Sci. Instrum. 2008. V. 79. 10E501. https://doi.org/10.1063/1.2989140

  9. Suvorov E.V., Holzhauer E., Kasparek W., Lubyako L.V., Burov A.B., Dryagin Y.A., Fil’chenkov S.E., Fraiman A.A., Kukin L.M., Kostrov A.V., Ryndyk D.A., Shtanyuk A.M., Skalyga N.K., Smolyakova O.B., Erckmann V., Geist T., Kick M., Laqua H., Rust M., W7-AS Team, ECRH Team and NBI Team // Plasma Phys. Control. Fusion. 1997. V. 39(12B). P. B337. https://doi.org/10.1088/0741-3335/39/12B/026

  10. Suvorov E.V., Holzhauer E., Kasparek W., Burov A.B., Dryagin Y.A., Fil’chenkov S.E., Fraiman A.A., Lubya-ko L.V., Ryndyk D.A., Skalyga N.K., Smolyakova O.B., Erckmann V., Geist T., Kick M., Rust N., W7-AS Team and ECRH Team // Nucl. Fusion. 1998. V. 38(5). P. 661. https://doi.org/10.1088/0029-5515/38/5/302

  11. Shalashov A.G., Suvorov E.V., Lubyako L.V., Maassberg H. and W7-AS Team // Plasma Phys. Control. Fusion. 2003. V. 45. P. 395. https://doi.org/10.1088/0741-3335/45/4/306

  12. Kubo S., Nishiura M., Tanaka K., Shimozuma T., Yoshimura Y., Igami H., Takahash H., Mutoh T., Tamura N., Tatematsu Y., Saito T., Notake T., Korsholm S.B., Meo F., Nielsen S.K., Salewski M., Stejner M. // Rev. Sci. Instrum. 2010. V. 81. 10D535. https://doi.org/10.1063/1.3481165

  13. Laqua H.P., Baldzuhn J., Braune H., Bozhenkov S., Brunner K.J., Kazakov Ye.O., Marsen S., Moseev D., Stange T., Wolf R.C., Zanini M., W7-X Team // European Physical Journal Conferences. 2018. V. 187(2). 01011. https://doi.org/10.1051/epjconf/201818701011

  14. Moseev D., Salewski M., Garcia-Muñoz M., Geiger B., Nocente M. // Rev. Mod. Plasma Phys. 2018. V. 2. 7. https://doi.org/10.1007/s41614-018-0019-4

  15. Bindslev H., Meo F., Tsakadze E.L., Korsholm S.B., Woskov P. // Rev. Sci. Instrum. A. 2004. V. 75(10). P. 3598. https://doi.org/10.1063/1.1779620

  16. Bagryansky P.A., Shalashov A.G., Gospodchikov E.D., Lizunov A.A., Maximov V.V., Prikhodko V.V., Soldatkina E.I., Solomakhin A.L., Yakovlev D.V. // Phys. Rev. Lett. 2015. V.114 (20). 205001. https://doi.org/10.1103/PhysRevLett.114.205001

  17. Bagryansky P.A., Anikeev A.V., Denisov G.G., Gospodchikov E.D., Ivanov A.A., Lizunov A.A., Kovalenko Yu.V., Malygin V.I., Maximov V.V., Korobeinikova O.A., Murakhtin S.V., Pinzhenin E.I., Prikhodko V.V., Savkin V.Ya., Shalashov A.G., Smolyakova O.B., Soldatkina E.I., Solomakhin A.L., Yakovlev D.V., Zaytsev K.V. // Nucl. Fusion. 2015. V. 55(5). 053009. https://doi.org/10.1088/0029-5515/55/5/053009

  18. Yakovlev D.V., Shalashov A.G., Gospodchikov E.D., Maximov V.V., Prikhodko V.V., Savkin V.Ya., Soldatkina E.I., Solomakhin A.L., Bagryansky P.A. // Nucl. Fusion. 2018. V. 58. 094001. https://doi.org/10.1088/1741-4326/aacb88

  19. Simonen T.C. // J. Fusion Energy. 2016. V. 35(1). P. 63. https://doi.org/10.1007/s10894-015-0017-2

  20. Gota H., Binderbauer M.W., Tajima T., Putvinski S., Tuszewski M., Dettrick S, Garate E., Korepanov S., Smirnov A., Thompson M.C., Trask E., Yang X., Schmitz L., Lin Z., Ivanov A.A., Asai T., Allfrey I., Andow R., Beall M., Bolte N., Bui D.Q., Cappello M., Ceccherini F., Clary R., Cheung A.H., Conroy K., Deng B.H., Douglass J., Dunaevsky A., Feng P., Fulton D., Galeotti L., Granstedt E., Griswold M., Gupta D., Gupta S., Hubbard K., Isakov I., Kinley J.S., Knapp K., Magee1 R., Matvienko V., Mendoza R., Mok Y., Necas A., Primavera S., Onofri M., Osin D., Rath N., Roche T., Romero J., Schindler T., Schroeder J.H., Sevier L., Sheftman D., Sibley A., Song Y., Steinhauer L.C., Valentine T., Van Drie A.D., Walters J.K., Waggoner W., Yushmanov P., Zhai K., the TAE Team // Nucl. Fusion. 2017. V. 57. 116021. https://doi.org/10.1088/1741-4326/aa7d7b

  21. Bagryansky P.A., Beklemishev A.D., Postupaev V.V. // J. Fusion Energy. 2019. V. 38. P. 162. https://doi.org/10.1007/s10894-018-0174-1

  22. Shalashov A.G., Gospodchikov E.D., Khusainov T.A., Lubyako L.V., Smolyakova O.B., Solomakhin A.L. // Plasma Phys. Control. Fusion. 2020. V. 62. 065010. https://doi.org/10.1088/1361-6587/ab83cc

  23. Ivanov A.A., Prikhodko V.V. // Plasma Phys. Control. Fusion. 2013. V. 55. 063001. https://doi.org/10.1088/0741-3335/55/6/063001

  24. Salewski M., Nielsen S.K., Bindslev H., Furtula V., Gorelenkov N.N., Korsholm S.B., Leipold F., Meo F., Michelsen P.K., Moseev D., Stejner M. // Nucl. Fusion. 2011. V. 51. 083014. https://doi.org/10.1088/0029-5515/51/8/083014

  25. Salewski M., Geiger B., Nielsen S.K., Bindslev H., García-Muñoz M., Heidbrink W.W., Korsholm S.B., Leipold F., Meo F., Michelsen P.K., Moseev D., Stejner M., Tardini G., the ASDEX Upgrade team // Nucl. Fusion. 2012. V. 52. 103008. https://doi.org/10.1088/0029-5515/52/10/103008

  26. Юшманов Е.Е. // ЖЭТФ. 1966. Т. 49. С. 588. [E.E.Yushmanov // J. Exp. Theoret. Phys. (U.S.S.R.). 1966. V. 22. P. 409]

  27. Краснов М.Л. Интегральные уравнения. Введение в теорию. М.: Наука, 1975.

  28. Бейтмен Г., Эрдейи А. Таблицы интегральных преобразований. Том 1. Преобразования Фурье, Лапласа, Меллина. М.: Наука, 1969.

  29. Бейтмен Г., Эрдейи А. Функции Бесселя, функции параболического цилиндра, ортогональные многочлены. Высшие трансцендентные функции. Т. 2. М.: Наука, 1966.

  30. Abramowitz M., Stegun I.A., Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables. Applied Mathematics Series. Washington, D.C., New York: United States Department of Commerce, National Bureau of Standards, Dover Publications, 1983.

  31. Bracewell R. The Fourier Transform and Its Applications, 3rd ed. New York: McGraw-Hill, 1999.

Дополнительные материалы отсутствуют.