Астрономический вестник, 2019, T. 53, № 1, стр. 34-60

Моделирование ускорения протонов в магнитном острове в складке гелиосферного токового слоя

О. В. Мингалев a*, О. В. Хабарова b**, Х. В. Малова cd, И. В. Мингалев a, Р. А. Кислов c, М. Н. Мельник a, П. В. Сецко a, Л. М. Зеленый c, G. P. Zank ef

a Полярный геофизический институт
Мурманская обл, Апатиты, Россия

b Институт земного магнетизма, ионосферы и распространения радиоволн им. Н.В. Пушкова РАН
Москва, Троицк, Россия

c Научно-исследовательский институт ядерной физики им. Д.В. Скобельцына МГУ
, Россия

d Институт космических исследований РАН
Москва, Россия

e Center for Space Plasma and Aeronomic Research (CSPAR), University of Alabama in Huntsville
Huntsville, USA

f Department of Space Science, University of Alabama in Huntsville
Huntsville, USA

* E-mail: mingalev_o@pgia.ru
** E-mail: habarova@izmiran.ru

Поступила в редакцию 21.12.2017
После доработки 25.07.2018
Принята к публикации 01.03.2018

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

Аннотация

Пересечения гелиосферного токового слоя (ГТС) на орбите Земли часто ассоциируются с наблюдением анизотропных пучков высокоэнергичных протонов, ускоренных до энергий от сотен кэВ до нескольких МэВ и выше. Недавно обнаружилась связь между этим явлением и присутствием маломасштабных магнитных островов (МО) вблизи пересоединяющихся токовых слоев. В работе предложен механизм набора энергии предускоренными протонами при колебаниях множественных магнитных островов внутри складки пересоединяющегося ГТС. Разработана модель электромагнитного поля колеблющегося трехмерного МО с характерным размером ~0.001 а. е., на который налетают протоны, ускоренные ранее в процессе магнитного пересоединения до энергий порядка от кэВ до десятков кэВ. Численные расчеты продемонстрировали, что возникающие продольные индукционные электрические поля способны дополнительно ускорять падающие на МО протоны. Показано, что существует локальная “ускорительная” область внутри острова, в которой приобретение энергии частицами наиболее эффективно, в результате чего их средние энергии вылета составляют от сотен кэВ до 2 МэВ и выше. Вне этой области ускорение частиц практически отсутствует. Показано, что набранные частицами энергии существенно зависят от начальных фаз и места попадания в МО, но слабо зависят от величин начальных энергий. Таким образом, низкоэнергичные частицы могут ускоряться эффективнее, чем высокоэнергичные, а на вылете из МО всеми частицами может достигаться общая предельная энергия. Обнаружено также, что направления скоростей вылета имеют сильную пространственную анизотропию. Полученные результаты согласуются с наблюдениями в плазме солнечного ветра.

Ключевые слова: гелиосферный токовый слой, магнитные острова, ускорение заряженных частиц

ВВЕДЕНИЕ

Гелиосферный токовый слой (ГТС) – это самое долгоживущее магнитное образование в гелиосфере. ГТС представляет собой условное продолжение солнечного магнитного экватора и является сравнительно тонкой токовой структурой толщиной порядка ~104 км, разделяющей магнитные сектора – крупномасштабные области с полями противоположного направления (от Солнца/к Солнцу). ГТС окружен более толстым плазменным слоем, заполненным вторичными токовыми слоями, магнитными островами и другими продуктами магнитного пересоединения (Bemporad, 2008; Cartwright и др., 2008; Khabarova и др., 2015a; 2015b; 2016). ГТС динамически реагирует на прохождение потоков солнечного ветра различной природы, вследствие чего даже в относительно спокойной гелиосфере его форма разнообразна и варьируется от крупномасштабной волны до мелкой складчатости (Khabarova и др., 2015a; 2015b; 2016; 2017). Данные наблюдений указывают на то, что в результате активных процессов в солнечной короне вблизи ГТС образуются магнитные острова (MO), называемые также магнитными пузырями, плазмоидами или баблами, которые затем растут, эволюционируют и движутся вместе с солнечным ветром вдоль складок ГТС на большие расстояния (Bemporad, 2008; Cartwright и др., 2008; Khabarova и др., 2015a; 2015b; 2016). На орбите Земли их размер в среднем имеет порядок 0.001–0.01 а. е., хотя разброс может быть большим (Cartwright и др., 2008). В то же время часть островов образуется непосредственно в солнечном ветре в результате магнитного пересоединения на ГТС или же как следствие развития различных неустойчивостей, которым подвержены токовые слои (Eastwood и др., 2002; Drake и др., 2006; Markidis и др., 2013; Greco и др., 2016; Zheng, Hu, 2018).

На рис. 1а показаны относительно мелкомасштабные складки гелиосферного токового слоя, профиль которого восстановлен по данным межпланетных сцинтилляций (см. Khabarova и др., 2015a; 2015b; 2016 и http://smei.ucsd.Edu/new_smei/data&images/data&images.html). В каждой такой складке может находиться от одного до нескольких магнитных островов, как проиллюстрировано на рис. 1б. Система “ГТС–магнитные острова” представлена на рис. 1б в вертикальном срезе. Тройная линия ГТС отражает его расслоение при пересоединении, что соответствует наблюдениям (Khabarova и др., 2015a; 2015b; 2016). Плиссированная форма ГТС часто возникает после прохождения мощных выбросов корональных масс, что провоцирует интенсификацию магнитного пересоединения в ГТС и возбуждение разнообразных МГД-волн, бегущих по ГТС и колеблющих магнитные острова. В данной работе мы рассмотрим простейший случай одного магнитного острова, зажатого в складке пересоединяющегося токового слоя. В работах (Khabarova и др., 2015a; 2015b; 2016) показано, что эти складки подобны токамаку, удерживающему МО и обеспечивающему эффективность любых механизмов ускорения. То же происходит при любом магнитном экранировании объема, занятого магнитными островами, например, при возникновении магнитных полостей между токовыми слоями на фронтах различных потоков солнечного ветра и ГТС. Присутствие МО в любой естественной магнитной ловушке, ограниченной сильными токовыми слоями, часто сопровождается нетипичными увеличениями потоков заряженных частиц, энергии которых могут достигать нескольких МэВ (обычно, 100 кэВ–1 МэВ) (Khabarova и др., 2015a; 2015b; 2016). Нетипичность этих всплесков потоков энергичных частиц заключается в том, что ранее предполагалось, что процесс ускорения частиц до таких энергий может происходить лишь далеко от земной орбиты. Считалось, будто частицы могут ускоряться до таких энергий либо на Солнце во время вспышек, либо на ударной волне выброса корональных масс (также вблизи Солнца), либо на 2–3 а. е. на фронтах ударных волн, сформированных потоками из корональных дыр. Однако наблюдения указывают на локальность ускорения энергичных частиц, поскольку области, в которых обнаруживаются ускоренные частицы, перемещаются вместе с окружающим солнечным ветром, что отчетливо трассируется несколькими пространственно разнесенными космическими аппаратами.

Рис. 1.

Складки ГТС и магнитные острова. (а) – Профиль ГТС, восстановленный по данным межпланетных сцинтилляций STELab. Орбита и положение Земли обозначены схематически. (б) – Вертикальный срез ГТС. Из-за плиссированной формы ГТС, магнитные острова размерами 0.01–0.001 а. е. оказываются захваченными внутри складок ГТС. Радиальная компонента магнитного поля вне системы “ГТС–магнитные острова” разнонаправлена (обозначено знаками +/–).

Возможность стохастического ускорения энергичных частиц в МО в связи с динамическими процессами была показана в результате моделирований (Matthaeus и др., 1984; Drake и др., 2006; Oka и др., 2010; Bian, Kontar, 2013; Zhou и др., 2015), а теория ускорения заряженных частиц в МО, испытывающих слияние или сжатие, была развита в работах (Zank и др., 2014; 2015a; 2015b; le Roux и др., 2015; 2016). В Zank и др. (2015a; 2015b) показано, что теоретические предсказания комбинированного ускорения частиц на ударных волнах с ускорением в МО за фронтом ударной волны соответствуют наблюдениям космического аппарата Voyager-2. В предложенном комбинированном механизме первичное ускорение обеспечивается фронтом ударной волны. В то же время, в механизме слияния или сжатия МО (Zank и др., 2014; le Roux и др., 2015; 2016) нет ограничений на источник предварительно ускоренных частиц, т.е. предускоренные частицы могут иметь различную природу.

При сравнении предсказаний теории (Zank и др., 2015b) с наблюдениями был впервые обнаружен эффект двух популяций ускоренных частиц, ассоциирующихся с магнитными островами: ниже определенной пороговой энергии относительное ускорение тем сильнее, чем меньшую энергию имели изначально поступившие в систему частицы, а выше пороговой энергии – наоборот, интенсивнее ускоряются более энергичные частицы. Механизм (Zank и др., 2015b) объясняет поведение частиц из последней популяции. Чем объясняется выраженная инверсия ускорения, и как образуется первая популяция – не было известно. Очевидно, здесь может играть роль иной механизм, не связанный со сжатием или слиянием МО.

Рассмотрим возможные источники дополнительного ускорения частиц в МО в солнечном ветре. Если ударные волны отсутствуют, предускоренные частицы с типичными энергиями от 1 кэВ до нескольких десятков кэВ поступают в области, занятые магнитными островами, преимущественно из квазирегулярно пересоединяющегося ГТС (Zharkova, Khabarova, 2012; 2015; Khabarova и др., 2017). При этом характерная температура тепловых протонов солнечного ветра на несколько порядков меньше и составляет в среднем около 10 эВ (может варьироваться в пределах 1–100 эВ). Одним из возможных механизмов дальнейшего ускорения предускоренных протонов рядом с ГТС представляется ускорение продольным электрическим полем в системе достаточно крупномасштабных колеблющихся магнитных островов, которые движутся вместе с солнечным ветром вдоль складки ГТС.

Коллективные эффекты ускорения частиц, возникающие при сжатии множества МО подробно изучены в (Zank и др., 2014; 2015a; 2015b; Zhou и др., 2015; le Roux и др., 2015; 2016), при этом механизм ускорения является стохастическим. Области сжатия могут наблюдаться между взаимодействующими потоками разной природы, а также между гелиосферным токовым слоем и любым набегающим потоком. Однако в условиях спокойного солнечного ветра наиболее распространенный и естественный механизм, вызывающий сжатие/разрежение в МО – это магнитогидродинамические волны с характерным периодом от минут до часов, которые распространяются в азимутальном направлении вдоль ГТС (Musielak, Suess, 1988; Wang и др., 1988; Ruderman, 1990; Dai и др., 2014). Численные эксперименты показывают, что наличие изгибов магнитного экватора Солнца на уровне фотосферы приводит к усилению этого эффекта и возникновению выраженных складок поверхности ГТС на расстояния порядка 1 а. е. (Merkin и др., 2011).

В связи с этим естественной кажется цепь событий, описанная ниже. Поверхностные волны создают колебания складок ГТС и стимулируют процессы пересоединения в складках ГТС. Предускоренные магнитным пересоединением частицы падают на магнитные острова. При этом МО колеблются вместе с колебанием складок ГТС, что порождает возникновение в них индукционного электрического поля. Для произвольно выбранного острова размера порядка $L \sim {{10}^{8}}$ м с периодом колебаний порядка 10 мин и амплитудой колебаний магнитного поля 5 нТл на расстоянии 1 а. е. наши оценки показывают, что продольная часть индукционного электрического поля может составлять ${{E}_{{v\parallel }}} \sim {{10}^{{ - 5}}} - {{10}^{{ - 4}}}$ В/м, в то время как потенциальная часть электрического поля за счет разделения зарядов ${{E}_{p}}\,$ может иметь значения на 1–2 порядка ниже: не более чем 10–6 В/м. Таким образом, достаточно сильное продольное индукционное электрическое поле не может быть экранировано полем за счет разделения заряда. В результате последовательного ускорения частиц на множестве рядом расположенных магнитных островов конечная энергия частиц может увеличиваться на несколько порядков.

Действительно, как показано в (Zank и др., 2014; 2015a; 2015b; Zhou и др., 2015; le Roux и др., 2015; 2016), коллективные эффекты играют решающую роль в процессе ускорения в обширных областях, заполненных динамическими магнитными островами, поскольку приводят к доминированию ускорения Ферми второго рода над всеми возможными потерями энергии в отдельных островах (le Roux и др., 2016) .

Моделирование этого крупномасштабного и сложного явления во всей полноте требует огромного объема вычислительных ресурсов. Поэтому нами предложен следующий подход. Вначале строится аналитическая модель электромагнитного поля в одном периодически колеблющемся МО с определенным набором входных параметров. В разработанной ниже модели электрическое поле формально считается чисто индукционным. В качестве второго этапа с помощью этой модели выбирается ряд вариантов полей модельного острова, каждый из которых согласуется с известными данными для МО на 1 а. е. Далее в этих полях рассчитываются траектории протонов для достаточно представительного набора начальных данных (начальное время/фаза колебания, начальная точка запуска частиц и начальная скорость).

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

ГЕОМЕТРИЯ СИСТЕМЫ И СЦЕНАРИЙ КОЛЕБАНИЙ ОСТРОВА

Рассмотрим периодически колеблющийся МО, имеющий форму эллипсоида с полуосями, колеблющимися вокруг своих средних значений. Заметим, что модель представляет собой достаточно упрощенную картину устройства МО, которые могут образовываться в результате пересоединения магнитных полей многочисленных токовых слоев, окружающих ГТС. Так, в модели остров считается автономным и не связанным с общей токовой системой внутри гелиосферного плазменного слоя, что, конечно же, является достаточно сильным упрощением. Однако разработанная модель не является самосогласованной, основное внимание в ней уделяется только лишь процессам ускорения частиц в той области острова, где возникающие электрические поля максимальны. Структура прилежащих к МО областей, равно как и процессы попадания частиц внутрь острова не рассматриваются. Будем использовать сферическую систему координат $\left( {r,\alpha ,\beta } \right)$ с началом отсчета в центре Солнца, физический базис которой обозначим через $\left\{ {{{{\mathbf{h}}}_{{{\kern 1pt} r}}},{{{\mathbf{h}}}_{\alpha }},{{{\mathbf{h}}}_{\beta }}} \right\},$ а физические компоненты магнитного поля в этой системе обозначим через ${{B}_{r}},{{B}_{\alpha }},{{B}_{\beta }}.$ Также будем использовать солнечно-магнитосферную декартову систему координат GSM, декартов базис которой обозначим через $\left\{ {{{{\mathbf{e}}}_{x}},{{{\mathbf{e}}}_{y}},{{{\mathbf{e}}}_{z}}} \right\},$ а компоненты магнитного поля в этой системе обозначим через ${{B}_{x}},{{B}_{y}},{{B}_{z}}.$ На расстоянии 1 а. е. от Солнца кривизной силовых линий можно пренебречь и считать, что

$\begin{gathered} {{{\mathbf{e}}}_{x}} = - {{{\mathbf{h}}}_{r}},\,\,\,\,{{{\mathbf{e}}}_{y}} = - {{{\mathbf{h}}}_{\alpha }},\,\,\,\,{{{\mathbf{e}}}_{z}} = {{{\mathbf{h}}}_{\beta }},\,\,\,\,{{B}_{x}} = - {{B}_{r}}, \\ {{B}_{y}} = - {{B}_{\alpha }},\,\,\,\,{{B}_{z}} = {{B}_{\beta }}. \\ \end{gathered} $

Будем считать, что крупномасштабная складка ГТС направлена вдоль ${{{\mathbf{e}}}_{z}},$ т.е. перпендикулярно плоскости эклиптики, что соответствует наблюдениям. В модели мы полагаем, что в складке находится только один остров, что также часто встречается вблизи гелиосферного токового слоя.

Будем считать, что в невозмущенном “стационарном” случае внутри этой складки фоновые поля солнечного ветра однородны и стационарны, и в системе координат GSM имеют вид

(1)
${{{\mathbf{B}}}^{{(sw)}}} = {{B}_{{x0}}}{{{\mathbf{e}}}_{x}} + {{B}_{{y0}}}{{{\mathbf{e}}}_{y}},\,\,\,\,{{{\mathbf{E}}}^{{(sw)}}} = {{E}_{{z0}}}{{{\mathbf{e}}}_{z}}.$

Тогда скорость солнечного ветра (используем систему СИ) в системе координат GSM имеет вид

(2)
${{{\mathbf{V}}}^{{(sw)}}} = {{{\mathbf{v}}}_{E}} = \frac{{\left[ {{{{\mathbf{E}}}^{{(sw)}}} \times {{{\mathbf{B}}}^{{(sw)}}}} \right]}}{{{{{\left| {{{{\mathbf{B}}}^{{(sw)}}}} \right|}}^{2}}}} = \frac{{{{E}_{{z0}}}\left( { - {{B}_{{y0}}}{{{\mathbf{e}}}_{x}} + {{B}_{{x0}}}{{{\mathbf{e}}}_{y}}} \right)}}{{{{{\left| {{{B}_{{x0}}}} \right|}}^{2}} + {{{\left| {{{B}_{{y0}}}} \right|}}^{2}}}}.$

Отметим, что, согласно данным спутниковых измерений, величина скорости солнечного ветра ${{V}^{{(sw)}}} = \left| {{\kern 1pt} {{{\mathbf{V}}}^{{(sw)}}}} \right|$ может изменяться в пределах примерно 250–1200 км/с, а ее среднее значение составляет около 400 км/с, то есть всегда выполнено условие ${{{{V}^{{(sw)}}}} \mathord{\left/ {\vphantom {{{{V}^{{(sw)}}}} c}} \right. \kern-0em} c} \ll 1.$ Из последнего условия и формулы (2) вытекает условие

$\frac{{\left| {{{{\mathbf{E}}}^{{(sw)}}}} \right|}}{{c\left| {{{{\mathbf{B}}}^{{(sw)}}}} \right|{\kern 1pt} }} = \frac{{\left| {{{E}_{{z0}}}} \right|}}{{c\sqrt {{\kern 1pt} \left| {{{B}_{{x0}}}} \right|{{{\kern 1pt} }^{2}} + {{{\left| {{{B}_{{y0}}}} \right|}}^{2}}} }} \ll 1.$

Далее мы перейдем в систему координат $K{\text{'}},$ оси которой параллельны осям системы координат GSM и которая движется вместе с солнечным ветром со скоростью ${{{\mathbf{V}}}^{{(sw)}}}.$ Согласно формулам преобразования Лоренца для полей в этой системе отсчета электрическое поле отсутствует: ${{{\mathbf{E}}}^{{(sw)}}}{\text{'}} = 0,$ а магнитное имеет вид

$\begin{gathered} {{{\mathbf{B}}}^{{(sw)}}}{\text{'}} = \frac{1}{{\sqrt {1 - {{{\left( {{{{{V}^{{(sw)}}}} \mathord{\left/ {\vphantom {{{{V}^{{(sw)}}}} c}} \right. \kern-0em} c}} \right)}}^{2}}} }} \times \\ \times \,\,\left( {1 - {{{\left( {{{\left| {{{{\mathbf{E}}}^{{(sw)}}}} \right|} \mathord{\left/ {\vphantom {{\left| {{{{\mathbf{E}}}^{{(sw)}}}} \right|} {\left( {c\left| {{{{\mathbf{B}}}^{{(sw)}}}} \right|} \right)}}} \right. \kern-0em} {\left( {c\left| {{{{\mathbf{B}}}^{{(sw)}}}} \right|} \right)}}} \right)}}^{2}}} \right){{{\mathbf{B}}}^{{(sw)}}} \approx {{{\mathbf{B}}}^{{(sw)}}}, \\ \end{gathered} $
т.е. практически не отличается от поля в системе координат GSM. Будем рассматривать в системе отсчета $K{\text{'}}$ МО, движущийся вместе с солнечным ветром. В этой системе отсчета центр острова неподвижен, и в невозмущенном “стационарном” случае электрическое поле отсутствует, а магнитное поле имеет вид

(3)
${\mathbf{B}}{\kern 1pt} \left( {\mathbf{x}} \right) = {{{\mathbf{B}}}^{{(sw)}}} + {{{\mathbf{B}}}^{{(i)}}}\left( {\mathbf{x}} \right).$

Здесь и далее штрихи у полей опускаем. Пространственно неоднородная часть магнитного поля ${{{\mathbf{B}}}^{{{\kern 1pt} (i)}}}\left( {\mathbf{x}} \right)$ создается находящейся внутри острова замкнутой токовой системой (Григоренко, 2016) с током ${{{\mathbf{j}}}^{{{\kern 1pt} (i)}}}\left( {\mathbf{x}} \right),$ удовлетворяющим уравнению Ампера

${\text{rot}}{\kern 1pt} {{{\mathbf{B}}}^{{(i)}}}\left( {\mathbf{x}} \right) = {{{\mu }}_{0}}{{{\mathbf{j}}}^{{(i)}}}\left( {\mathbf{x}} \right).$

Геометрия задачи схематически показана на рис. 2. Панель (а) демонстрирует положение острова в складке гофрированного ГТС. На панели (б) изображен общий вид острова и система координат. Ось Х направлена к Солнцу, против направления распространения солнечного ветра. На панели (в) показано поперечное сечение острова с периодом колебаний ${\Theta }$плоскостью $\left\{ {x = 0} \right\}$ в моменты времени $t = 0,{{\Theta } \mathord{\left/ {\vphantom {{\Theta } 4}} \right. \kern-0em} 4},3{{\Theta } \mathord{\left/ {\vphantom {{\Theta } 4}} \right. \kern-0em} 4},$ соответствующие начальному, максимально сжатому и максимально широкому положениям. Из-за наличия доминирующего направления движения солнечного ветра остров будет вытянут вдоль координаты (Х), а колебания будут происходить в поперечном направлении (YZ). В общем случае бабл можно представить зажатым токовым слоем с двух сторон, предполагая, что сверху и сбоку он соседствует с другими магнитными островами (см. геометрию на рис. 2а). В таком случае ток течет по ГТС с двух сторон бабла в противоположную сторону, как показано на рис. 2а и 2в.

Рис. 2.

Топология численного эксперимента. (а) – Схематическое положение магнитного острова в складке гелиосферного токового слоя (ГТС). (б) – Общий вид острова и система координат. (в) – Сечение острова плоскостью $\left\{ {x = 0} \right\}$ при различных фазах колебаний $\tau = 0\,$ – сплошная линия, $\tau = {{\Theta } \mathord{\left/ {\vphantom {{\Theta } 4}} \right. \kern-0em} 4}\,$ – штрих-пунктирная линия, $\tau = {{3{\Theta }} \mathord{\left/ {\vphantom {{3{\Theta }} 4}} \right. \kern-0em} 4}$ – пунктирная линия.

На рис. 3 показано распределение модуля нормированного магнитного поля острова ${{\left| {{{{\mathbf{B}}}^{{{\text{(}}i{\text{)}}}}}\left( {x = 0,y,z,t} \right)} \right|} \mathord{\left/ {\vphantom {{\left| {{{{\mathbf{B}}}^{{{\text{(}}i{\text{)}}}}}\left( {x = 0,y,z,t} \right)} \right|} {{{B}_{{x0}}}}}} \right. \kern-0em} {{{B}_{{x0}}}}}$ в плоскости $\left\{ {x = 0} \right\}$ при фазах колебаний $\tau = 0,{{\Theta } \mathord{\left/ {\vphantom {{\Theta } 8}} \right. \kern-0em} 8},{{\Theta } \mathord{\left/ {\vphantom {{\Theta } 4}} \right. \kern-0em} 4}.$

Рис. 3.

Модуль нормированного магнитного поля острова ${{\left| {{{{\mathbf{B}}}^{{{\text{(}}i{\text{)}}}}}\left( {x = 0,y,z,t} \right)} \right|} \mathord{\left/ {\vphantom {{\left| {{{{\mathbf{B}}}^{{{\text{(}}i{\text{)}}}}}\left( {x = 0,y,z,t} \right)} \right|} {{{B}_{{x0}}}}}} \right. \kern-0em} {{{B}_{{x0}}}}}$ в плоскости $\left\{ {x = 0} \right\}$ при различных фазах колебаний: (1) $\tau = 0,$ 2) $\tau = {\Theta \mathord{\left/ {\vphantom {\Theta 8}} \right. \kern-0em} 8},$ (3) $\tau = {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4}.$

Рассмотрим кратко систему уравнений, которая описывает медленные (с характерным временем изменения порядка 10 мин) крупномасштабные процессы в бесстолкновительной плазме. С ее помощью проанализируем возможный сценарий колебаний магнитного острова в солнечном ветре вблизи орбиты Земли. Как известно, в солнечном ветре имеется два вида основных ионов: протоны Н+ и альфа-частицы Не2+, причем характерное значение доли альфа-частиц ${{{{n}_{\alpha }}} \mathord{\left/ {\vphantom {{{{n}_{\alpha }}} {{{n}_{e}}}}} \right. \kern-0em} {{{n}_{e}}}}$ составляет 5%, но вблизи ГТС доля альфа-частиц существенно понижена, поэтому ими в приближенных оценках можно пренебречь.

На малых пространственных масштабах порядка дебаевского расстояния электронов бесстолкновительная плазма описывается нерелятивистской системой уравнений Власова–Максвелла. Будем обозначать через $\left( {{\mathbf{a}},{\mathbf{b}}} \right)$ и $\left[ {{\mathbf{a}} \times {\mathbf{b}}} \right]$ скалярное и векторное произведения векторов в пространстве ${{\mathbb{R}}^{3}},$ а через ${\mathbf{a}} \otimes {\mathbf{b}}$ – образованный этими векторами диадный тензор с декартовыми компонентами ${{\left( {{\mathbf{a}} \otimes {\mathbf{b}}} \right)}_{{k,l}}} = {{a}_{k}}{{b}_{l}}.$ Для каждой из плазменных компонент (протоны p+ и электроны e) $\beta = p,e$ через ${{f}_{\beta }}\left( {t,{\mathbf{x}},{\mathbf{v}}} \right)$ обозначим функцию распределения, зависящую от времени t, пространственной координаты ${\mathbf{x}} = {{\left( {{{x}_{1}},{{x}_{2}},{{x}_{3}}} \right)}^{T}} \in {{\mathbb{R}}^{3}}$ и скорости ${\mathbf{v}} = {{\left( {{{v}_{1}},{{v}_{2}},{{v}_{3}}} \right)}^{T}} \in {{\mathbb{R}}^{3}},$ через ${{e}_{{{\kern 1pt} \beta }}}$ обозначим заряд частицы (для электронов ${{e}_{e}} = - e,$ где e – заряд протона), через ${{m}_{\beta }}$ – массу частицы, через ${{n}_{{{\kern 1pt} \beta }}}\left( {{\mathbf{x}},t} \right)$ – концентрацию, через ${{{\mathbf{j}}}_{{{\kern 1pt} \beta }}}\left( {{\mathbf{x}},t} \right)$ – плотность тока. Через ${{{\varepsilon }}_{0}}$ и ${{{\mu }}_{0}}$ обозначим электрическую и магнитную постоянные, а через $c = {1 \mathord{\left/ {\vphantom {1 {\sqrt {{{{\varepsilon }}_{0}}{{{\mu }}_{0}}{\kern 1pt} } }}} \right. \kern-0em} {\sqrt {{{{\varepsilon }}_{0}}{{{\mu }}_{0}}{\kern 1pt} } }}$ – скорость света в вакууме. В системе СИ систему уравнений Власова–Максвелла можно представить в следующем виде:

(4)
$\begin{gathered} \frac{{\partial {{f}_{\beta }}}}{{\partial {\kern 1pt} t}} + \left( {{\mathbf{v}},\frac{{\partial {{f}_{\beta }}}}{{\partial {\mathbf{x}}}}} \right) + \frac{{{{e}_{{{\kern 1pt} \beta }}}}}{{{\kern 1pt} {{m}_{\beta }}}}\left( {{\mathbf{E}}\left( {{\mathbf{x}},t} \right) + \left[ {{\mathbf{v}} \times {\mathbf{B}}\left( {x,t} \right)} \right]{\kern 1pt} ,\frac{{\partial {{f}_{\beta }}}}{{\partial {\mathbf{v}}}}} \right) = 0, \\ \beta = p,e, \\ \end{gathered} $
(5)
$\begin{gathered} {{n}_{{{\kern 1pt} \beta }}}\left( {{\mathbf{x}},t} \right) = \int\limits_{{{\mathbb{R}}^{3}}} {{{f}_{\beta }}\left( {t,{\mathbf{x}},{\mathbf{v}}} \right){{d}^{3}}{\mathbf{v}}} , \\ \rho \left( {{\mathbf{x}},t} \right) = e{\kern 1pt} \left( {{{n}_{p}}\left( {{\mathbf{x}},t} \right) - {{n}_{e}}\left( {{\mathbf{x}},t} \right)} \right), \\ \end{gathered} $
(6)
$\begin{gathered} {{{\mathbf{j}}}_{{{\kern 1pt} \beta }}}\left( {{\mathbf{x}},t} \right) = {{e}_{\beta }}\int\limits_{{{\mathbb{R}}^{3}}} {{\mathbf{v}}{{f}_{\beta }}\left( {t,{\mathbf{x}},{\mathbf{v}}} \right){{d}^{3}}{\mathbf{v}}} , \\ {\mathbf{j}}\left( {{\mathbf{x}},t} \right) = {{{\mathbf{j}}}_{p}}\left( {{\mathbf{x}},t} \right) + {{{\mathbf{j}}}_{e}}\left( {{\mathbf{x}},t} \right), \\ \end{gathered} $
(7)
${\text{div}}{\mathbf{B}}\left( {{\mathbf{x}},t} \right) = 0,$
(8)
$\frac{{\partial {\mathbf{B}}\left( {{\mathbf{x}},t} \right)}}{{\partial {\kern 1pt} t}} = - {\text{rot}}{\mathbf{E}}\left( {{\mathbf{x}},t} \right),$
(9)
${\text{div}}{\mathbf{E}}\left( {{\mathbf{x}},t} \right) = {{\rho \left( {{\mathbf{x}},t} \right)} \mathord{\left/ {\vphantom {{\rho \left( {{\mathbf{x}},t} \right)} {{{\varepsilon }_{0}}}}} \right. \kern-0em} {{{\varepsilon }_{0}}}},$
(10)
${{\varepsilon }_{0}}\frac{{\partial {\mathbf{E}}\left( {{\mathbf{x}},t} \right)}}{{\partial t}} = \frac{1}{{{{\mu }_{0}}}}{\text{rot}}{\kern 1pt} {\mathbf{B}}\left( {{\mathbf{x}},t} \right) - {\mathbf{j}}\left( {{\mathbf{x}},t} \right){\kern 1pt} .$

Здесь ${\mathbf{E}}\left( {{\mathbf{x}},t} \right)$ – вектор напряженности электрического поля, ${\mathbf{B}}{\kern 1pt} \left( {{\mathbf{x}},t} \right)$ – вектор магнитной индукции, ${{{\rho }}_{i}}\left( {{\mathbf{x}},t} \right)$ и ${{{\mathbf{j}}}_{{{\kern 1pt} i}}}\left( {{\mathbf{x}},t} \right)$ – суммарные плотности заряда и тока ионов, $\rho \left( {{\mathbf{x}},t} \right)$ и ${\mathbf{j}}\left( {{\mathbf{x}},t} \right)$ – полные плотности заряда и тока. Из уравнений Власова (4) вытекают следующие гидродинамические уравнения потока импульса для каждой компоненты плазмы $\beta = p,e\,:$

(11)
$\begin{gathered} \frac{{\partial {{{\mathbf{j}}}_{\beta }}}}{{\partial t}} = \frac{{{{e}_{\beta }}}}{{{{m}_{\beta }}}}\left( {{{e}_{\beta }}{{n}_{\beta }}{\mathbf{E}} + \left[ {{{{\mathbf{j}}}_{\beta }} \times {\mathbf{B}}} \right]{\kern 1pt} - {\text{div}}{{{{\mathbf{\hat {\Pi }}}}}_{\beta }}{\kern 1pt} } \right), \\ {{{{\mathbf{\hat {\Pi }}}}}_{\beta }}\left( {{\mathbf{x}},t} \right) = {{m}_{\beta }}\int\limits_{{{\mathbb{R}}^{3}}} {{\mathbf{v}} \otimes {\mathbf{v}}{\kern 1pt} {{f}_{\beta }}\left( {t,{\mathbf{x}},{\mathbf{v}}} \right){{d}^{3}}{\mathbf{v}}} , \\ \end{gathered} $
где ${{{\mathbf{\hat {\Pi }}}}_{\beta }}\left( {x,t} \right)$ – полный тензор напряжений для частиц сорта β, а ${\text{div}}{{{\mathbf{\hat {\Pi }}}}_{\beta }}\left( {{\mathbf{x}},t} \right)$ – вектор дивергенции этого тензора. Декартовы компоненты ${{\Pi }_{{\beta ,k,l}}}\left( {{\mathbf{x}},t} \right)$ тензора ${{{\mathbf{\hat {\Pi }}}}_{\beta }}\left( {{\mathbf{x}},t} \right)$ и компоненты ${{\left( {{\text{div}}{{{{\mathbf{\hat {\Pi }}}}}_{\beta }}\left( {{\mathbf{x}},t} \right)} \right)}_{k}}$ вектора ${\text{div}}{{{\mathbf{\hat {\Pi }}}}_{\beta }}\left( {{\mathbf{x}},t} \right)$ определяются формулами

$\begin{gathered} {{\Pi }_{{\beta ,k,l}}}\left( {{\mathbf{x}},t} \right) = {{m}_{\beta }}\int\limits_{{{\mathbb{R}}^{3}}} {{{v}_{k}}{{v}_{l}}{{f}_{\beta }}\left( {t,{\mathbf{x}},{\mathbf{v}}} \right){{d}^{3}}{\mathbf{v}}} , \\ {{\left( {{\text{div}}{{{{\mathbf{\hat {\Pi }}}}}_{\beta }}\left( {{\mathbf{x}},t} \right)} \right)}_{k}} = \sum\limits_{l = 1}^3 {\frac{{\partial {{\Pi }_{{\beta ,k,l}}}\left( {{\mathbf{x}},t} \right)}}{{\partial {{x}_{l}}}}} . \\ \end{gathered} $

Отметим, что полный тензор напряжений ${{{\mathbf{\hat {\Pi }}}}_{\beta }}\left( {{\mathbf{x}},t} \right)$ принято представлять в виде суммы тензора инерции ${{m}_{{{\kern 1pt} \beta }}}{{n}_{\beta }}{{{\mathbf{u}}}_{\beta }} \otimes {{{\mathbf{u}}}_{\beta }}$ и тензора давления ${{{\mathbf{\hat {{\rm P}}}}}_{\beta }}\,:$

(12)
${{{\mathbf{\hat {\Pi }}}}_{\beta }}\left( {{\mathbf{x}},t} \right) = {{m}_{\beta }}{{n}_{\beta }}\left( {{\mathbf{x}},t} \right){{{\mathbf{u}}}_{\beta }}\left( {{\mathbf{x}},t} \right) \otimes {{{\mathbf{u}}}_{\beta }}\left( {{\mathbf{x}},t} \right) + {{{\mathbf{\hat {{\rm P}}}}}_{\beta }}\left( {{\mathbf{x}},t} \right),$
где гидродинамическая скорость ${{{\mathbf{u}}}_{{{\kern 1pt} \beta }}}\left( {{\mathbf{x}},t} \right)$ и тензор давления определяются следующими формулами:

(13)
${{{\mathbf{u}}}_{{{\kern 1pt} \beta }}}\left( {{\mathbf{x}},t} \right) = \frac{{{{{\mathbf{j}}}_{\beta }}\left( {{\mathbf{x}},t} \right)}}{{{{e}_{{{\kern 1pt} \beta }}}{{n}_{\beta }}\left( {{\mathbf{x}},t} \right)}} = \frac{1}{{{{n}_{{{\kern 1pt} \beta }}}\left( {{\mathbf{x}},t} \right)}}\int\limits_{{{\mathbb{R}}^{3}}} {{\mathbf{v}}{\kern 1pt} {{f}_{{{\kern 1pt} \beta }}}\left( {t,{\mathbf{x}},{\mathbf{v}}} \right){{d}^{3}}{\mathbf{v}}} ,$
(14)
$\begin{gathered} {{{{\mathbf{\hat {{\rm P}}}}}}_{\beta }}\left( {{\mathbf{x}},t} \right) = {{m}_{\beta }}\int\limits_{{{\mathbb{R}}^{3}}} {\left( {{\mathbf{v}} - {{{\mathbf{u}}}_{\beta }}\left( {{\mathbf{x}},t} \right){\kern 1pt} } \right)} \otimes \\ \otimes \left( {{\mathbf{v}} - {{{\mathbf{u}}}_{\beta }}\left( {{\mathbf{x}},t} \right)} \right){{f}_{{{\kern 1pt} \beta }}}\left( {t,{\mathbf{x}},{\mathbf{v}}} \right){{d}^{3}}{\mathbf{v}}. \\ \end{gathered} $

Из формулы (12) следует, что дивергенция полного тензора напряжений ${\text{div}}{\kern 1pt} {{{\mathbf{\hat {\Pi }}}}_{\beta }}\left( {{\mathbf{x}},t} \right)$ является суммой объемной плотности силы инерции ${{m}_{{{\kern 1pt} \beta }}}{\text{div}}\left( {{{n}_{\beta }}{{{\mathbf{u}}}_{\beta }} \otimes {{{\mathbf{u}}}_{\beta }}} \right)$ и дивергенции тензора давления ${\text{div}}{{{\mathbf{\hat {{\rm P}}}}}_{\beta }}\,:$

(15)
${\text{div}}{\kern 1pt} {{{\mathbf{\hat {\Pi }}}}_{\beta }} = {{m}_{\beta }}{\text{div}}\left( {{\kern 1pt} {{n}_{\beta }}{{{\mathbf{u}}}_{\beta }} \otimes {{{\mathbf{u}}}_{\beta }}} \right) + {\text{div}}{{{\mathbf{\hat {{\rm P}}}}}_{\beta }}.$

Суммирование уравнений (11) с учетом (15) приводит к уравнению

(16)
$\left. \begin{gathered} \frac{{\partial {\mathbf{j}}}}{{\partial {\kern 1pt} t}} = \frac{{{{e}^{2}}{{n}_{p}}}}{{{{m}_{p}}}}\left( {{\mathbf{E}} + \left[ {{{{\mathbf{u}}}_{p}} \times {\mathbf{B}}} \right] - \frac{{{{m}_{p}}{\text{div}}\left( {{{n}_{p}}{{{\mathbf{u}}}_{p}} \otimes {{{\mathbf{u}}}_{p}}} \right)}}{{e{{n}_{p}}}} - \frac{{{\text{div}}{{{{\mathbf{\hat {{\rm P}}}}}}_{p}}}}{{e{{n}_{p}}}}} \right) + \hfill \\ + \,\,\frac{{{{e}^{2}}{{n}_{e}}}}{{{{m}_{e}}}}\left( {{\kern 1pt} {\mathbf{E}} + \left[ {{{{\mathbf{u}}}_{e}} \times {\mathbf{B}}} \right]{\kern 1pt} + \frac{{{{m}_{e}}{\text{div}}\left( {{{n}_{e}}{{{\mathbf{u}}}_{e}} \otimes {{{\mathbf{u}}}_{e}}} \right)}}{{e{{n}_{e}}}} + \frac{{{\text{div}}{{{{\mathbf{\hat {{\rm P}}}}}}_{e}}}}{{e{{n}_{e}}}}} \right), \hfill \\ \end{gathered} \right\}{\kern 1pt} $
которое часто называют обобщенным законом Ома. Для удобства дальнейшего анализа введем следующие обозначения:
(17)
$\left. \begin{gathered} {{\left( {\frac{{\partial {\mathbf{j}}}}{{\partial t}}} \right)}_{E}} = \left( {\frac{{{{e}^{2}}{{n}_{p}}}}{{{{m}_{p}}}} + \frac{{{{e}^{2}}{{n}_{e}}}}{{{{m}_{e}}}}} \right){\mathbf{E}},\,\,\,\,{{\left( {\frac{{\partial {\mathbf{j}}}}{{\partial t}}} \right)}_{B}} = \frac{{{{e}^{2}}{{n}_{e}}}}{{{{m}_{e}}}}\left[ {{{{\mathbf{u}}}_{e}} \times {\mathbf{B}}} \right] + \frac{{{\kern 1pt} {{e}^{2}}{{n}_{p}}}}{{{{m}_{p}}}}\left[ {{{{\mathbf{u}}}_{p}} \times {\mathbf{B}}} \right], \hfill \\ {{\left( {\frac{{\partial {\mathbf{j}}}}{{\partial t}}} \right)}_{{{\text{in}}}}} = - e{\text{div}}\left( {{{n}_{p}}{{{\mathbf{u}}}_{p}} \otimes {{{\mathbf{u}}}_{p}}} \right) + e{\text{div}}\left( {{{n}_{e}}{{{\mathbf{u}}}_{e}} \otimes {{{\mathbf{u}}}_{e}}} \right),\,\,\,\,{{\left( {\frac{{\partial {\mathbf{j}}}}{{\partial t}}} \right)}_{P}} = \frac{e}{{{{m}_{e}}}}{\text{div}}{{{{\mathbf{\hat {{\rm P}}}}}}_{e}} - \frac{{{\kern 1pt} e}}{{{{m}_{p}}}}{\text{div}}{{{{\mathbf{\hat {{\rm P}}}}}}_{p}}, \hfill \\ \end{gathered} \right\}$
с учетом которых уравнение (16) примет вид

(18)
$\frac{{\partial {\mathbf{j}}}}{{\partial {\kern 1pt} t}} = {{\left( {\frac{{\partial {\mathbf{j}}}}{{\partial t}}} \right)}_{E}} + {{\left( {\frac{{\partial {\mathbf{j}}}}{{\partial {\kern 1pt} t}}} \right)}_{B}} + {{\left( {\frac{{\partial {\mathbf{j}}}}{{\partial t}}} \right)}_{{{\text{in}}}}} + {{\left( {\frac{{\partial {\mathbf{j}}}}{{\partial t}}} \right)}_{P}}.$

Удобно использовать разложение электрического поля на потенциальную ${{{\mathbf{E}}}_{p}}\left( {{\mathbf{x}},t} \right)$ и соленоидальную ${{{\mathbf{E}}}_{v}}\left( {{\mathbf{x}},t} \right)$ части в рассматриваемой области моделирования Ω:

(19)
$\begin{gathered} {\mathbf{E}}{\kern 1pt} \left( {{\mathbf{x}},t} \right) = {{{\mathbf{E}}}_{{p}}}\left( {{\mathbf{x}},t} \right) + {{{\mathbf{E}}}_{v}}{\kern 1pt} \left( {{\mathbf{x}},t} \right),\,\,\,\,{{{\mathbf{E}}}_{p}}\left( {{\mathbf{x}},t} \right) = - \nabla {\varphi }\left( {{\mathbf{x}},t} \right){\text{,}} \\ {\text{div}}{{{\mathbf{E}}}_{v}}\left( {{\mathbf{x}},t} \right) = 0,\,\,\,\,{\mathbf{x}} \in {\Omega }. \\ \end{gathered} $

Такое разложение определено с точностью до градиента гармонической в области Ω функции:

$\begin{gathered} {\mathbf{E}}_{v}^{'}\left( {{\mathbf{x}},t} \right) = {{{\mathbf{E}}}_{v}}\left( {{\mathbf{x}},t} \right) + \nabla w\left( {{\mathbf{x}},t} \right),\,\,\,{\mathbf{E}}_{p}^{'}\left( {{\mathbf{x}},t} \right) = \\ = - \nabla {\varphi }\left( {{\mathbf{x}},t} \right) - \nabla w\left( {{\mathbf{x}},t} \right),\,\,\,\,\Delta w\left( {{\mathbf{x}},t} \right){\kern 1pt} \, = 0,\,\,\,\,{\mathbf{x}} \in {\Omega }{\kern 1pt} , \\ \end{gathered} $
где $\Delta {\kern 1pt} w = {\text{div}}\nabla w$ – оператор Лапласа. Для единственности такого разложения необходимо использовать какое-либо граничное условие, которое обычно вытекает из специфики каждой конкретной задачи.

Следует отметить, что магнитные острова, наблюдающиеся вблизи токовых слоев в гелиосфере, образуют сильно неоднородную плазменную систему, в которой магнитное поле и плотность колеблются в противофазе. Такая система может включать как подобласти с увеличенным по сравнению с окружающим солнечным ветром магнитным полем и пониженным давлением плазмы, так и аналогичные подобласти с ослабленным магнитным полем и повышенным давлением плазмы. Для того, чтобы показать, что создающиеся в системе колеблющихся островов индукционные электрические поля могут быть велики и не компенсируются амбиполярным электрическим полем разделения зарядов, оценим слагаемые в обобщенном законе Ома (16), используя известные по данным измерений на КА на расстоянии 1 а. е. от Солнца характерные значения амплитуды колебаний магнитного поля в острове $\Delta B \approx 5$ нТл, периода колебаний ${\Theta } \sim {\text{300}}$ с, размера острова ${{L}_{0}} \sim {{10}^{8}}$ м, концентрации ${{n}_{e}} \approx {{n}_{p}} \sim 1 - 10\,\,{\text{с }}{{{\text{м }}}^{{ - 3}}}$ = ${{10}^{6}} - {{10}^{7}}\,\,{{{\text{м }}}^{{ - 3}}},$ температуры электронов ${{T}_{e}} \approx 12$ эВ и протонов ${{T}_{p}} \approx 8$ эВ. Из условия квазинейтральности ${{n}_{e}} \approx {{n}_{p}}$ следует, что

(20)
${{\left( {\frac{{\partial {\mathbf{j}}}}{{\partial t}}} \right)}_{E}} = \left( {\frac{{{{e}^{2}}{{n}_{p}}}}{{{{m}_{p}}}} + \frac{{{{e}^{2}}{{n}_{e}}}}{{{{m}_{e}}}}} \right){\mathbf{E}} \approx \frac{{{{e}^{2}}{{n}_{e}}}}{{{{m}_{e}}}}\left( {1 + \frac{{{{m}_{e}}}}{{{{m}_{p}}}}} \right){\mathbf{E}} \approx \frac{{{{e}^{2}}{{n}_{e}}}}{{{{m}_{e}}}}{\mathbf{E}}.$

В рассматриваемых медленных процессах ток смещения в уравнении Максвелла (10) пренебрежимо мал, и оно переходит в уравнение Ампера

(21)
${\kern 1pt} {\text{rot}}{\mathbf{B}}\left( {{\mathbf{x}},t} \right) = {{\mu }_{0}}{\mathbf{j}}\left( {{\mathbf{x}},t} \right).$

Упрощенно можно считать, что магнитный остров образуется двумя токовыми системами с одинаковым по величине, но противоположно направленным дипольным моментом, и с характерным размером неоднородности порядка $L \sim 0.1{{L}_{0}} \sim {{10}^{7}}$ м. Из уравнения (21) вытекает оценка

(22)
${\mathbf{j}}\left( {{\mathbf{x}},t} \right) = \frac{{{\text{rot}}{\mathbf{B}}\left( {{\mathbf{x}},t} \right)}}{{{{\mu }_{0}}}} \sim \frac{{\Delta B}}{{{{\mu }_{0}}L}} \sim 5 \times {{10}^{{ - 9}}}\,\,\,\,{{\text{А }} \mathord{\left/ {\vphantom {{\text{А }} {{{{\text{м }}}^{2}}}}} \right. \kern-0em} {{{{\text{м }}}^{2}}}}.$

Дифференцирование по времени уравнения Ампера (21) и учет (22) позволяют получить оценку для производной тока по времени

(23)
$\frac{{\partial {\mathbf{j}}}}{{\partial {\kern 1pt} t}} = \frac{1}{{{{\mu }_{0}}}}{\text{rot}}\frac{{\partial {\mathbf{B}}}}{{\partial t}} \sim \frac{{\Delta B}}{{{{\mu }_{0}}L\Delta t}} \sim {{10}^{{ - 11}}}\,\,{{\text{А }} \mathord{\left/ {\vphantom {{\text{А }} {\left( {{{{\text{м }}}^{2}}{\text{с }}} \right)}}} \right. \kern-0em} {\left( {{{{\text{м }}}^{2}}{\text{с }}} \right)}}$

Из условия квазинейтральности ${{n}_{e}} \approx {{n}_{p}},$ формул (13) и оценки (22) вытекает, что

(24)
${{u}_{p}} - {{u}_{e}} \sim \frac{j}{{e{{n}_{e}}}} \sim 0.1 - 1\,\,{{{\text{к м }}} \mathord{\left/ {\vphantom {{{\text{к м }}} {\text{с }}}} \right. \kern-0em} {\text{с }}}.$

Эта оценка вместе условием квазинейтральности и оценкой ${{u}_{p}} \sim {{u}_{e}} \sim {{V}^{{(sw)}}} \sim 400$ км/с означают, что скорость протонов в движущейся системе координат составляет

(25)
${{u}_{p}} \approx {{u}_{e}} \sim {{V}^{{(sw)}}} \sim 100\,\,{{{\text{к м }}} \mathord{\left/ {\vphantom {{{\text{к м }}} {\text{с }}}} \right. \kern-0em} {\text{с }}}.$

Отсюда следует, что

(26)
${{\left( {\frac{{\partial {\mathbf{j}}}}{{\partial t}}} \right)}_{B}} \sim 1.5 \times {{10}^{{ - 4}}}\,\,{{\text{А }} \mathord{\left/ {\vphantom {{\text{А }} {({{{\text{м }}}^{2}}\,{\text{с }})}}} \right. \kern-0em} {({{{\text{м }}}^{2}}\,{\text{с }})}} \gg \frac{{\partial {\mathbf{j}}}}{{\partial {\kern 1pt} t}}.$

Из (24) и (25) следует, что третьи слагаемые в каждой скобке в (16) близки по величине и противоположны по знаку, то есть их сумма очень мала по сравнению с каждым из них:

(27)
$\left. \begin{gathered} {{\left( {\frac{{\partial {\mathbf{j}}}}{{\partial t}}} \right)}_{{{\text{in}}}}} = - e{\text{div}}\left( {{{n}_{p}}{{{\mathbf{u}}}_{p}} \otimes {{{\mathbf{u}}}_{p}}} \right) + e{\text{div}}\left( {{{n}_{e}}{{{\mathbf{u}}}_{e}} \otimes {{{\mathbf{u}}}_{e}}} \right) \sim \hfill \\ \sim \,\,\frac{{2{{u}_{e}}j}}{L} \sim {{10}^{{ - 13}}} - {{10}^{{ - 12}}}\,\,{{\text{А }} \mathord{\left/ {\vphantom {{\text{А }} {({{{\text{м }}}^{2}}\,{\text{с }})}}} \right. \kern-0em} {({{{\text{м }}}^{2}}\,{\text{с }})}} \sim \frac{{\partial {\mathbf{j}}}}{{\partial t}}. \hfill \\ \end{gathered} \right\}$

Если предположить, что в спокойном солнечном ветре и в магнитном острове величина магнитного поля $B \sim 2 - 10$ нТл, а температуры электронов и протонов составляют ${{T}_{e}} \approx 12$ эВ и ${{T}_{p}} \approx 8$ эВ, то их характерные гирорадиусы (${{r}_{{c{\beta }}}} = {{{{V}_{{T{\beta }}}}} \mathord{\left/ {\vphantom {{{{V}_{{T{\beta }}}}} {{{{\omega }}_{{c{\beta }}}}}}} \right. \kern-0em} {{{{\omega }}_{{c{\beta }}}}}}$) и гиропериоды (${{\theta }_{{c\beta }}} = {{2\pi } \mathord{\left/ {\vphantom {{2\pi } {{{\omega }_{{c\beta }}}}}} \right. \kern-0em} {{{\omega }_{{c\beta }}}}},$ где ${{{{\omega }_{{c\beta }}} = \left| {{{e}_{\beta }}} \right|B} \mathord{\left/ {\vphantom {{{{\omega }_{{c\beta }}} = \left| {{{e}_{\beta }}} \right|B} {{{m}_{\beta }}}}} \right. \kern-0em} {{{m}_{\beta }}}}$ – гирочастота, а ${{V}_{{T\alpha }}} = \sqrt {{{e{{T}_{\alpha }}} \mathord{\left/ {\vphantom {{e{{T}_{\alpha }}} {{{m}_{\alpha }}}}} \right. \kern-0em} {{{m}_{\alpha }}}}} $ – тепловая скорость частиц сорта β) имеют значения, соответственно, ${{r}_{{cp}}} \approx 30 - 100$км, ${{r}_{{ce}}} \approx 1 - 4$км и ${{\theta }_{{cp}}} \approx 2\pi - 10\pi $с, ${{\theta }_{{ce}}} \approx 0.003 - 0.02$ с. Отметим, что в этом случае тепловые электроны являются сильно замагниченными, и их тензор давления должен быть гиротропным: ${{{\mathbf{\hat {{\rm P}}}}}_{e}}\left( {{\mathbf{x}},t} \right)$ = ${{p}_{{e \bot }}}\left( {{\mathbf{x}},t} \right){\mathbf{\hat {I}}}$ + + $\left( {{{p}_{{e\parallel }}}\left( {{\mathbf{x}},t} \right) - {{p}_{e}}_{ \bot }\left( {{\mathbf{x}},t} \right){\kern 1pt} } \right){\mathbf{b}}\left( {{\mathbf{x}},t} \right) \otimes {\mathbf{b}}\left( {{\mathbf{x}},t} \right).$

Здесь ${{p}_{{e \bot }}}$и ${{p}_{{e\parallel }}}$– ортогональное и продольное давление электронов, ${\mathbf{\hat {I}}}$ – единичный тензор, ${\mathbf{b}}({\mathbf{x}},t) = {{{\mathbf{B}}({\mathbf{x}},t)} \mathord{\left/ {\vphantom {{{\mathbf{B}}({\mathbf{x}},t)} {\left| {{\mathbf{B}}({\mathbf{x}},t)} \right|}}} \right. \kern-0em} {\left| {{\mathbf{B}}({\mathbf{x}},t)} \right|}}$ – единичный вектор вдоль магнитного поля. Дивергенция этого тензора дается формулой

(28)
$\begin{gathered} {\text{div}}{{{{\mathbf{\hat {{\rm P}}}}}}_{e}} = \nabla {{p}_{{e \bot }}} + {\mathbf{b}}\left( {{\mathbf{b}}{\kern 1pt} ,\nabla {\kern 1pt} \left( {{{p}_{{e\parallel }}} - {{p}_{e}}_{ \bot }} \right)} \right) + \\ + \,\,\left( {{{p}_{{e\parallel }}} - {{p}_{e}}_{ \bot }} \right)\left( {\left( {{\mathbf{b}},\nabla } \right){\mathbf{b}} + {\mathbf{b}}{\kern 1pt} {\text{div}}{\mathbf{b}}} \right). \\ \end{gathered} $

Последние слагаемые в каждой скобке в правой части (16) для электронов и протонов можно оценить в приближении изотропного давления ${{{\mathbf{\hat {{\rm P}}}}}_{\beta }} \approx {{p}_{\beta }}{\hat {I}}\,:$

(29)
$\left. \begin{gathered} \frac{{{\text{div}}{{{{\mathbf{\hat {{\rm P}}}}}}_{e}}}}{{e{{n}_{e}}}} \sim \frac{{\nabla {{p}_{e}}}}{{e{{n}_{e}}}} = \frac{{\nabla \left( {e{{n}_{e}}{{T}_{e}}} \right)}}{{e{{n}_{e}}}} \sim \frac{{{{T}_{e}}}}{L} \sim {{10}^{{ - 6}}}\,\,{{\text{В }} \mathord{\left/ {\vphantom {{\text{В }} {\text{м }}}} \right. \kern-0em} {\text{м }}} \hfill \\ \frac{{{\text{div}}{{{{\mathbf{\hat {{\rm P}}}}}}_{p}}}}{{e{{n}_{p}}}} \sim \frac{{\nabla {{p}_{p}}}}{{e{{n}_{p}}}} = \frac{{\nabla \left( {e{{n}_{p}}{{T}_{p}}} \right)}}{{e{{n}_{p}}}} \sim \frac{{{{T}_{p}}}}{L} \sim {{10}^{{ - 6}}}\,\,{{\text{В }} \mathord{\left/ {\vphantom {{\text{В }} {{\text{м }}{\text{.}}}}} \right. \kern-0em} {{\text{м }}{\text{.}}}} \hfill \\ \end{gathered} \right\}$

Отсюда следует, что

(30)
$\begin{gathered} {{\left( {\frac{{\partial {\mathbf{j}}}}{{\partial t}}} \right)}_{P}} = \frac{e}{{{{m}_{e}}}}\left( {{\text{div}}{{{{\mathbf{\hat {{\rm P}}}}}}_{e}} - \frac{{{{m}_{e}}}}{{{{m}_{p}}}}{\text{div}}{{{{\mathbf{\hat {{\rm P}}}}}}_{p}}} \right) \approx \\ \approx \frac{e}{{{{m}_{e}}}}{\text{div}}{{{{\mathbf{\hat {{\rm P}}}}}}_{e}} \sim {{10}^{{ - 6}}}\,\,{{\text{А }} \mathord{\left/ {\vphantom {{\text{А }} {({{{\text{м }}}^{2}}\,\,{\text{с }})}}} \right. \kern-0em} {({{{\text{м }}}^{2}}\,\,{\text{с }})}} \gg \frac{{\partial {\mathbf{j}}}}{{\partial t}}. \\ \end{gathered} $

Из (23), (26), (27), (29) и (30) вытекает, что для магнитного острова как в стационарном случае, так и в случае колебаний верны оценки

(31)
$\begin{gathered} \frac{{\partial {\mathbf{j}}}}{{\partial {\kern 1pt} t}} \approx {{\left( {\frac{{\partial {\mathbf{j}}}}{{\partial t}}} \right)}_{E}} + {{\left( {\frac{{\partial {\mathbf{j}}}}{{\partial {\kern 1pt} t}}} \right)}_{B}} + {{\left( {\frac{{\partial {\mathbf{j}}}}{{\partial t}}} \right)}_{P}} \approx \\ \approx \frac{{{{e}^{2}}{{n}_{e}}}}{{{{m}_{e}}}}\left( {{\kern 1pt} {\mathbf{E}} + \left[ {{{{\mathbf{u}}}_{e}} \times {\mathbf{B}}} \right]{\kern 1pt} + \frac{{{\text{div}}{{{{\mathbf{\hat {{\rm P}}}}}}_{e}}}}{{e{{n}_{e}}}}} \right) \approx 0, \\ \end{gathered} $
которые означают, что имеет место приближенное равновесие (гидростатика) электронов. В результате получаем значение электрического поля:

(32)
${\mathbf{E}} \approx - \left[ {{{{\mathbf{u}}}_{e}} \times {\mathbf{B}}} \right]{\kern 1pt} \, - \frac{{{\text{div}}{{{{\mathbf{\hat {{\rm P}}}}}}_{e}}}}{{e{{n}_{e}}}}.$

Отсюда, используя также (28), можно оценить продольную, экранирующую, компоненту электрического поля ${{E}_{{{\kern 1pt} \parallel }}} = \left( {{\mathbf{E}},{\mathbf{b}}} \right),$ обусловленную разделением зарядов:

(33)
$\begin{gathered} {{E}_{{{\kern 1pt} \parallel }}} \approx \frac{1}{{e{{n}_{e}}}}\left( {{\text{div}}{{{{\mathbf{\hat {{\rm P}}}}}}_{e}}{\kern 1pt} ,{\mathbf{b}}} \right) = \\ = \frac{{\left( {{\mathbf{b}},\nabla {{p}_{{e\parallel }}}} \right) + \left( {{{p}_{{e\parallel }}} - {{p}_{e}}_{ \bot }} \right){\text{div}}{\mathbf{b}}}}{{e{{n}_{e}}}} \sim {{10}^{{ - 6}}}\,\,{{\text{В }} \mathord{\left/ {\vphantom {{\text{В }} {{\text{м }}{\text{.}}}}} \right. \kern-0em} {{\text{м }}{\text{.}}}} \\ \end{gathered} $

С другой стороны, известна оценка максимального отклонения от электронейтральности в неоднородной плазме (см., например, Лифшиц и Питаевский, стр. 186), которая в рассматриваемом случае с учетом обозначений в (5) имеет вид:

(34)
${\kern 1pt} \frac{1}{{{{n}_{e}}}}\left( {{{n}_{e}} - {{n}_{p}}} \right) \sim {{\left( {\frac{{{{{\lambda }}_{{{\text{De}}}}}}}{L}} \right)}^{2}},$
где L – размер неоднородности, ${{{\lambda }}_{{{\text{De}}}}} = {{{{V}_{{{\text{Te}}}}}} \mathord{\left/ {\vphantom {{{{V}_{{{\text{Te}}}}}} {{{{\omega }}_{{pe}}}}}} \right. \kern-0em} {{{{\omega }}_{{pe}}}}}$ – дебаевское расстояние, ${{{\omega }}_{{pe}}} = e\sqrt {{{{{n}_{e}}} \mathord{\left/ {\vphantom {{{{n}_{e}}} {\left( {{{{\varepsilon }}_{0}}{{m}_{e}}} \right)}}} \right. \kern-0em} {\left( {{{{\varepsilon }}_{0}}{{m}_{e}}} \right)}}} $ – плазменная частота, ${{V}_{{{\text{Te}}}}} = \sqrt {{{e{{T}_{e}}} \mathord{\left/ {\vphantom {{e{{T}_{e}}} {{{m}_{e}}}}} \right. \kern-0em} {{{m}_{e}}}}} $ – тепловая скорость электронов. В солнечном ветре указанные выше перед (20) и (22) характерные значения параметров дают значение дебаевского расстояния ${{{\lambda }}_{{{\text{De}}}}} \sim 10$ м. Для потенциальной части электрического поля ${{{\mathbf{E}}}_{p}}\left( {{\mathbf{x}},t} \right) = - \nabla {\varphi }\left( {{\mathbf{x}},t} \right)$ из уравнения Пуассона (9) вытекает уравнение

${\text{div}}{{{\mathbf{E}}}_{p}}\left( {{\mathbf{x}},t} \right) = {{{{\rho \left( {{\mathbf{x}},t} \right)} \mathord{\left/ {\vphantom {{\rho \left( {{\mathbf{x}},t} \right)} {\varepsilon }}} \right. \kern-0em} {\varepsilon }}}_{0}}.$

Из этого уравнения и (34) получаем следующую оценку:

(35)
$\begin{gathered} {\kern 1pt} {{E}_{p}} \sim L\frac{\rho }{{{{\varepsilon }_{0}}}} = L\frac{{e{{n}_{e}}}}{{{{\varepsilon }_{0}}}}\frac{1}{{{{n}_{e}}}}\left( {{{n}_{e}} - {{n}_{p}}} \right) \sim L\frac{{e{{n}_{e}}}}{{{{\varepsilon }_{0}}}}{{\left( {\frac{{{{\lambda }_{{{\text{De}}}}}}}{L}} \right)}^{2}} = \\ = \frac{{e{{n}_{e}}}}{{{{\varepsilon }_{0}}L}}{{\left( {{{\lambda }_{{{\text{De}}}}}} \right)}^{2}} = \frac{{{{T}_{e}}}}{L} \sim {{10}^{{ - 6}}}\,\,{{\text{В }} \mathord{\left/ {\vphantom {{\text{В }} {\text{м }}}} \right. \kern-0em} {\text{м }}}, \\ \end{gathered} $
которая, как мы видим, совпадает с оценкой электрического поля за счет дивергенции тензора давления электронов. Также отметим, что из (35) и (32) следует, что в стационарном случае ортогональная часть гидродинамической скорости электронов ${{{\mathbf{u}}}_{{e \bot }}}$ должна быть достаточно малой.

Таким образом, из оценок (33) и (35) вытекает, что продольное электрическое поле экранировки за счет разделения заряда в магнитном острове по порядку величины составляет ${{E}_{{{\kern 1pt} \parallel }}} \sim {{10}^{{ - 6}}}\,\,{{\text{В }} \mathord{\left/ {\vphantom {{\text{В }} {\text{м }}}} \right. \kern-0em} {\text{м }}}.$ В случае колебаний острова, согласно уравнению Фарадея (8), появляется индукционное поле

(36)
${\kern 1pt} {{E}_{v}} \sim L\frac{{{\Delta }B}}{{\Theta }} \sim {{10}^{7}} \times \frac{{5 \times {{{10}}^{{ - 9}}}}}{{100}} \sim 5 \times {{10}^{{ - 4}}}\,\,{{\text{В }} \mathord{\left/ {\vphantom {{\text{В }} {\text{м }}}} \right. \kern-0em} {\text{м }}},$
которое, как следует из оценки (36), на два порядка больше поля экранировки. При этом ортогональная к B часть индукционного электрического поля в приближенном уравнении (32) компенсируется слагаемым – $\left[ {{{{\mathbf{u}}}_{e}} \times {\mathbf{B}}} \right].$ Итак, мы видим, что индукционное электрическое поле, возникающее в колеблющемся магнитном острове, оказывается достаточно большим, чтобы поле разделения зарядов (экранировки) не могло его скомпенсировать. Как результат, индуктивное поле может эффективно ускорять заряженные частицы, попадающие внутрь магнитного острова.

МОДЕЛЬ ПОЛЕЙ В МАГНИТНОМ ОСТРОВЕ

Самосогласованное моделирование распределения полей и параметров плазмы в МО, находящемся в солнечном ветре на расстоянии около 1 а. е. от Солнца в районе орбиты Земли, представляет собой большую и сложную задачу, которая лежит далеко за рамками данной работы. Построить разумную модель потенциального электрического поля весьма сложно. Как показывают приведенные выше оценки (35), (36) это поле примерно на два порядка меньше индукционного электрического поля, возникающего из-за колебаний магнитного острова. Поэтому в приближенной аналитической модели потенциальное электрическое поле ${{{\mathbf{E}}}_{p}}\left( {{\mathbf{x}},t} \right)$ мы не учитываем и считаем, что в острове электрическое поле является чисто индукционным. Таким образом, поля определяются векторным потенциалом ${\mathbf{A}}({\mathbf{x}},t)$ по формулам

(37)
$\begin{gathered} {\mathbf{B}}{\kern 1pt} \left( {{\mathbf{x}},t} \right) = {{B}_{{x0}}}{{{\mathbf{e}}}_{x}} + {{B}_{{y0}}}{{{\mathbf{e}}}_{y}} + {\text{rot}}{\mathbf{A}}\left( {{\mathbf{x}},t} \right), \\ {\mathbf{E}}\left( {{\mathbf{x}},t} \right) = - \frac{{\partial {\mathbf{A}}\left( {{\mathbf{x}},t} \right)}}{{\partial t}}. \\ \end{gathered} $

В этом случае формально выполнено условие электронейтральности

(38)
${\text{div}}{\mathbf{E}}{\kern 1pt} \left( {{\mathbf{x}},t} \right) = 0,$
которое равносильно условию

(39)
${\text{div}}{\mathbf{A}}\left( {{\mathbf{x}},t} \right) = 0.$

Для автоматического выполнения условий (38) и (39) удобно задавать векторный потенциал в форме

(40)
$\begin{gathered} {\mathbf{A}}{\kern 1pt} \left( {{\mathbf{x}},t} \right) = L{{B}_{m}}{\text{rot}}\left( {{{\psi }_{1}}\left( {{\mathbf{x}},t} \right)\nabla {{\psi }_{2}}\left( {\mathbf{x}} \right)} \right) = \\ = L{\kern 1pt} {{B}_{m}}\left[ {\nabla {{\psi }_{1}}\left( {{\mathbf{x}},t} \right) \times \nabla {{\psi }_{2}}\left( {\mathbf{x}} \right)} \right], \\ \end{gathered} $
где L – модельный параметр, имеющий размерность расстояния, и ${{B}_{m}}$ – модельный параметр, имеющий размерность магнитного поля. Формула (40) означает, что векторный потенциал определяется своими эйлеровыми потенциалами: зависящим от времени первым потенциалом ${{\psi }_{1}}\left( {{\mathbf{x}},t} \right)$ и стационарным вторым потенциалом ${{\psi }_{2}}\left( {\mathbf{x}} \right),$ которые должны быть независимы друг от друга, то есть их градиенты должны быть не параллельны. Из формул (40) и (37) вытекают следующие выражения для полей, выраженные через функции ${{\psi }_{1}}\left( {{\mathbf{x}},t} \right)$ и ${{\psi }_{2}}\left( {\mathbf{x}} \right)\,:$

(41)
$\left. \begin{gathered} {\mathbf{B}}{\kern 1pt} \left( {{\mathbf{x}},t} \right) = {{B}_{{x0}}}{{{\mathbf{e}}}_{x}} + {{B}_{{y0}}}{{{\mathbf{e}}}_{y}} + {\text{rot}}{\mathbf{A}}\left( {{\mathbf{x}},t} \right) = {{B}_{{x0}}}{{{\mathbf{e}}}_{x}} + {{B}_{{y0}}}{{{\mathbf{e}}}_{y}} + L{{B}_{m}}{\text{rotrot}}\left( {{\kern 1pt} {{\psi }_{1}}\left( {{\mathbf{x}},t} \right)\nabla {{\psi }_{2}}\left( {\mathbf{x}} \right)} \right), \hfill \\ {\mathbf{E}}{\kern 1pt} \left( {{\mathbf{x}},t} \right) = - \frac{{\partial {\mathbf{A}}\left( {{\mathbf{x}},t} \right)}}{{\partial t}} = L{{B}_{m}}\left[ {\nabla {{\psi }_{2}}\left( {\mathbf{x}} \right) \times \nabla \frac{{\partial {{\psi }_{1}}\left( {{\mathbf{x}},t} \right)}}{{\partial t}}} \right]. \hfill \\ \end{gathered} \right\}$

Будем считать, что магнитный остров начал колебания в момент времени ${{T}_{0}}$ с периодом Θ. Удобно представить время в виде $t = {{T}_{0}} + {{n}_{t}}\Theta + \tau ,$ где τ – текущая фаза колебаний острова, а ${{n}_{t}} = \left[ {{{\left( {t - {{T}_{0}}} \right)} \mathord{\left/ {\vphantom {{\left( {t - {{T}_{0}}} \right)} \Theta }} \right. \kern-0em} \Theta }} \right] \in \mathbb{Z}$ – целое число. В качестве первого потенциала ${{\psi }_{1}}\left( {{\mathbf{x}},t} \right)$ удобно выбрать параметризацию эллипсоида, у которого полуоси колеблются вокруг своих средних значений ${{L}_{x}},{{L}_{y}},{{L}_{z}}$ с периодом $\Theta = {{2\pi } \mathord{\left/ {\vphantom {{2\pi } \omega }} \right. \kern-0em} \omega }\,:$

(42)
$\left. \begin{gathered} {{\psi }_{1}}\left( {{\mathbf{x}},t} \right) = \hfill \\ = {\text{exp}}\left( { - {{{\left| {{{S}_{x}}\left( {x,\tau } \right)} \right|}}^{2}} - {{{\left| {{{S}_{y}}\left( {y,\tau } \right)} \right|}}^{2}} - {{{\left| {{{S}_{z}}\left( {z,\tau } \right)} \right|}}^{2}}} \right), \hfill \\ {\text{г д е }} \hfill \\ {{S}_{x}}\left( {x,\tau } \right) = \frac{{x - {{x}_{0}}}}{{{{L}_{x}} - {{\lambda }_{x}}\sin \omega \tau }},\,\,\,\, \hfill \\ {{S}_{y}}\left( {y,\tau } \right) = \frac{{y - {{y}_{0}}}}{{{{L}_{y}} - {{\lambda }_{y}}\sin \omega \tau }},\,\, \hfill \\ \,\,{{S}_{z}}\left( {z,\tau } \right) = \frac{{z - {{z}_{0}}}}{{{{L}_{z}} - {{\lambda }_{z}}\sin \omega \tau }}. \hfill \\ \end{gathered} \right\}$

Здесь ${{{\mathbf{x}}}_{0}} = \left( {{{x}_{0}},{{y}_{0}},{{z}_{0}}} \right)$ – координаты центра магнитного острова, ${{\lambda }_{x}},{{\lambda }_{y}},{{\lambda }_{z}}$ – амплитуды колебаний полуосей, которые должны удовлетворять условию $0 < {{\lambda }_{k}} < {{{{L}_{k}}} \mathord{\left/ {\vphantom {{{{L}_{k}}} 2}} \right. \kern-0em} 2},$ $k = x,y,z.$

В качестве 2-го потенциала можно выбрать функцию

(43)
${{\psi }_{2}}\left( {\mathbf{x}} \right) = \frac{1}{2}{{\left( {{\mathbf{x}} - {{{\mathbf{x}}}_{0}}} \right)}^{2}},\,\,\,{\text{т о е с т ь }}\,\,\,\,\nabla {{\psi }_{2}}\left( {\mathbf{x}} \right) = {\mathbf{x}} - {{{\mathbf{x}}}_{0}}.$

В представленных расчетах мы рассматривали три варианта пространственного размера $L = {{L}_{x}}\,:$ $L = {{10}^{8}};\,\,1.5 \times {{10}^{8}};\,\,2 \times {{10}^{8}}\,\,{\text{м ,}}$ а остальные пространственные модельные параметры определяли по нему следующим образом: ${{L}_{y}} = {{L}_{z}} = {L \mathord{\left/ {\vphantom {L 2}} \right. \kern-0em} 2},$ ${{\lambda }_{x}} = {{\lambda }_{z}} = 0,$ ${{\lambda }_{y}} = {{{{L}_{y}}} \mathord{\left/ {\vphantom {{{{L}_{y}}} 2}} \right. \kern-0em} 2} = {L \mathord{\left/ {\vphantom {L 4}} \right. \kern-0em} 4}.$ То есть рассматривался остров в виде эллипсоида, вытянутого в направлении оси Х, у которого колеблется только полуось в направлении оси Y. Компоненты магнитного поля солнечного ветра имели значения ${{B}_{{x0}}} = {{B}_{{y0}}} = 5\,\,{\text{н Т л }}{\text{.}}$ Рассматривалось три варианта для амплитуды колебаний магнитного поля острова ${{B}_{m}} = 0.5{{B}_{{x0}}};$ $0.75{{B}_{{x0}}};\,\,{{B}_{{x0}}},$ а также два варианта периода колебаний $\Theta = {{2\pi } \mathord{\left/ {\vphantom {{2\pi } \omega }} \right. \kern-0em} \omega } = 300;\,\,600\,\,{\text{с }}$.

Начало координат системы отсчета $K{\text{'}}$ удобно поместить в центр острова, тогда ${{x}_{0}} = 0.$ Также удобно ввести безразмерные переменные $\tilde {x} = {x \mathord{\left/ {\vphantom {x {L\,,}}} \right. \kern-0em} {L\,,}}\;\;\tilde {y} = {y \mathord{\left/ {\vphantom {y {L\,,}}} \right. \kern-0em} {L\,,}}\;\;\tilde {z} = {z \mathord{\left/ {\vphantom {z L}} \right. \kern-0em} L}$ и параметр $\tilde {\lambda } = {{{{\lambda }_{y}}} \mathord{\left/ {\vphantom {{{{\lambda }_{y}}} L}} \right. \kern-0em} L} = {1 \mathord{\left/ {\vphantom {1 4}} \right. \kern-0em} 4}$.

В этом случае из формул (3), (41)–(43) получаются следующие выражения для компонент полей:

(44)
$\begin{gathered} {{\psi }_{0}}\left( t \right) = {{\left( {1 - 2\tilde {\lambda }\sin \omega \tau } \right)}^{{ - 2}}}, \\ {{\psi }_{1}}\left( {{\mathbf{x}},t} \right) = {\text{exp}}\left( { - {{{\tilde {x}}}^{2}} - 4{{{\tilde {z}}}^{2}} - 4{{{\tilde {y}}}^{2}}{{\psi }_{0}}\left( t \right)} \right), \\ \end{gathered} $
(45)
$\left. \begin{gathered} {{B}_{x}}\left( {{\mathbf{x}},t} \right) = {{B}_{{x0}}} + {\text{4}}{{B}_{m}}{{\psi }_{1}}\left( {{\mathbf{x}},t} \right) \times \hfill \\ \times \,\,\tilde {x}\left( {1 - 12{{{\tilde {z}}}^{2}} + 2{{\psi }_{0}}\left( t \right)\left( {1 + 2{{{\tilde {y}}}^{2}}\left( {1 - 4{{\psi }_{0}}\left( t \right)} \right)} \right)} \right), \hfill \\ {{B}_{y}}\left( {{\mathbf{x}},t} \right) = {\kern 1pt} {{B}_{{y0}}} + {\text{2}}{{B}_{m}}{{\psi }_{1}}\left( {{\mathbf{x}},t} \right) \times \hfill \\ \times \,\,\tilde {y}\left( {5 - 2{{{\tilde {x}}}^{2}} - 32{{{\tilde {z}}}^{2}} + 8{{\psi }_{0}}\left( t \right)\left( {4{{{\tilde {z}}}^{2}} + {{{\tilde {x}}}^{2}} - 1} \right)} \right), \hfill \\ {{B}_{z}}\left( {{\mathbf{x}},t} \right) = 2{{B}_{m}}{{\psi }_{1}}\left( {{\mathbf{x}},t} \right)\tilde {z}\left( {6{{{\tilde {x}}}^{2}} - 7 + 4{{\psi }_{0}}\left( t \right)} \right. \times \hfill \\ \times \,\,\left. {\left( {1 - 32\tilde {\lambda }{{{\tilde {y}}}^{2}}{{\psi }_{0}}\left( t \right)\left( {1 - \tilde {\lambda }\sin \omega \tau } \right)\sin \omega \tau } \right)} \right) \hfill \\ \end{gathered} \right\}$
(46)
${{E}_{0}}\left( {{\mathbf{x}},t} \right) = \frac{{32\omega {\kern 1pt} L{{B}_{m}}{\tilde {\lambda }}\tilde {y}\cos {\omega }\tau }}{{{{{\left( {1 - 2{\tilde {\lambda }}\sin \omega {\kern 1pt} \tau } \right)}}^{3}}}}{{\psi }_{1}}\left( {{\mathbf{x}},t} \right),$
(47)
$\left. \begin{gathered} {{E}_{x}}\left( {{\mathbf{x}},t} \right) = {{E}_{0}}\left( {{\mathbf{x}},t} \right) \times \hfill \\ \times \,\,\tilde {z}\left( {1 - 32\tilde {\lambda }{{{\tilde {y}}}^{2}}{{\psi }_{0}}\left( t \right)\left( {1 - \tilde {\lambda }\sin \omega \tau } \right)\sin \omega \tau } \right), \hfill \\ {{E}_{y}}\left( {{\mathbf{x}},t} \right) = - 3{{E}_{0}}\left( {{\mathbf{x}},t} \right)\tilde {x}\tilde {z}\tilde {y}, \hfill \\ {{E}_{z}}\left( {{\mathbf{x}},t} \right) = {{E}_{0}}\left( {{\mathbf{x}},t} \right)\tilde {x}\left( {{{{\tilde {y}}}^{2}}\left( {4{{\psi }_{0}}\left( t \right) - 1} \right) - 1} \right). \hfill \\ \end{gathered} \right\}$

Отметим, что из формул (45)–(47) вытекает наличие в предложенной модели полей в течение большей части периода колебаний в основной области магнитного острова – эллипсоиде

$\left\{ {{\mathbf{x}}\,:{{\left( {{{x}^{2}} + 4{{y}^{2}} + 4{{z}^{2}}} \right)} \mathord{\left/ {\vphantom {{\left( {{{x}^{2}} + 4{{y}^{2}} + 4{{z}^{2}}} \right)} {{{L}^{2}} \leqslant 1}}} \right. \kern-0em} {{{L}^{2}} \leqslant 1}}} \right\}$
существенного продольного электрического поля.

МЕТОДИКА МОДЕЛИРОВАНИЯ

Для расчета траектории частицы с зарядом e и массой покоя ${{m}_{0}}$ необходимо численно решать задачу Коши для системы уравнений движения заряда в электромагнитном поле (систему Ньютона–Лоренца). В системе единиц СИ в релятивистском случае она имеет вид:

(48)
$\begin{gathered} \frac{{d{\mathbf{x}}\left( t \right)}}{{dt}} = {\mathbf{v}}\left( t \right),\,\,\,\,{\mathbf{p}}\left( {\mathbf{v}} \right) = \frac{{{{m}_{0}}{\mathbf{v}}}}{{\sqrt {{\kern 1pt} 1 - {{{{{\left| {\mathbf{v}} \right|}}^{2}}} \mathord{\left/ {\vphantom {{{{{\left| {\mathbf{v}} \right|}}^{2}}} {{{c}^{2}}{\kern 1pt} }}} \right. \kern-0em} {{{c}^{2}}{\kern 1pt} }}} }}, \\ \frac{{d{\mathbf{p}}\left( {{\mathbf{v}}\left( t \right)} \right)}}{{dt}} = e\left( {{\mathbf{E}}\left( {{\mathbf{x}}\left( t \right),t} \right) + \left[ {{\mathbf{v}}\left( t \right) \times {\mathbf{B}}\left( {{\mathbf{x}}\left( t \right),t} \right)} \right]} \right), \\ \end{gathered} $
где $c$ – скорость света в вакууме, p – импульс. В классическом пределе ${{\left| {\mathbf{v}} \right|} \mathord{\left/ {\vphantom {{\left| {\mathbf{v}} \right|} c}} \right. \kern-0em} c} \ll {\text{1}}$ эта система принимает вид

(49)
$\begin{gathered} \frac{{d{\mathbf{x}}\left( t \right)}}{{dt}} = {\mathbf{v}}\left( t \right), \\ \frac{{d{\mathbf{v}}\left( t \right)}}{{d{\kern 1pt} t}} = \frac{e}{{{{m}_{0}}}}\left( {{\mathbf{E}}\left( {{\mathbf{x}}\left( t \right),t} \right) + \left[ {{\mathbf{v}}\left( t \right){\mathbf{ \times B}}\left( {{\mathbf{x}}\left( t \right),t} \right)} \right]{\kern 1pt} } \right). \\ \end{gathered} $

Начальные условия следующие:

(50)
${\mathbf{x}}\left( {{{t}^{0}}} \right) = {{{\mathbf{x}}}^{0}},\,\,\,\,{\mathbf{v}}\left( {{{t}^{0}}} \right) = {{{\mathbf{v}}}^{0}}.$

Кинетическая энергия в электрон-вольтах $W\left( v \right)$ связана с модулем скорости $v = \left| {\mathbf{v}} \right|$ формулами:

(51)
$\begin{gathered} W(\left| {\mathbf{v}} \right|) = \frac{{{{m}_{0}}{{c}^{2}}}}{e}\left( {\frac{1}{{\sqrt {{\kern 1pt} 1 - {{{{{\left| {\mathbf{v}} \right|}}^{2}}} \mathord{\left/ {\vphantom {{{{{\left| {\mathbf{v}} \right|}}^{2}}} {{{c}^{2}}{\kern 1pt} }}} \right. \kern-0em} {{{c}^{2}}{\kern 1pt} }}} }} - 1} \right), \\ \frac{{\left| {\mathbf{v}} \right|}}{с } = {{\sqrt {{\kern 1pt} \frac{{eW}}{{{{m}_{0}}{{c}^{2}}}}\left( {\frac{{eW}}{{{{m}_{0}}{{c}^{2}}}} + 2} \right)} } \mathord{\left/ {\vphantom {{\sqrt {{\kern 1pt} \frac{{eW}}{{{{m}_{0}}{{c}^{2}}}}\left( {\frac{{eW}}{{{{m}_{0}}{{c}^{2}}}} + 2} \right)} } {\left( {\frac{{eW}}{{{{m}_{0}}{{c}^{2}}}} + 1} \right)}}} \right. \kern-0em} {\left( {\frac{{eW}}{{{{m}_{0}}{{c}^{2}}}} + 1} \right)}}. \\ \end{gathered} $

В классическом пределе она переходит в кинетическую энергию для классического приближения: $W\left( v \right) \to \frac{1}{{2e}}{{m}_{0}}{{v}^{2}}$ при $\frac{v}{с } \to 0.$ Отметим, что из уравнений (48) и (51) вытекает уравнение

(52)
$\frac{{dW\left( {\left| {{\mathbf{v}}\left( t \right)} \right|} \right)}}{{dt}} = \left( {{\mathbf{E}}\left( {{\mathbf{x}}\left( t \right),t} \right),{\mathbf{v}}\left( t \right)} \right).$

В ходе расчета траектории для экономии вычислительных ресурсов при значении параметра $\beta = {{\left| {\mathbf{v}} \right|} \mathord{\left/ {\vphantom {{\left| {\mathbf{v}} \right|} c}} \right. \kern-0em} c} < {{\beta }_{0}} = 0.02$ (которому соответствует кинетическая энергия ${{W}_{0}} \approx 188$ кэВ) выполняется временной шаг в рамках классической системы уравнений (49), а при $\beta \geqslant {{\beta }_{0}} = 0.02$ выполняется временной шаг в рамках релятивистской системы уравнений (48). Детальное описание численной методики расчета траекторий протонов в рамках классической и релятивистской систем уравнений представлено в Приложении.

Решение задачи Коши (48), (50) для последующего изложения удобно обозначить как

(53)
${\mathbf{x}}\left( t \right) = {\mathbf{X}}\left( {t,{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right),\,\,\,{\mathbf{v}}\left( t \right) = {\mathbf{V}}\left( {t,{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right),$
где ${{t}^{0}},$ ${{{\mathbf{x}}}^{0}}$ и ${{{\mathbf{v}}}^{0}}$ – соответственно начальные время, точка и скорость (вообще далее верхний индекс нуль обозначает начальные данные), чтобы анализировать зависимость результатов моделирования от начальных данных частиц. То есть эти векторные функции удовлетворяют следующим уравнениям и начальным условиям:

$\begin{gathered} \frac{{\partial {\mathbf{X}}\left( {t,{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right)}}{{\partial {\kern 1pt} t}} = {\mathbf{V}}\left( {t,{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right),\,\,\,\,\frac{{\partial {\mathbf{V}}\left( {t,{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right)}}{{\partial t}} = \\ = e{\text{'}}\left( {{\mathbf{E}}\left( {{\mathbf{X}}\left( {t,{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right),t} \right)} \right. + \\ + \,\,\left. {\left[ {{\mathbf{V}}\left( {t,{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right) \times {\mathbf{B}}\left( {{\mathbf{X}}\left( {t,{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right),t} \right)} \right]{\kern 1pt} } \right), \\ {\mathbf{X}}\left( {{{t}^{0}},{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right) = {{{\mathbf{x}}}^{0}},\,\,\,\,{\mathbf{V}}\left( {{{t}^{0}},{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right) = {{{\mathbf{v}}}^{0}}. \\ \end{gathered} $

Отметим, что, поскольку дополнительное к полю в солнечном ветре электромагнитное поле острова на его краях плавно нарастает от нуля до максимального значения на краях, то падающие из складок ГТС потоки предускоренных протонов, в принципе, могут попасть в любую часть острова, причем в различные моменты времени, соответствующие различным фазам колебаний острова. Поэтому для моделирования возможного ускорения протонов в работе используется следующий подход. Было выбрано достаточно широкое множество начальных точек $\left\{ {{{{\mathbf{x}}}^{0}}\left[ {{{k}_{x}}} \right]} \right\}$ внутри магнитного острова – эллипсоида:

(54)
${{\Omega }_{0}} = \left\{ {{\mathbf{x}}\,:{{\left( {{{x}^{2}} + 4{{y}^{2}} + 4{{z}^{2}}} \right)} \mathord{\left/ {\vphantom {{\left( {{{x}^{2}} + 4{{y}^{2}} + 4{{z}^{2}}} \right)} {{{L}^{2}} \leqslant 1}}} \right. \kern-0em} {{{L}^{2}} \leqslant 1}}} \right\},$
а также набор начальных фаз колебаний острова ${{{\tau }}^{0}}\left[ {{{k}_{t}}} \right]$ от 0 до Θ с шагом $\frac{{\Theta }}{{16}},$ которому соответствует набор начальных моментов времени $\left\{ {{{t}^{0}}\left[ {{{k}_{t}}} \right] = {{T}_{0}} + {{\tau }^{0}}\left[ {{{k}_{t}}} \right]} \right\}.$ Далее для простоты будем отсчитывать время от начала колебаний острова, то есть положим ${{T}_{0}} = 0$ и отождествим начальное время и начальную фазу: ${{t}^{0}}\left[ {{{k}_{t}}} \right] = {{{\tau }}^{0}}\left[ {{{k}_{t}}} \right].$

Набор вариантов начальной скорости ${{{\mathbf{v}}}^{0}}$ строится следующим образом. Был выбран набор начальных кинетических энергий ${{W}^{{0{\kern 1pt} }}}\left[ {{{k}_{w}}} \right]$ в диапазоне от 10 эВ до 100 кэВ, причем шаг в диапазоне от 10 эВ до 0.5 кэВ составляет 10 эВ, а в диапазоне от 0.5 до 100 кэВ составляет 0.5 кэВ. Каждой начальной энергии в эВ соответствует модуль начальной скорости, определяемый по 2-й формуле в (51):

(55)
${{v}^{0}}\left( {{{W}^{0}}} \right) = c{{\sqrt {{\kern 1pt} \frac{{e{{W}^{0}}}}{{{{m}_{0}}{{c}^{2}}}}\left( {\frac{{e{{W}^{0}}}}{{{{m}_{0}}{{c}^{2}}}} + 2} \right)} } \mathord{\left/ {\vphantom {{\sqrt {{\kern 1pt} \frac{{e{{W}^{0}}}}{{{{m}_{0}}{{c}^{2}}}}\left( {\frac{{e{{W}^{0}}}}{{{{m}_{0}}{{c}^{2}}}} + 2} \right)} } {\left( {\frac{{e{{W}^{0}}}}{{{{m}_{0}}{{c}^{2}}}} + 1} \right)}}} \right. \kern-0em} {\left( {\frac{{e{{W}^{0}}}}{{{{m}_{0}}{{c}^{2}}}} + 1} \right)}}.$

Набор из ${{N}_{{\left( {\alpha \beta } \right)}}} = 2{{N}_{\alpha }}\left( {{{N}_{\beta }} + 1} \right) = 360 \times 181$ начальных скоростей с данной начальной энергией ${{W}^{{0{\kern 1pt} }}}$ и с шагом в 10 по углам направления скорости

(56)
$\begin{gathered} \left\{ {{{{\mathbf{v}}}^{0}}\left[ {{{W}^{0}},{{k}_{\alpha }},{{k}_{\beta }}} \right]} \right.,\,\,\,{{k}_{\alpha }} = - {{N}_{\alpha }} + 1, \ldots ,{{N}_{\alpha }}, \\ \left. {^{{}}{{k}_{\beta }} = 0, \ldots ,{{N}_{\beta }},\,\;\;{{N}_{\alpha }} = {{N}_{\beta }} = 180} \right\} \\ \end{gathered} $
определяется следующими формулами для компонент начальных скоростей:

(57)
$\left. \begin{gathered} v_{x}^{0}\left[ {{{W}^{0}},{{k}_{\alpha }},{{k}_{\beta }}} \right] = {{v}^{0}}\left( {{{W}^{0}}} \right){\text{sin}}\left( {{{\pi {{k}_{\beta }}} \mathord{\left/ {\vphantom {{\pi {{k}_{\beta }}} {{{N}_{\beta }}}}} \right. \kern-0em} {{{N}_{\beta }}}}} \right){\text{cos}}\left( {{{\pi {{k}_{\alpha }}} \mathord{\left/ {\vphantom {{\pi {{k}_{\alpha }}} {{{N}_{\alpha }}}}} \right. \kern-0em} {{{N}_{\alpha }}}}} \right), \\ v_{y}^{0}\left[ {{{W}^{0}},{{k}_{\alpha }},{{k}_{\beta }}} \right] = {{v}^{0}}\left( {{{W}^{0}}} \right){\text{sin}}\left( {{{\pi {{k}_{\beta }}} \mathord{\left/ {\vphantom {{\pi {{k}_{\beta }}} {{{N}_{\beta }}}}} \right. \kern-0em} {{{N}_{\beta }}}}} \right){\text{sin}}\left( {{{\pi {{k}_{\alpha }}} \mathord{\left/ {\vphantom {{\pi {{k}_{\alpha }}} {{{N}_{\alpha }}}}} \right. \kern-0em} {{{N}_{\alpha }}}}} \right), \\ v_{z}^{0}\left[ {{{W}^{0}},{{k}_{\alpha }},{{k}_{\beta }}} \right] = {{v}^{0}}\left( {{{W}^{0}}} \right){\text{cos}}\left( {{{\pi {{k}_{\beta }}} \mathord{\left/ {\vphantom {{\pi {{k}_{\beta }}} {{{N}_{\beta }}}}} \right. \kern-0em} {{{N}_{\beta }}}}} \right). \\ \end{gathered} \right\}$

Таким образом, формируется множество начальных условий $\left\{ {{{t}^{{0}}}\left[ {{{k}_{t}}} \right],{{{\mathbf{x}}}^{0}}\left[ {{{k}_{x}}} \right],{{{\mathbf{v}}}^{0}}\left[ {{{W}^{0}}\left[ {{{k}_{w}}} \right],{{k}_{\alpha }},{{k}_{\beta }}} \right]} \right\}.$ Используется описанная в разделе “Геометрия системы…” движущаяся вместе с солнечным ветром система отсчета $K{\text{',}}$ в которой центр острова неподвижен, а магнитное и электрическое поля определяются формулами (44)(47). Область контроля траекторий (область моделирования) задана в виде прямоугольного параллелепипеда

(58)
$\Pi = \left\{ {{\mathbf{x}}\,:\left| {{x \mathord{\left/ {\vphantom {x L}} \right. \kern-0em} L}} \right| \leqslant 1.25,\,\,\,\left| {{y \mathord{\left/ {\vphantom {y L}} \right. \kern-0em} L}} \right| \leqslant 1,\,\,\,\left| {{z \mathord{\left/ {\vphantom {z L}} \right. \kern-0em} L}} \right| \leqslant 0.75} \right\},$
содержащий эллипсоид ${{\Omega }_{0}},$ заданный согласно (4.7).

Для каждого варианта поля ${{{\mathbf{F}}}^{{(i)}}}\left[ {{{k}_{f}}} \right]\left( {{\mathbf{x}},t} \right)$ модельного острова (где ${{{\mathbf{F}}}^{{(i)}}} = \left\{ {{{{\mathbf{B}}}^{{(i)}}},{{{\mathbf{E}}}^{{(i)}}}} \right\}$), которые описаны в разделе “Модель полей…”, и каждого начального условия $\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right)$ были рассчитаны траектории протонов, т.е. находилось решение задачи Коши (47), (50) вплоть до вылета частицы из области контроля Π. В результате рассчитывается время нахождения траектории в области контроля Π, которое обозначим через ${{T}_{{{\text{end}}}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right),$ координата точки вылета ${{{\mathbf{x}}}_{{{\text{end}}}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right),$ и скорость вылета ${{{\mathbf{v}}}_{{{\text{end}}}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right),$ которая по первой формуле в (4.4) дает энергию вылета (в эВ):

(59)
$\begin{gathered} {{W}_{{{\text{end}}}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right) = \\ = \frac{{{{m}_{0}}{{c}^{2}}}}{e}\left( {\frac{1}{{\sqrt {{\kern 1pt} 1 - |{{{{{\mathbf{v}}}_{{{\text{end}}}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right)} \mathord{\left/ {\vphantom {{{{{\mathbf{v}}}_{{{\text{end}}}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right)} {{{{\left. c \right|}}^{2}}}}} \right. \kern-0em} {{{{\left. c \right|}}^{2}}}}} }} - 1} \right). \\ \end{gathered} $

Для начального момента времени ${{t}^{0}},$ начальной точки ${{{\mathbf{x}}}^{0}}$ и набора (57) начальных скоростей с данной начальной энергией ${{W}^{{0{\kern 1pt} }}}$ будем исследовать как наиболее показательные величины среднюю энергию частиц $\left\langle {{{W}_{{{\text{end}}}}}} \right\rangle \left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{{0{\kern 1pt} }}}} \right)$ и среднее время $\left\langle {{{T}_{{{\text{end}}}}}} \right\rangle \left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right)$ их нахождения в области контроля Π. Эти величины можно определить через параметры каждой траектории следующим образом:

(60)
$\begin{gathered} \left\langle {{{W}_{{{\text{end}}}}}} \right\rangle \left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{{0{\kern 1pt} }}}} \right) = \\ = \frac{1}{{{{N}_{{(\alpha \beta )}}}}}\sum\limits_{{{k}_{\beta }} = 0}^{{{N}_{\beta }}} {\sum\limits_{{{k}_{\alpha }} = - {{N}_{\alpha }} + 1}^{{{N}_{\alpha }}} {{{W}_{{{\text{end}}}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}\left[ {{{W}^{0}},{{k}_{\alpha }},{{k}_{\beta }}} \right]} \right),} } \\ \end{gathered} $
(61)
$\begin{gathered} \left\langle {{{T}_{{{\text{end}}}}}} \right\rangle \left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{{0{\kern 1pt} }}}} \right) = \\ = \frac{1}{{{{N}_{{(\alpha \beta )}}}}}\sum\limits_{{{k}_{\beta }} = 0}^{{{N}_{\beta }}} {\sum\limits_{{{k}_{\alpha }} = - {{N}_{\alpha }} + 1}^{{{N}_{\alpha }}} {{{T}_{{{\text{end}}}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{{0{\kern 1pt} }}}\left[ {{{W}^{0}},{{k}_{\alpha }},{{k}_{\beta }}} \right]} \right)} } . \\ \end{gathered} $

Для каждого начального набора начальных моментов времени и координат $\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}}} \right)$ введенные в (60) и (61) величины являются функциями одной переменной: начальной энергии ${{W}^{0}}\left[ {{{k}_{w}}} \right].$ Графики этих функций наглядно демонстрируют наличие и характер ускорения.

Также представляет интерес анизотропия вылетевших из области контроля $\Pi $ протонов. Для данного набора $\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{{0{\kern 1pt} }}}} \right)$ ее можно описать при помощи функции распределения $f\left[ {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right]{\kern 1pt} \left( {\alpha ,\beta } \right)$ направления скоростей вылета ${{{\mathbf{v}}}_{{{\text{end}}}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}\left[ {{{W}^{0}},k_{\alpha }^{0},k_{\beta }^{0}} \right]} \right)$ по углам вылета в сферической системе координат $\alpha \in \left( { - \pi ;\pi } \right]$ и $\beta \in \left[ {0;\pi } \right].$ Эти углы определяются формулами

(62)
$\left. \begin{gathered} \beta \left[ {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right]\left( {k_{\alpha }^{0},k_{\beta }^{0}} \right) = {\text{arccos}}\left( {{{{{{\mathbf{v}}}_{{{\text{end}},z}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}\left[ {{{W}^{0}},k_{\alpha }^{0},k_{\beta }^{0}} \right]} \right)} \mathord{\left/ {\vphantom {{{{{\mathbf{v}}}_{{{\text{end}},z}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}\left[ {{{W}^{0}},k_{\alpha }^{0},k_{\beta }^{0}} \right]} \right)} {\left| {{{{\mathbf{v}}}_{{{\text{end}}}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}\left[ {{{W}^{0}},k_{\alpha }^{0},k_{\beta }^{0}} \right]{\kern 1pt} } \right)} \right|}}} \right. \kern-0em} {\left| {{{{\mathbf{v}}}_{{{\text{end}}}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}\left[ {{{W}^{0}},k_{\alpha }^{0},k_{\beta }^{0}} \right]{\kern 1pt} } \right)} \right|}}} \right), \\ \alpha \left[ {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right]\left( {k_{\alpha }^{0},k_{\beta }^{0}} \right) = 2{\text{arctg}}\left( {\frac{{{{{\mathbf{v}}}_{{{\text{end}},y}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}\left[ {{{W}^{0}},k_{\alpha }^{0},k_{\beta }^{0}} \right]} \right)}}{{\left| {{{{\mathbf{v}}}_{{{\text{end}}}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{{0{\kern 1pt} }}}\left[ {{{W}^{0}},k_{\alpha }^{0},k_{\beta }^{0}} \right]} \right)} \right| + {{v}_{{{\text{end}},x}}}\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}\left[ {{{W}^{0}},k_{\alpha }^{0},k_{\beta }^{0}} \right]} \right)}}} \right), \\ \end{gathered} \right\}$
которые равносильны соотношениям
$\begin{gathered} {{v}_{{{\text{end}},x}}} = \left| {{\kern 1pt} {{{\mathbf{v}}}_{{{\text{end}}}}}} \right|{\text{sin}}\beta {\text{cos}}\alpha ,\,\,\,{{v}_{{{\text{end}},y}}} = \left| {{\kern 1pt} {{{\mathbf{v}}}_{{{\text{end}}}}}} \right|{\text{sin}}\beta {\text{sin}}\alpha , \\ {{v}_{{{\text{end}},z}}} = \left| {{{{\mathbf{v}}}_{{{\text{end}}}}}} \right|{\text{cos}}\beta . \\ \end{gathered} $
Значения функции $f\left[ {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right]{\kern 1pt} \left( {\alpha ,\beta } \right)$ в узлах сетки
$\begin{gathered} \beta \left[ {{{k}_{\beta }}} \right] = {{\pi {{k}_{{\beta }}}} \mathord{\left/ {\vphantom {{\pi {{k}_{{\beta }}}} {{{N}_{\beta }}}}} \right. \kern-0em} {{{N}_{\beta }}}},\,\,\,\,\alpha \left[ {{{k}_{\alpha }}} \right] = {{\pi {{k}_{\alpha }}} \mathord{\left/ {\vphantom {{\pi {{k}_{\alpha }}} {{{N}_{\alpha }}}}} \right. \kern-0em} {{{N}_{\alpha }}}}, \\ {{k}_{\alpha }} = - {{N}_{\alpha }} + 1, \ldots ,{{N}_{\alpha }},\,\,\,{{k}_{{\beta }}} = 0, \ldots ,{{N}_{\beta }}, \\ {{N}_{\alpha }} = {{N}_{\beta }} = 180 \\ \end{gathered} $
определяются при помощи линейного взвешивания по формулам
(63)
$\begin{gathered} f\left[ {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right]{\kern 1pt} \left( {\alpha \left[ {{{k}_{\alpha }}} \right],\beta \left[ {{{k}_{\beta }}} \right]} \right) = \\ = \sum\limits_{k_{\beta }^{0} = 0}^{{{N}_{\beta }}} {\sum\limits_{k_{\alpha }^{0} = - {{N}_{\alpha }} + 1}^{{{N}_{\alpha }}} {{{W}_{1}}\left( {\frac{{\alpha \left[ {{{k}_{\alpha }}} \right] - \alpha \left[ {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right]\left( {k_{\alpha }^{0},k_{\beta }^{0}} \right)}}{{\Delta \alpha }}} \right)} } \times \\ \times \,\,{{W}_{1}}\left( {\frac{{\beta \left[ {{{k}_{\beta }}} \right] - \beta \left[ {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right]\left( {k_{\alpha }^{0},k_{\beta }^{0}} \right)}}{{\Delta \beta }}} \right), \\ \end{gathered} $
где $\Delta \alpha = {\pi \mathord{\left/ {\vphantom {\pi {{{N}_{\alpha }}}}} \right. \kern-0em} {{{N}_{\alpha }}}},$ $\Delta \beta = {\pi \mathord{\left/ {\vphantom {\pi {{{N}_{\beta }}}}} \right. \kern-0em} {{{N}_{\beta }}}}$ и W1(λ) = ${\text{ = }}\,{\text{max}}\left( {1 - \left| \lambda \right|,0} \right).$

РЕЗУЛЬТАТЫ МОДЕЛИРОВАНИЯ

Рассмотрим результаты моделирования ускорения частиц в магнитном острове размером $L = {{10}^{8}}\,\,{\text{м ,}}$ с периодом колебаний $\Theta = 600\,\,{\text{с }}$ и амплитудой колебаний магнитного поля ${{B}_{m}} = 0.75{{B}_{{x0}}} = 3.75\,\,{\text{н Т л }}{\text{.}}$ На рис. 4 показаны средние энергии вылетевших из плазмоида частиц в зависимости от их начальных энергий W0 в кэВ, отложенных по оси абсцисс. Кривые 1–6 на трех панелях соответствуют шести разным начальным точкам запуска частиц внутри МО. Можно видеть, что внутри острова имеются точки, в которых ускорение частиц большое, им соответствуют кривые 13, ускорение слабое или отсутствует (кривые 4–6). Множество начальных данных $\left[ {t_{1}^{0};t_{2}^{0}} \right] \times \Omega \left( {{{t}^{0}}} \right) \subset \mathbb{R}_{{t,{\mathbf{x}}}}^{4}$ (область $\Omega ({{t}^{0}})$ характеризует окрестность начальной точки ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.5;{\kern 1pt} 0.2{\kern 1pt} ;{\kern 1pt} 0.25} \right)$) в совокупности представляет собой “ускорительное” множество точек. В ускорительных точках частицы получают большой прирост энергии по сравнению с начальными значениями. Средние энергии вылета у них составляют $\left\langle {{{W}_{{{\text{end}}}}}} \right\rangle \left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{{0{\kern 1pt} }}}} \right) > 150\,\,{\text{к э В }}$ во всем рассматриваемом диапазоне начальных энергий. При этом центральной на промежутке “ускорительных” начальных фаз $\left[ {t_{1}^{0};t_{2}^{0}} \right] \approx \left[ {{{\Theta } \mathord{\left/ {\vphantom {{\Theta } 8}} \right. \kern-0em} 8};3{{\Theta } \mathord{\left/ {\vphantom {{\Theta } 8}} \right. \kern-0em} 8}} \right]$ является фаза ${{t}^{0}} = {{\Theta } \mathord{\left/ {\vphantom {{\Theta } 4}} \right. \kern-0em} 4}.$ Вне “ускорительного” множества начальных точек $\left( {{{t}^{0}},{{{\mathbf{x}}}^{0}}} \right) \notin \left[ {t_{1}^{0};t_{2}^{0}} \right] \times \Omega ({{t}^{0}}),$ ускорение либо отсутствовало, либо имело место только в ограниченном диапазоне малых начальных энергий и было слабым. Ускорительная область ${\Omega (}{{t}^{0}}{\text{)}}$ может менять размер во времени – она максимальна в момент ${{t}^{0}} = {{\Theta } \mathord{\left/ {\vphantom {{\Theta } 4}} \right. \kern-0em} 4}$ и уменьшается для других значений ${{t}^{0}}.$ На рис. 4 кривая 1 соответствует центру “ускорительной” области. Ход кривой показывает, что для начальных энергий ${{W}^{{0{\kern 1pt} }}} \geqslant 10$ кэВ средняя энергия вылета превышает 700 кэВ и почти не зависит от начальной энергии ${{W}^{0}}$. Кривые 2 и 3 соответствуют центральной части “ускорительной” области ${\Omega }\left( {{{t}^{0}} = {{\Theta } \mathord{\left/ {\vphantom {{\Theta } 4}} \right. \kern-0em} 4}} \right).$ Мы видим, что для частиц в этих областях средняя энергия вылета превышает, соответственно, 600 кэВ (кривая 2) и 400 кэВ (кривая 3), и также очень слабо зависит от начальных энергий ${{W}^{0}}.$ Получается, что частицы разных энергий (различающихся в разы и даже на порядки), происходящие из одной начальной точки, выходят из области контроля с близкими конечными энергиями.

Рис. 4.

Средняя энергия вылета $\left\langle {{{W}_{{{\text{end}}}}}} \right\rangle \left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{{0{\kern 1pt} }}}} \right)$ в кэВ для нескольких начальных точек в магнитном острове размером $L = {{10}^{8}}$ м, магнитным полем ${{B}_{m}} = 0.75{{B}_{{x0}}} = 3.75$ нТл и периодом колебаний $\Theta = 600$ с и начальным временем ${{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4}.$ По оси абсцисс отложены начальные энергии частиц ${{W}^{{0{\kern 1pt} }}}$ в кэВ. Кривые 1–3 соответствуют начальным точкам ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.5;0.2;0.25} \right);$ ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.4;{\kern 1pt} 0.16;{\kern 1pt} 0.2} \right)$ и ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.6;{\kern 1pt} 0.24;{\kern 1pt} 0.3} \right)$ в ускорительной области. Кривые 4–6 соответствуют начальным точкам ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.2;{\kern 1pt} 0.1{\kern 1pt} ;0.1} \right),$ ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.7;{\kern 1pt} 0.28;{\kern 1pt} 0.35} \right)$ и ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( { - 0.2; - 0.1;{\kern 1pt} - 0.1} \right).$

Рассмотрим точки, лежащие вне “ускорительной” области $\Omega \left( {{{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4}} \right).$ Кривая (5) показывает зависимость энергий вылета частиц, запущенных близко от границы “ускорительной” области из начальной точки ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.7;{\kern 1pt} 0.28;{\kern 1pt} 0.35} \right).$ Для них средние энергии вылета превышают 150 кэВ и слабо зависят от ее начальных энергий, как и для ускорительной области. Кривые 4 и 6 для двух выбранных начальных точек ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.2;{\kern 1pt} 0.1{\kern 1pt} ;0.1} \right)$ и ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( { - 0.2{\kern 1pt} ; - 0.1;{\kern 1pt} - 0.1} \right)$ и частиц с малыми энергиями ${{W}^{{0{\kern 1pt} }}} \leqslant 10$ кэВ ускоряются до средних энергий, соответственно, около 100 и 50 кэВ; при бóльших начальных энергиях ${{W}^{{0{\kern 1pt} }}}$ ускорение практически не происходит.

На рис. 5 показано в минутах среднее время $\left\langle {{{T}_{{{\text{end}}}}}} \right\rangle \left( {{{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right)$ нахождения траектории в области моделирования для тех же траекторий, что и на рис. 4. Кривые 1, 2, 3 показывают, что траекторий из внутренней части “ускорительной” области $\Omega \left( {{{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4}} \right)$ зависимость $\left\langle {{{T}_{{{\text{end}}}}}} \right\rangle \left( {{{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right)$ от начальной энергии ${{W}^{{0{\kern 1pt} }}}$ в диапазоне ${{W}^{{0{\kern 1pt} }}} \geqslant 20$ кэВ также является очень слабой, а среднее время нахождения траектории в области контроля составляет менее 2.7 мин. В то же время кривая 6 демонстрирует достаточно большое время пребывания частиц внутри области моделирования – от 8 до 13 мин. Кривые 4 и 5 показывают, что основные отличия времен пребывания частиц в системе наблюдаются в основном на малых начальных энергиях (менее 20 кэВ), при больших же они сходятся вблизи одних и тех же значений.

Рис. 5.

Среднее время $\left\langle {{{T}_{{{\text{end}}}}}} \right\rangle \left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{{0{\kern 1pt} }}}} \right)$ (в минутах) нахождения траектории в области контроля $\Pi $ для нескольких начальных точек для варианта $L = {{10}^{8}}$ м, ${{B}_{m}} = 0.75{{B}_{{x0}}} = 3.75$ нТл, периода колебаний $\Theta = 600$ с и начального времени ${{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4}.$ Обозначения кривых такие же, как на рис. 4.

Если сравнивать различные варианты колеблющихся магнитных островов, то картина следующая: по мере усиления продольного индукционного электрического поля “ускорительный” промежуток начальных фаз, при котором частицы набирают энергии выше 150 кэВ, заметно увеличивается. На рис. 6 показано моделирование ускорения частиц в МО с характерным масштабом $L = 1.5 \times {{10}^{8}}$ м, магнитным полем ${{B}_{m}} = \,3.75$ нТл и периодом колебаний $\Theta = 300$ с. По сравнению с первым вариантом просчета, описанном на рис. 4, амплитуда электрического поля увеличилась в три раза, а временной интервал ускорительной области составил $\left[ {t_{1}^{0};t_{2}^{0}} \right] \approx \left[ {{\Theta \mathord{\left/ {\vphantom {\Theta {16}}} \right. \kern-0em} {16}};{{\Theta 7} \mathord{\left/ {\vphantom {{\Theta 7} {16}}} \right. \kern-0em} {16}}} \right].$ Центр “ускорительной” области остается прежним в пространственных координатах, но сама область расширилась, а средняя энергия вылета частиц увеличилась в несколько раз по сравнению с предыдущим рассмотренным случаем. Кривая 1 на рис. 6 соответствует начальной точке ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.5;{\kern 1pt} 0.2{\kern 1pt} ;0.25} \right)$ – центру “ускорительной” области, в которой средняя энергия вылета максимальна и превышает 2 МэВ, т.е. примерно в три раза больше, чем для первой конфигурации, показанной на рис. 4. Кривая 2 соответствует ускорительной области (аналог кривой 3 на рис. 4). По сравнению с предыдущим случаем, средняя энергия вылетающих частиц из этой области увеличилась примерно в четыре раза. Кривая 3 с начальной точкой ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.4;{\kern 1pt} 0.16;{\kern 1pt} 0.2} \right)$ (аналог кривой 2 на рис. 4) демонстрирует увеличение средней энергии вылета приблизительно в 2.3 раза. Кривые 4–8, взятые вне ускорительной области острова, демонстрируют общее увеличение прироста энергии частиц примерно в 2–4 раза по сравнению с конфигурацией на рис. 4. Отметим, что в данном случае частицы существенно ускоряются даже в тех точках, где в конфигурации на рис. 4 ускорение либо отсутствовало, либо было очень слабым.

Рис. 6.

Средняя энергия вылета $\left\langle {{{W}_{{{\text{end}}}}}} \right\rangle \left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right)$ (в кэВ) для нескольких начальных точек для МО с параметрами: $L = 1.5 \times {{10}^{8}}$ м, ${{B}_{m}} = 3.75$ нТл, $\Theta = 300$ с и ${{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4}.$ Кривая 1 соответствует центру “ускорительной” области – начальной точке ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.5;{\kern 1pt} 0.2;0.25} \right),$ кривая 2 соответствует начальной точке ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.6;0.24;0.3} \right),$ кривая 3 соответствует начальной точке ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.4;{\kern 1pt} 0.16;0.2} \right),$ кривая 4 соответствует начальной точке ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.7;{\kern 1pt} 0.28;0.35} \right),$ кривая 5 соответствует начальной точке ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.2;{\kern 1pt} 0.1;0.1} \right),$ кривая 6 соответствует начальной точке ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = - \left( {0.5;0.2;0.25} \right),$ которая симметрична центру ускорительной области, кривая 7 соответствует начальной точке ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0;{\kern 1pt} 0;0} \right)$ в центре острова, а кривая 8 соответствует начальной точке ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.7;0.3;0.4} \right).$

На рис. 7 показано в минутах среднее время $\left\langle {{{T}_{{{\text{end}}}}}} \right\rangle \left( {{{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right)$ нахождения траектории в моделируемой области для тех же начальных точек, что и на рис. 6. Сравнение с рис. 5 показывает, что для начальных точек из центральной части “ускорительной” области $\Omega \left( {{{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4}} \right)$ среднее время нахождения траектории в области контроля $\Pi $ сократилось примерно в 1.6 раза, что говорит о более сильном ускорении протонов по сравнению с конфигурацией на рис. 4.

Рис. 7.

Среднее время $\left\langle {{{T}_{{{\text{end}}}}}} \right\rangle \left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right)$ (в минутах) нахождения траектории в области контроля $\Pi $ для нескольких начальных точек для варианта $L = 1.5 \times {{10}^{8}}$ м, ${{B}_{m}} = 0.75{{B}_{{x0}}} = 3.75$ нТл, периода колебаний $\Theta = 300$ с и начального времени ${{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4}.$ Обозначения кривых такие же, как на рис. 6.

В конфигурации с МО, время осцилляций которого составляет $\Theta = 300$ с (остальные параметры указаны в подписи к рис. 6) исследованы процессы ускорения в одной и той же начальной точке внутри МО, в зависимости от того, с какими фазами (начальными точками во времени), частицы начинают движение в системе. На рис. 8 показаны средние энергии вылета $\left\langle {{{W}_{{{\text{end}}}}}} \right\rangle \left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{{0{\kern 1pt} }}}} \right)$ для различных начальных моментов времени в центре “ускорительной” области Ω. В момент времени ${{t}^{0}} = 0$ (кривая 1) частицы слабо ускоряются в диапазоне начальных энергий ${{W}^{{0{\kern 1pt} }}} \leqslant 65$ кэВ. Для момента времени ${{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta {16}}} \right. \kern-0em} {16}}$ (кривая 2) происходит существенное ускорение при малых начальных энергиях ${{W}^{{0{\kern 1pt} }}} \leqslant 20$ кэВ и более слабое ускорение в остальном диапазоне. Кривые 3 и 7 демонстрируют эффективное ускорение до 700–800 кэВ и выше. Наибольший прирост энергии частицы получают при ${{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4}$(кривая 5) и ${{t}^{0}} = {{3\Theta } \mathord{\left/ {\vphantom {{3\Theta } {16}}} \right. \kern-0em} {16}};{{5\Theta } \mathord{\left/ {\vphantom {{5\Theta } {16}}} \right. \kern-0em} {16}};$ (кривые 4, 6). Кривые 8–11 (${{t}^{0}} = {{3\Theta } \mathord{\left/ {\vphantom {{3\Theta } 8}} \right. \kern-0em} 8},$ ${{7\Theta } \mathord{\left/ {\vphantom {{7\Theta } {16}}} \right. \kern-0em} {16}},$ ${\Theta \mathord{\left/ {\vphantom {\Theta 2}} \right. \kern-0em} 2},$ $3{\kern 1pt} {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4},$ $7{\Theta \mathord{\left/ {\vphantom {\Theta 8}} \right. \kern-0em} 8}$) снова показывают относительно слабое ускорение, с некоторыми различиями, особенно для малых начальных энергий до 50 кэВ.

Рис. 8.

Средняя энергия вылета $\left\langle {{{W}_{{{\text{end}}}}}} \right\rangle \left( {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right)$ (в кэВ) для различных начальных моментов времени ${{t}^{0}}$ для варианта $L = 1.5 \times {{10}^{8}}$ м, ${{B}_{m}} = 0.75{{B}_{{x0}}} = 3.75$ нТл с периодом колебаний $\Theta = 300$ с для начальной точки ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.5;{\kern 1pt} 0.2{\kern 1pt} ;0.25} \right)$ – центра “ускорительной” области. Кривые (111) соответствуют начальным моментам времени: ${{t}^{0}} = 0$1, ${{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta {16}}} \right. \kern-0em} {16}}$ 2, ${{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 8}} \right. \kern-0em} 8}$ 3, ${{t}^{0}} = {{3\Theta } \mathord{\left/ {\vphantom {{3\Theta } {16}}} \right. \kern-0em} {16}}$ 4, ${{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4}$ 5, ${{t}^{0}} = {{5\Theta } \mathord{\left/ {\vphantom {{5\Theta } {16}}} \right. \kern-0em} {16}}$ 6, ${{t}^{0}} = {{3\Theta } \mathord{\left/ {\vphantom {{3\Theta } 8}} \right. \kern-0em} 8}$ 7, ${{t}^{0}} = {{7\Theta } \mathord{\left/ {\vphantom {{7\Theta } {16}}} \right. \kern-0em} {16}}$ 8, ${{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 2}} \right. \kern-0em} 2}$ 9, ${{t}^{0}} = 3{\kern 1pt} {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4}$ 10, ${{t}^{0}} = 7{\kern 1pt} {\Theta \mathord{\left/ {\vphantom {\Theta 8}} \right. \kern-0em} 8}$ 11.

Отметим, что для других начальных моментов ${{t}^{0}} > {\Theta \mathord{\left/ {\vphantom {\Theta 2}} \right. \kern-0em} 2}$ ускорение либо отсутствует, либо является слабым. Для второй половины периода колебаний острова $\left( {{\Theta \mathord{\left/ {\vphantom {\Theta {2;\Theta }}} \right. \kern-0em} {2;\Theta }}} \right)$ также можно выделить аналогичный интервал времени, когда ускорение частиц наиболее эффективно.

Таким образом, результаты нашего моделирования показывают, что существует не только пространственный, но и временной интервал ($\left[ {t_{1}^{0};t_{2}^{0}} \right] \approx \left[ {{\Theta \mathord{\left/ {\vphantom {\Theta {16}}} \right. \kern-0em} {16}};7{\Theta \mathord{\left/ {\vphantom {\Theta {16}}} \right. \kern-0em} {16}}} \right]$), когда взаимодействие индуктивного электрического поля колеблющегося МО с частицами плазмы оказывается очень эффективным, т.е. этот механизм может играть существенную роль в дополнительной энергизации предускоренных протонов, движущихся в магнитно неоднородной среде вблизи ГТС. В некотором смысле можно говорить о резонансном характере ускорения частиц внутри колеблющегося острова, сходным с явлением параметрического резонанса.

Далее рассмотрим распределение в пространстве скоростей потоков вылетающих из острова частиц. На рис. 9 и 10 показаны функции распределения плазмы в зависимости от направляющих углов скоростей $\alpha ,\beta $ (см. формулы (62) и (63)) в сферической системе координат.

Рис. 9.

Логарифм функции распределения $\lg \left( {1 + f{\kern 1pt} \left[ {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right]{\kern 1pt} \left( {\alpha ,\beta } \right)} \right)$ (см. раздел “Методика…”) для шести начальных энергий ${{W}^{0}}$ в зависимости от направлений скоростей вылета частиц для “медленно” ($\Theta = 600$ с) колеблющегося плазмоида, соответствующего рис. 4. Начальная точка ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.5;{\kern 1pt} 0.2{\kern 1pt} ;0.25} \right)$выбрана в центре ускорительной области; начальный момент времени ${{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4}$. Функции распределения показаны для следующих энергий: ${{W}^{0}} = 5$ кэВ (панель 1), ${{W}^{0}} = 15$ кэВ (2), 25 кэВ (3), 50 кэВ (4), 75 кэВ (5), 95 кэВ (6).

Рис. 10.

Функция $\lg \left( {1 + f{\kern 1pt} \left[ {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right]{\kern 1pt} \left( {\alpha ,\beta } \right)} \right)$от функции распределения $f{\kern 1pt} \left[ {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right]{\kern 1pt} \left( {\alpha ,\beta } \right)$ направления скоростей вылета по углам вылета для варианта $L = 1.5 \times {{10}^{8}}$м, ${{B}_{m}} = 0.75{{B}_{{x0}}} = 3.75$ нТл, периода колебаний $\Theta = 300$ с для начальной точки ${{{{{\mathbf{x}}}^{0}}} \mathord{\left/ {\vphantom {{{{{\mathbf{x}}}^{0}}} L}} \right. \kern-0em} L} = \left( {0.5;{\kern 1pt} 0.2{\kern 1pt} ;0.25} \right)$ центра “ускорительной” области начальных точек, и наилучшего для ускорения начального времени ${{t}^{0}} = {\Theta \mathord{\left/ {\vphantom {\Theta 4}} \right. \kern-0em} 4}$ для шести начальных энергий ${{W}^{0}}.$ На панели (1) ${{W}^{0}} = 5$ кэВ, на панели (2) ${{W}^{0}} = 15$ кэВ, на панели (3) ${{W}^{0}} = 25$ кэВ, на панели (4) ${{W}^{0}} = 50$ кэВ, на панели (5) ${{W}^{0}} = 75$ кэВ, на панели (6) ${{W}^{0}} = 95$ кэВ.

Рис. 9 приведен для первого рассмотренного острова, осциллирующего с периодом $\Theta = 600$ с. Рис. 10 соответствует второй конфигурации с $\Theta = 300$ с. Сравнение двух рисунков демонстрирует сходство общей структуры функции распределения независимо от начальных энергий частиц в плазме. Видно, что вытекающий поток плазмы имеет явную пространственную анизотропию и имеет максимумы на определенных азимутальных и широтных углах. Если характеризовать направления вылета плазмы в целом, то частицы выбрасываются из острова в области азимутальных углов $60^\circ \leqslant \beta \leqslant 120^\circ ,$ концентрируясь вокруг $\beta \sim 90^\circ .$ Это означает, что основная зона вылета частиц сконцентрирована вблизи экваториальной области магнитного пузыря, изображенной на рис. 2. С другой стороны, распределение вылетающих частиц концентрируется в двух симметричных областях относительно угла $\alpha \sim - 30^\circ ,$ а именно $\alpha \sim 60^\circ $ и $\alpha \sim - 120^\circ .$ Это означает, что вылет частиц из острова неоднороден и происходит вблизи экватора не диффузно во все стороны, а направлен в виде двух пучков под углом 180° друг к другу. Сравнение рис. 9-1 и 9-6 показывает, что, с ростом энергий частиц картина распределений вылета становится все более сложной, границы ее размываются, и она приобретает сложный характер в пространстве скоростей.

Характер функции распределения вылетающих частиц, показанный на рис. 10 для варианта более быстро колеблющегося МО, сходен с рис. 9. Сравнение рис. 9 и 10 показывает, что для данного набора начальных точек $\left[ {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right]$ функция распределения $f\left[ {{{t}^{0}},{{{\mathbf{x}}}^{0}},{{W}^{0}}} \right]{\kern 1pt} \left( {\alpha ,\beta } \right)$ обладает ярко выраженной анизотропией. Поскольку предускоренные частицы, которые попадают в колеблющийся магнитный остров, имеют достаточно ограниченный диапазон энергий, следует ожидать, что в результате их ускорения получится также анизотропное по направлениям распределение скоростей, которое получится в результате сложения разных анизотропных распределений. Отметим, что этот тонкий вопрос заслуживает дальнейшего более детального исследования.

ОБСУЖДЕНИЕ РЕЗУЛЬТАТОВ

Результаты нашего моделирования показывают, что рассматриваемая в работе модельная система демонстрирует разнообразное и сложное поведение, схожее с параметрическим резонансом, которое может резко изменяться при небольших изменениях как сценария, так и входных параметров. Из результатов моделирования следует, что предускоренные в ГТС протоны получат дальнейшее ускорение в колеблющемся магнитном острове, если они попадут в ускорительное множество начальных положений $\left( {{{\tau }^{0}},{{{\mathbf{x}}}^{0}}} \right) \in \left[ {t_{1}^{0};t_{2}^{0}} \right] \times {\Omega }\left( {{{\tau }^{0}}} \right),$ где $\left[ {t_{1}^{0};t_{2}^{0}} \right]$ – отрезок начальных фаз, в котором происходит эффективное ускорение предускоренных протонов, а ${\Omega (}{{\tau }^{0}}{\text{)}}$ – “ускорительная” область начальных точек. При этом для достижения максимально возможной энергии вылета предускоренные протоны должны попасть в район центра ${{{\mathbf{x}}}^{{(с )}}}$ “ускорительной” области начальных точек ${\Omega (}{{\tau }^{0}}{\text{)}},$ причем в момент времени, близкий к фазе колебаний острова в четверть периода, т.е. в положение $\left( {{{\tau }^{0}},{{{\mathbf{x}}}^{0}}} \right) \approx \left( {\frac{\Theta }{4},{{{\mathbf{x}}}^{{(с )}}}} \right).$ Как видно на рис. 5 и 7, по мере отклонения от этого наилучшего для ускорения начального положения, средняя энергия вылета уменьшается в разы, и далее на порядок.

Как показывают приведенные выше оценки на основе обобщенного закона Ома (см. раздел “Геометия системы и …”), в колеблющемся магнитном острове на расстоянии 1 а. е. от Солнца для компенсации вызванных неоднородностью плазмы напряжений может возникать потенциальное электрическое поле порядка ${{E}_{p}} \sim {{10}^{{ - 6}}}$ В/м, обусловленное малым по относительной величине разделением зарядов. Максимальные значения продольного индукционного электрического поля оказываются порядка ${{E}_{\parallel }} \sim {{10}^{{ - 4}}}$ В/м, при этом пространственное распределение максимумов ${{E}_{\parallel }}$ оказывается существенно неравномерным. Взятое по верхней границе возможных значений индукционное электрическое поле, используемое в модели, позволяет имитировать одним расчетом однократного взаимодействия предускоренного протона с модельным магнитным островом процесс многократного взаимодействия частицы с островами всей системы. При взаимодействии предускоренного пересоединением в ГТС протона с одним островом, он рассеивается на неоднородностях магнитного поля и ускоряется продольным электрическим полем. Значительная часть предускоренных протонов может испытывать многократное, от 10 до 100 взаимодействий с различными островами системы, в результате чего они могут получать значительное ускорение до энергий, на несколько порядков больше начальной, аналогично коллективному процессу ускорений на магнитных островах, описанному в (Zank и др., 2014; 2015a; 2015b; Zhou и др., 2015; le Roux и др., 2015; 2016). Этот упрощенный подход позволяет выявить основные качественные особенности процесса ускорения продольным электрическим полем в магнитном острове и дает возможность сравнения результатов с наблюдениями КА.

На данный момент становится ясным, что такое поведение на качественном и количественном уровне хорошо согласуется с данными in situ измерений для потоков протонов, локально ускоренных до энергий в диапазоне от сотен кэВ до нескольких МэВ, которые были получены при пролетах космических аппаратов через магнитные острова вблизи ГТС (Khabarova и др., 2015a; 2015b; 2016). Поставленная и решенная в настоящей работе задача представляется важной, поскольку не требует никаких дополнительных источников возмущения, кроме наблюдающихся в спокойное время поверхностных волн ГТС, хорошо известных уже несколько десятилетий (Musielak и др., 1988; Ruderman, 1990; Yamauchi и др., 1997). Общность процессов распространения таких волн с гармоническим осциллятором особо отмечена в (Ruderman, 1998). Топология ГТС довольно часто предполагает наличие относительно небольших складок, в которые помещается от одного до нескольких МО с размерами порядка 0.001–0.01 а. е. Первичное ускорение энергичных частиц магнитным пересоединением на сильных токовых слоях в солнечном ветре может достигать сотен кэВ, при этом максимальные полученные энергии и траектории ускоренных частиц существенно зависят от начальных энергий и топологии магнитных полей вокруг нейтральной плоскости (Zharkova, Khabarova, 2012; 2015). В целом, ускоренные таким образом частицы создают облака различных форм вокруг токового слоя. Выявить углы и энергии влета предускоренных протонов в магнитный остров можно только в результате моделирования конкретной ситуации, наблюдающейся в солнечном ветре, что на данный момент выходит за рамки поставленной задачи, но планируется сделать в будущем. Таким образом, система “токовый слой–магнитный остров” исключительно чувствительна к топологическим факторам, чем отчасти может объясняться нерегулярность обнаружений частиц, ускоренных до МэВ рядом с токовыми слоями в солнечном ветре.

Отметим, что нами также были проведены расчеты для модели полей в колеблющемся магнитном острове, в которой электрическое поле строго ортогонально магнитному. Как и следовало ожидать, для таких полей никакого среднего ускорения или замедления не получилось. При этом у подавляющего числа траекторий энергия вылета очень мало отличается от начальной. Лишь у ничтожно малой доли траекторий энергия вылета отличается от начальной более чем на 5%, но лежит в пределах $0.7 < {{{{W}_{{{\text{end}}}}}} \mathord{\left/ {\vphantom {{{{W}_{{{\text{end}}}}}} {{{W}^{{0{\kern 1pt} }}} < 1.5}}} \right. \kern-0em} {{{W}^{{0{\kern 1pt} }}} < 1.5}}.$ Этим подчеркивается определяющая роль продольного электрического поля в ускорении заряженных частиц.

ВЫВОДЫ

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

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

Интересным фактом является то, что ускоренные частицы формируют анизотропное по направлению облако, вылетая преимущественно в районе экватора в стороны ГТС. Физически это означает, что эти частицы будут впоследствии захвачены ГТС и могут участвовать в процессе магнитного пересоединения. Этот сценарий будет рассмотрен нами в дальнейшем. Отметим, что наблюдения также показывают наличие существенной анизотропии энергичных протонов, ускоренных локально в областях, занятых магнитными островами (Khabarova и др., 2015a; 2015b; 2016).

Существенным результатом является также неодинаковая эффективность ускорения частиц с изначально низкими и высокими энергиями. Как было отмечено во введении, механизм ускорения сжатием/слиянием магнитных островов (Zank и др., 2014) предполагает более эффективное ускорение высокоэнергичных частиц, чем низкоэнергичных. Этот эффект действительно наблюдается в солнечном ветре, но ниже определенной энергии отсечки происходит инверсия, и наиболее эффективно оказываются ускоренными самые низкоэнергичные частицы. Наши результаты показывают, что, наиболее вероятно, это является проявлением другого механизма ускорения, а именно резонансной накачки, изученной в данной работе.

Таким образом, в результате исследования поведения частиц в колеблющемся магнитном острове нами получены следующие основные результаты.

1. Построена модель магнитного острова, колеблющегося в складке пересоединяющегося ГТС с периодом в несколько минут. Протестированы различные варианты условий попадания в него протонов, предускоренных до кэВ магнитным пересоединением, включая различные фазы колебаний и различное множество точек внутри острова. Для моделирования взяты максимально возможные значения продольного электрического поля (~10–4–10–5 В/м), что позволило с минимальными техническими затратами описать движение частиц в одиночном острове с конфигурацией, наиболее благоприятной для ускорения частиц, имитируя ускорение частиц в системе из достаточно большого числа магнитных островов с меньшими значениями продольного электрического поля.

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

3. Предположено, что эффекты резонансной накачки, возникающие в колеблющемся острове, ответственны за обнаруженное ранее явление так называемого инверсивного ускорения (Zank и др., 2015b), когда частицы с изначально меньшими энергиями ускоряются эффективнее, чем частицы с большими энергиями. В результате всеми частицами достигается некое пороговое значение энергии, выше которого ускорения нет. В случае системы магнитных островов в складках ГТС, эта энергия отсечки, по-видимому, находится на уровне 1–2 МэВ, что соответствует наблюдениям (Khabarova и др., 2015a; 2015b; 2016; 2017). Выше этой энергии начинает работать стохастический механизм (Zank и др., 2014; 2015a; 2015b) при котором сильнее ускоряются частицы с большими энергиями.

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

Работа выполнена при поддержке подпрограммы V.15 “Динамика разреженной плазмы в космосе и лаборатории” Комплексной программы фундаментальных исследований ОФН РАН. Работа Х.В. Маловой выполнена при поддержке гранта РНФ 14-12-00824. Работа О.В. и И.В. Мингалева и О.В. Хабаровой поддержана грантами РФФИ 16-02-00479, 17-02-00300 и 17-02-01328. Авторы выражают благодарность группе обработки данных SMEI и STELab за 3-D восстановление картин плотности, скорости и магнитного поля в гелиосфере (использованы данные сайта http:// smei.ucsd.edu/new_smei/data&images/data&images. html).

ПРИЛОЖЕНИЕ. МЕТОД РАСЧЕТА ТРАЕКТОРИЙ ПРОТОНОВ

Для сокращения вычислительных затрат и удобства системы уравнений (47) и (48) нужно представить в оптимальной безразмерной (нормализованной) форме. Для этого надо выбрать масштабы ${\kern 1pt} \hat {e} = \left| e \right|,$ $\hat {m} = {{m}_{0}},$ $\hat {v} = c,$ а также произвольное значение какого-либо одного из следующих масштабов: ${\kern 1pt} \hat {t},$ ${\kern 1pt} \hat {x},$ ${\kern 1pt} \hat {B},$ ${\kern 1pt} \hat {E},$ а остальные из них определить из соотношений ${\kern 1pt} \hat {x} = с \hat {t},$ $\hat {t} = {{{{m}_{0}}} \mathord{\left/ {\vphantom {{{{m}_{0}}} {\left( {e\hat {B}} \right)}}} \right. \kern-0em} {\left( {e\hat {B}} \right)}},$ $\hat {E} = c\hat {B}.$ В результате, если обозначить безразмерный заряд как $e{\text{'}} = {\text{sign}}\left( e \right) = \pm 1$ и по традиции обозначить безразмерный импульс через u:

${\mathbf{u}} = {\mathbf{p}}{\text{'}} = {{\mathbf{p}} \mathord{\left/ {\vphantom {{\mathbf{p}} {\left( {{{m}_{0}}c} \right)}}} \right. \kern-0em} {\left( {{{m}_{0}}c} \right)}} = {{\mathbf{v}} \mathord{\left/ {\vphantom {{\mathbf{v}} {\sqrt {{\kern 1pt} 1 - {{{{{\left| {\mathbf{v}} \right|}}^{2}}} \mathord{\left/ {\vphantom {{{{{\left| {\mathbf{v}} \right|}}^{2}}} {{{c}^{2}}{\kern 1pt} }}} \right. \kern-0em} {{{c}^{2}}{\kern 1pt} }}} }}} \right. \kern-0em} {\sqrt {{\kern 1pt} 1 - {{{{{\left| {\mathbf{v}} \right|}}^{2}}} \mathord{\left/ {\vphantom {{{{{\left| {\mathbf{v}} \right|}}^{2}}} {{{c}^{2}}{\kern 1pt} }}} \right. \kern-0em} {{{c}^{2}}{\kern 1pt} }}} }} = {{{\mathbf{v}}{\text{'}}} \mathord{\left/ {\vphantom {{{\mathbf{v}}{\text{'}}} {\sqrt {1 - {{{\left| {{\mathbf{v}}{\text{'}}} \right|}}^{2}}} }}} \right. \kern-0em} {\sqrt {1 - {{{\left| {{\mathbf{v}}{\text{'}}} \right|}}^{2}}} }},$
то получим следующие безразмерные формы систем (48) и (49) (далее штрихи опущены):

(64)
$\begin{gathered} \frac{{d{\mathbf{x}}\left( t \right)}}{{dt}} = {\mathbf{v}}\left( t \right),\,\,\,\,{\mathbf{u}}\left( {\mathbf{v}} \right) = \frac{{\mathbf{v}}}{{\sqrt {{\kern 1pt} 1 - {{{\left| {\mathbf{v}} \right|}}^{2}}} }}, \\ \frac{{d{\mathbf{u}}\left( {{\mathbf{v}}\left( t \right)} \right)}}{{d{\kern 1pt} t}} = e{\text{'}}\left( {{\mathbf{E}}\left( {{\mathbf{x}}\left( t \right),t} \right) + \left[ {{\mathbf{v}}\left( t \right) \times {\mathbf{B}}\left( {{\mathbf{x}}\left( t \right),t} \right)} \right]{\kern 1pt} } \right), \\ \end{gathered} $
(65)
$\begin{gathered} \frac{{d{\mathbf{x}}\left( t \right)}}{{dt{\kern 1pt} }} = {\mathbf{v}}\left( t \right),\,\,\,\,\frac{{d{\mathbf{v}}\left( t \right)}}{{d{\kern 1pt} t}} = \\ = e{\text{'}}\left( {{\mathbf{E}}\left( {{\mathbf{x}}\left( t \right),t} \right) + \left[ {{\mathbf{v}}\left( t \right){\mathbf{ \times B}}\left( {{\mathbf{x}}\left( t \right),t} \right)} \right]{\kern 1pt} } \right). \\ \end{gathered} $

В ходе расчета траектории для экономии вычислительных ресурсов при значении параметра $\beta = {{\left| {\mathbf{v}} \right|} \mathord{\left/ {\vphantom {{\left| {\mathbf{v}} \right|} c}} \right. \kern-0em} c} < {{\beta }_{0}} = 0.02$ (которому соответствует кинетическая энергия ${{W}_{0}} \approx 188$ кэВ) выполняется временной шаг в рамках классической системы уравнений (65), а в случае ускорения частицы при $\beta \geqslant {{\beta }_{0}} = 0.02$ выполняется временной шаг в рамках релятивистской системы уравнений (64).

Для численного решения задачи Коши (65), (50) на начальной части траектории используется разработанный авторами и показавший очень высокую эффективность алгоритм расчета траектории заряда, основанный на точном решении для нерелятивистского случая в постоянных (по пространству и времени) электрическом и магнитном полях. Этот алгоритм аналитически учитывает вращение в магнитном поле, имеет нулевую фазовую ошибку и снимает условие необходимости достаточно большого числа шагов по времени на период циклотронного вращения частиц. В случае плавно изменяющихся полей и достаточно сильного магнитного поля алгоритм допускает шаг по времени в десятки или более периодов циклотронного вращения, что позволяет получить значительную экономию вычислительных ресурсов. Разработанный нами алгоритм имеет второй порядок точности, допускает переменный шаг по времени, а также явный и неявный варианты, и основан на точном решении задачи Коши (65), (50) в случае постоянных полей ${\mathbf{E}} \equiv {\text{const}},$ ${\mathbf{B}} \equiv {\text{const}}{\text{.}}$ Это решение как функцию времени и начальных данных в обозначениях (53), используя стандартные обозначения:

(66)
$\begin{gathered} {\mathbf{b}} = \frac{{\mathbf{B}}}{{\left| {\mathbf{B}} \right|}},\,\,\,\,{{{\mathbf{v}}}_{E}} = \frac{{\left[ {{\mathbf{E}} \times {\mathbf{b}}} \right]}}{{\left| {\mathbf{B}} \right|}},\,\,\,{{\omega }_{c}} = \left| {\mathbf{B}} \right|,\,\,\,{{{\mathbf{E}}}_{{||}}} = \left( {{\mathbf{E}},{\mathbf{b}}} \right){\mathbf{b}}, \\ {{{\mathbf{v}}}_{{||}}} = \left( {{\mathbf{v}},{\mathbf{b}}} \right){\mathbf{b}},\,\,\,{{{\mathbf{v}}}_{L}} = {\mathbf{v}} - {{{\mathbf{v}}}_{{||}}} - {{{\mathbf{v}}}_{E}},\,\,\,\,\tau = t - {{t}^{0}}, \\ \end{gathered} $
можно представить в виде

(67)
$\left. \begin{gathered} {\mathbf{X}}\left( {t,{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right) = {{{\mathbf{x}}}^{0}} + \left( {{{{\mathbf{v}}}_{E}} + {\mathbf{v}}_{\parallel }^{0}} \right)\tau + \frac{{{{\tau }^{2}}}}{2}e{\text{'}}{{{\mathbf{E}}}_{\parallel }} + \hfill \\ + \,\,\frac{{{\mathbf{v}}_{L}^{0}}}{{{{\omega }_{c}}}}{\text{sin}}\left( {{{\omega }_{c}}\tau } \right) - e{\text{'}}\frac{{\left[ {{\mathbf{b}} \times {\mathbf{v}}_{L}^{0}} \right]}}{{{{\omega }_{c}}}}\left( {1 - {\text{cos}}\left( {{{\omega }_{c}}\tau } \right)} \right), \hfill \\ {\mathbf{V}}\left( {t,{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right) = {{{\mathbf{v}}}_{E}} + {\mathbf{v}}_{\parallel }^{0} + \tau e{\text{'}}{{{\mathbf{E}}}_{\parallel }} + \hfill \\ + \,\,{\mathbf{v}}_{L}^{0}{\text{cos}}\left( {{{\omega }_{c}}\tau } \right) - e{\text{'}}\left[ {{\mathbf{b}} \times {\mathbf{v}}_{L}^{0}} \right]{\kern 1pt} {\kern 1pt} {\text{sin}}\left( {{{\omega }_{c}}\tau } \right). \hfill \\ \end{gathered} \right\}$

Рассмотрим алгоритм в терминах послойного перехода ${{t}^{0}} \to {{t}^{1}} = {{t}^{0}} + \tau $ в полях ${\mathbf{E}}\left( {{\mathbf{x}},t} \right)$ и ${\mathbf{B}}{\kern 1pt} \left( {{\mathbf{x}},t} \right),$ заданных в любой точке во все моменты времени $t \geqslant {{t}^{0}}.$ Обозначим ${{t}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}} = {{t}^{0}} + {\tau \mathord{\left/ {\vphantom {\tau 2}} \right. \kern-0em} 2},$ и для функций вида $f\left( t \right)$ и $F{\kern 1pt} \left( {{\mathbf{x}},t} \right)\;$ введем обозначения ${{f}^{\alpha }} = f({{t}^{\alpha }}){\kern 1pt} $ и ${{F}^{\alpha }} = F{\kern 1pt} \left( {{\mathbf{x}}({{t}^{\alpha }}),{{t}^{\alpha }}} \right),$ $\alpha = 0,\frac{1}{2},1.$

I. СТАРТОВАЯ ИТЕРАЦИЯ

I. 1. Находим поля ${{{\mathbf{E}}}^{0}} = {\mathbf{E}}\left( {{{{\mathbf{x}}}^{0}},{{t}^{0}}} \right),$ ${{{\mathbf{B}}}^{0}} = {\mathbf{B}}\left( {{{{\mathbf{x}}}^{0}},{{t}^{0}}} \right),$ затем по ним и формулам (66), находим ${{{\mathbf{b}}}^{0}},$ ${\mathbf{v}}_{E}^{{{\kern 1pt} 0}},$ $\omega _{c}^{{{\kern 1pt} 0}},$ ${\mathbf{E}}_{{||}}^{0} = \left( {{{{\mathbf{E}}}^{0}};{{{\mathbf{b}}}^{0}}} \right){{{\mathbf{b}}}^{0}},$ ${\mathbf{v}}_{\parallel }^{0} = \left( {{{{\mathbf{v}}}^{0}};{{{\mathbf{b}}}^{0}}} \right){{{\mathbf{b}}}^{0}},$ ${\mathbf{v}}_{L}^{0} = {{{\mathbf{v}}}^{0}} - {\mathbf{v}}_{{||}}^{0} - {\mathbf{v}}_{E}^{{{\kern 1pt} 0}}.$

I. 2. По первой формуле в (67) находим ${{{\mathbf{x}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}}} = {\mathbf{X}}\left( {{{t}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}},{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right)$ в постоянных полях ${{{\mathbf{E}}}^{0}},{{{\mathbf{B}}}^{0}}\,:$

$\begin{gathered} {{{\mathbf{x}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}}} = {{{\mathbf{x}}}^{0}} + \left( {{\mathbf{v}}_{E}^{0} + {\mathbf{v}}_{\parallel }^{0}} \right)\frac{\tau }{2} + \frac{{{{\tau }^{2}}}}{8}e{\text{'}}{\mathbf{E}}_{\parallel }^{0} + \\ + \,\,\frac{{{\mathbf{v}}_{L}^{0}}}{{{{\omega }_{c}}}}{\text{sin}}\left( {\frac{{\omega _{c}^{{{\kern 1pt} 0}}\tau }}{2}} \right) - e{\text{'}}\frac{{\left[ {{{{\mathbf{b}}}^{0}} \times {\mathbf{v}}_{L}^{0}} \right]}}{{\omega _{c}^{{{\kern 1pt} 0}}}}\left( {1 - {\text{cos}}\left( {\frac{{\omega _{c}^{{{\kern 1pt} 0}}\tau }}{2}} \right)} \right). \\ \end{gathered} $

Находим поля ${{{\mathbf{E}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}}} = {\mathbf{E}}\left( {{{{\mathbf{x}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}}},{{t}^{{{\kern 1pt} {1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}} \right),$ B1/2(0) = $ = {\mathbf{B}}\left( {{{{\mathbf{x}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}}},{{t}^{{{\kern 1pt} {1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}} \right)$ а затем по ним и формулам (66), находим ${{{\mathbf{b}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}}},$ ${\mathbf{v}}_{E}^{{{\kern 1pt} {1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}},$ $\omega _{c}^{{{\kern 1pt} {\kern 1pt} {1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}},$ ${\mathbf{E}}_{{||}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}} = \left( {{{{\mathbf{E}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}}};{{{\mathbf{b}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}}}} \right){{{\mathbf{b}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}}},$ ${\mathbf{v}}_{{||}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}}$ = $\left( {{{{\mathbf{v}}}^{0}};{{{\mathbf{b}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}}}} \right){{{\mathbf{b}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}}},$ ${\mathbf{v}}_{L}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}} = {{{\mathbf{v}}}^{0}} - {\mathbf{v}}_{{||}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}} - {\mathbf{v}}_{E}^{{{\kern 1pt} {1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}}.$

I. 3. По первой формуле в (67) находим ${{{\mathbf{x}}}^{{1\left( 0 \right)}}} = {\mathbf{X}}\left( {{{t}^{1}},{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right)$ в постоянных полях ${{{\mathbf{E}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}}},$ ${{{\mathbf{B}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 0 \right)}}} \right. \kern-0em} {2\left( 0 \right)}}}}}\,:$

(68)
$\begin{gathered} {{{\mathbf{x}}}^{{1(\kappa )}}} = {{{\mathbf{x}}}^{0}} + \left( {{\mathbf{v}}_{E}^{{{1 \mathord{\left/ {\vphantom {1 {2(\kappa )}}} \right. \kern-0em} {2(\kappa )}}}} + {\mathbf{v}}_{\parallel }^{{{1 \mathord{\left/ {\vphantom {1 {2(\kappa )}}} \right. \kern-0em} {2(\kappa )}}}}} \right)\tau + \frac{{{{\tau }^{2}}}}{2}e{\text{'}}{\mathbf{E}}_{\parallel }^{{{1 \mathord{\left/ {\vphantom {1 {2(\kappa )}}} \right. \kern-0em} {2(\kappa )}}}} + \\ + \,\,\frac{{{\mathbf{v}}_{L}^{{{1 \mathord{\left/ {\vphantom {1 {2(\kappa )}}} \right. \kern-0em} {2(\kappa )}}}}}}{{\omega _{c}^{{{1 \mathord{\left/ {\vphantom {1 {2(\kappa )}}} \right. \kern-0em} {2(\kappa )}}}}}}{\text{sin}}\left( {\omega _{c}^{{{1 \mathord{\left/ {\vphantom {1 {2(\kappa )}}} \right. \kern-0em} {2(\kappa )}}}}\tau } \right) - \\ - \,\,e{\text{'}}\frac{{\left[ {{{{\mathbf{b}}}^{{{1 \mathord{\left/ {\vphantom {1 {2(\kappa )}}} \right. \kern-0em} {2(\kappa )}}}}} \times {\mathbf{v}}_{L}^{{{1 \mathord{\left/ {\vphantom {1 {2(\kappa )}}} \right. \kern-0em} {2(\kappa )}}}}} \right]}}{{\omega _{c}^{{{1 \mathord{\left/ {\vphantom {1 {2(\kappa )}}} \right. \kern-0em} {2(\kappa )}}}}}}\left( {1 - {\text{cos}}\left( {\omega _{c}^{{{1 \mathord{\left/ {\vphantom {1 {2(\kappa )}}} \right. \kern-0em} {2(\kappa )}}}}\tau } \right)} \right), \\ \end{gathered} $
где $\kappa = 0.$ Находим поля ${{{\mathbf{E}}}^{{1\left( 0 \right)}}} = {\mathbf{E}}\left( {{{{\mathbf{x}}}^{{1\left( 0 \right)}}},{{t}^{1}}} \right),$ ${{{\mathbf{B}}}^{{1\left( 0 \right)}}} = {\mathbf{B}}\left( {{{{\mathbf{x}}}^{{1\left( 0 \right)}}},{{t}^{1}}} \right).$

I. 4. Находим при $\kappa = 0$ поля ${{{\mathbf{E}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( {\kappa + 1} \right)}}} \right. \kern-0em} {2\left( {\kappa + 1} \right)}}}}} = \frac{1}{2}\left( {{{{\mathbf{E}}}^{0}} + {{{\mathbf{E}}}^{{1\left( \kappa \right)}}}} \right),$ ${{{\mathbf{B}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( {\kappa + 1} \right)}}} \right. \kern-0em} {2\left( {\kappa + 1} \right)}}}}}$ = $\frac{1}{2}\left( {{{{\mathbf{B}}}^{0}} + {{{\mathbf{B}}}^{{1\left( \kappa \right)}}}} \right),$ далее по ним и формулам в (66) находим ${{{\mathbf{b}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( {\kappa + 1} \right)}}} \right. \kern-0em} {2\left( {\kappa + 1} \right)}}}}},$ ${\mathbf{v}}_{E}^{{{\kern 1pt} {1 \mathord{\left/ {\vphantom {1 {2\left( {\kappa + 1} \right)}}} \right. \kern-0em} {2\left( {\kappa + 1} \right)}}}},$ $\omega _{c}^{{{\kern 1pt} {\kern 1pt} {1 \mathord{\left/ {\vphantom {1 {2\left( {\kappa + 1} \right)}}} \right. \kern-0em} {2\left( {\kappa + 1} \right)}}}},$ ${\mathbf{E}}_{{||}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( {\kappa + 1} \right)}}} \right. \kern-0em} {2\left( {\kappa + 1} \right)}}}},$ ${\mathbf{v}}_{{||}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( {\kappa + 1} \right)}}} \right. \kern-0em} {2\left( {\kappa + 1} \right)}}}}$ = = $\left( {{{{\mathbf{v}}}^{0}};{{{\mathbf{b}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( {\kappa + 1} \right)}}} \right. \kern-0em} {2\left( {\kappa + 1} \right)}}}}}} \right){{{\mathbf{b}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( {\kappa + 1} \right)}}} \right. \kern-0em} {2\left( {\kappa + 1} \right)}}}}},$ ${\mathbf{v}}_{L}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( {\kappa + 1} \right)}}} \right. \kern-0em} {2\left( {\kappa + 1} \right)}}}}$ = ${{{\mathbf{v}}}^{0}} - {\mathbf{v}}_{{||}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( {\kappa + 1} \right)}}} \right. \kern-0em} {2\left( {\kappa + 1} \right)}}}}$$ - \,{\mathbf{v}}_{E}^{{{\kern 1pt} {1 \mathord{\left/ {\vphantom {1 {2\left( {\kappa + 1} \right)}}} \right. \kern-0em} {2\left( {\kappa + 1} \right)}}}}.$ Затем по формуле (68) при $\kappa = 1$ находим ${{{\mathbf{x}}}^{{1(1)}}} = {\mathbf{X}}\left( {{{t}^{1}},{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right)$ в постоянных полях ${{{\mathbf{E}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 1 \right)}}} \right. \kern-0em} {2\left( 1 \right)}}}}},$ ${{{\mathbf{B}}}^{{{1 \mathord{\left/ {\vphantom {1 {2\left( 1 \right)}}} \right. \kern-0em} {2\left( 1 \right)}}}}}.$ Далее находим поля ${{{\mathbf{E}}}^{{1\left( 1 \right)}}} = {\mathbf{E}}\left( {{{{\mathbf{x}}}^{{1\left( 1 \right)}}},{{t}^{1}}} \right),$ ${{{\mathbf{B}}}^{{1\left( 1 \right)}}} = {\mathbf{B}}\left( {{{{\mathbf{x}}}^{{1\left( 1 \right)}}},{{t}^{1}}} \right),$ а также параметр точности ${{\delta }_{{{\kern 1pt} 1}}} = \left\| {{{{\mathbf{x}}}^{{1\left( 1 \right)}}} - {{{\mathbf{x}}}^{{1\left( 0 \right)}}}} \right\|.$ Если достигнута заданная точность δ: ${{\delta }_{1}} \leqslant \delta ,$ то процесс закончен. В противном случае выполняется обычная итерация.

II. Обычная итерация $\kappa \to \kappa + 1.$ По уже найденным ${{{\mathbf{x}}}^{{1(\kappa )}}},$ ${{{\mathbf{E}}}^{{1\left( \kappa \right)}}},$ ${{{\mathbf{B}}}^{{1\left( \kappa \right)}}}$ выполняя действия, описанные в пункте I. 4, находим ${{{\mathbf{x}}}^{{1(\kappa + 1)}}},$ ${{{\mathbf{E}}}^{{1\left( {\kappa + 1} \right)}}},$ ${{{\mathbf{B}}}^{{1\left( {\kappa + 1} \right)}}},$ а также параметр относительной точности ${{\delta }_{{{\kern 1pt} \kappa + 1}}} = \left\| {{{{\mathbf{x}}}^{{1(\kappa + 1)}}} - {{{\mathbf{x}}}^{{1(\kappa )}}}} \right\|$ и параметр сходимости ${{\gamma }_{{\kappa + 1}}} = {{\delta }_{{\kappa + 1}}} - {{\delta }_{\kappa }}.$ Итерации выполняются до достижения заданной относительной точности ${{\delta }_{{\kappa + 1}}} \leqslant \delta .$ При этом, если условие сходимости итерационного процесса ${{\gamma }_{{\kappa + 1}}} = {{\delta }_{{\kappa + 1}}} - {{\delta }_{\kappa }} < 0$ на некоторой итерации нарушается, или количество итераций превышает заданное максимально допустимое число ${{\kappa }_{0}}$, то текущий шаг по времени τ уменьшается в два раза: $\tau \to {\tau \mathord{\left/ {\vphantom {\tau 2}} \right. \kern-0em} 2},$ и описанный процесс начинается сначала.

III. Блок завершения шага. Выполняем присваивание ${{{\mathbf{x}}}^{1}} = {{{\mathbf{x}}}^{{1\left( {\kappa + 1} \right)}}},$ ${{{\mathbf{E}}}^{1}} = {{{\mathbf{E}}}^{{1\left( {\kappa + 1} \right)}}},$ ${{{\mathbf{B}}}^{1}} = {{{\mathbf{B}}}^{{1\left( {\kappa + 1} \right)}}}.$Затем находим поля ${{{\mathbf{E}}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}} = \frac{1}{2}\left( {{{{\mathbf{E}}}^{0}} + {{{\mathbf{E}}}^{1}}} \right),$ ${{{\mathbf{B}}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}} = \frac{1}{2}\left( {{{{\mathbf{B}}}^{0}} + {{{\mathbf{B}}}^{1}}} \right),$ далее по ним и формулам в (66) находим ${{{\mathbf{b}}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}},$ ${\mathbf{v}}_{E}^{{{\kern 1pt} {1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}},$ $\omega _{c}^{{{\kern 1pt} {1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}},$ ${\mathbf{E}}_{{||}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}},$ ${\mathbf{v}}_{\parallel }^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}} = \left( {{{{\mathbf{v}}}^{0}};{{{\mathbf{b}}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}} \right){{{\mathbf{b}}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}},$ ${\mathbf{v}}_{L}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}} = {{{\mathbf{v}}}^{0}} - {\mathbf{v}}_{\parallel }^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}} - {\mathbf{v}}_{E}^{{{\kern 1pt} {1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}.$ Затем по второй формуле в (67) находим ${{{\mathbf{v}}}^{1}} = {\mathbf{V}}\left( {t,{{t}^{0}},{{{\mathbf{x}}}^{0}},{{{\mathbf{v}}}^{0}}} \right)$ в постоянных полях ${{{\mathbf{E}}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}},$ ${{{\mathbf{B}}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}\,:$

$\begin{gathered} {{{\mathbf{v}}}^{1}} = {\mathbf{v}}_{E}^{{{\kern 1pt} {1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}} + {\mathbf{v}}_{{||}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}} + \tau e{\text{'}}{\mathbf{E}}_{{||}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}} + {\mathbf{v}}_{L}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}{\kern 1pt} {\text{cos}}\left( {\omega _{c}^{{{\kern 1pt} {1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}\tau } \right) - \\ - \,\,e{\text{'}}\left[ {{{{\mathbf{b}}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}} \times {\mathbf{v}}_{L}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}} \right]{\kern 1pt} {\kern 1pt} {\text{sin}}\left( {\omega _{c}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}\tau } \right). \\ \end{gathered} $

На этом шаг по времени закончен.

Для численного решения задачи Коши (64), (50) в релятивистском случае авторами разработана неявная схема 2 -го порядка точности. Для ее формулировки в терминах послойного перехода ${{t}^{0}} \to {{t}^{1}} = {{t}^{0}} + \tau $ необходимо ввести функцию

(69)
$\Gamma \left( {\mathbf{u}} \right) = \frac{1}{{\gamma \left( {\mathbf{v}} \right)}} = \frac{1}{{\sqrt {1 + {{{\left| {\mathbf{u}} \right|}}^{2}}} }} = \sqrt {{\kern 1pt} 1 - {{{\left| {\mathbf{v}} \right|}}^{2}}} ,$
и представить систему (64) в виде
(70)
$\begin{gathered} \frac{{d{\mathbf{x}}\left( t \right)}}{{dt}} = \Gamma \left( {{\mathbf{u}}\left( t \right)} \right){\mathbf{u}}\left( t \right),\,\,\,\,\frac{{d{\mathbf{u}}\left( t \right)}}{{dt}} = \\ = e{\text{'}}\left( {{\mathbf{E}}\left( {{\mathbf{x}}\left( t \right),t} \right) + \Gamma \left( {{\mathbf{u}}\left( t \right)} \right)\left[ {{\mathbf{u}}\left( t \right){\mathbf{ \times B}}\left( {{\mathbf{x}}\left( t \right),t} \right)} \right]{\kern 1pt} } \right), \\ \end{gathered} $
а также использовать вытекающие из нее следующие уравнения:

(71)
$\frac{{d\Gamma \left( {{\mathbf{u}}\left( t \right)} \right)}}{{d{\kern 1pt} t}} = - e{\text{'}}{{\left( {\Gamma \left( {{\mathbf{u}}\left( t \right)} \right)} \right)}^{3}}\left( {{\mathbf{u}}\left( t \right);{\mathbf{E}}\left( {{\mathbf{x}}\left( t \right),t} \right)} \right),$
(72)
$\begin{gathered} \frac{{{{d}^{2}}{\mathbf{x}}\left( t \right)}}{{d{{t}^{2}}}} = e{\text{'}}\Gamma \left( {{\mathbf{u}}\left( t \right)} \right)\left( {{\mathbf{E}}\left( {{\mathbf{x}}\left( t \right),t} \right) + } \right. \\ \left. { + \,\,\Gamma \left( {{\mathbf{u}}\left( t \right)} \right)\left[ {{\mathbf{u}}\left( t \right) \times {\mathbf{B}}\left( {{\mathbf{x}}\left( t \right),t} \right)} \right]} \right) + \,\,\frac{{d\Gamma \left( {{\mathbf{u}}\left( t \right)} \right)}}{{dt}}{\mathbf{u}}\left( t \right). \\ \end{gathered} $

Схема основана на следующих конечноразностных соотношениях:

(73)
$\begin{gathered} \frac{{{{{\mathbf{x}}}^{1}} - {{{\mathbf{x}}}^{0}}}}{\tau } = \frac{1}{2}\left( {{{\Gamma }^{0}}{{{\mathbf{u}}}^{0}} + {{\Gamma }^{1}}{{{\mathbf{u}}}^{1}}} \right) - \\ - \,\,\frac{\tau }{{12}}\left( {{{{\left( {\frac{{{{d}^{2}}{\mathbf{x}}}}{{d{{t}^{2}}}}} \right)}}^{1}} - {{{\left( {\frac{{{{d}^{2}}{\mathbf{x}}}}{{d{{t}^{2}}}}} \right)}}^{0}}} \right) + O\left( {{{\tau }^{4}}} \right), \\ \end{gathered} $
(74)
$\frac{{{{{\mathbf{u}}}^{1}} - {{{\mathbf{u}}}^{0}}}}{\tau } = \frac{1}{2}\left( {{{{\left( {\frac{{d{\mathbf{u}}}}{{dt}}} \right)}}^{0}} + {{{\left( {\frac{{d{\mathbf{u}}}}{{dt}}} \right)}}^{1}}} \right) + O\left( {{{\tau }^{2}}} \right).$

I. СТАРТОВАЯ ИТЕРАЦИЯ

I. 1. Находим поля ${{{\mathbf{E}}}^{0}} = {\mathbf{E}}\left( {{{{\mathbf{x}}}^{0}},{{t}^{0}}} \right),$ ${{{\mathbf{B}}}^{0}} = {\mathbf{B}}\left( {{{{\mathbf{x}}}^{0}},{{t}^{0}}} \right),$ затем по ним и формулам (69)(72) находим ${{\left( {\frac{{d{\mathbf{u}}}}{{dt}}} \right)}^{0}},$ ${{\left( {\frac{{d\Gamma }}{{dt}}} \right)}^{0}},$ ${{\left( {\frac{{{{d}^{2}}{\mathbf{x}}}}{{d{{t}^{2}}}}} \right)}^{0}}.$

I. 2. Находим ${{{\mathbf{x}}}^{{1\left( 0 \right)}}}$ с точностью $O({{\tau }^{3}})$ и ${{\Gamma }^{{1( - )}}}$ с точностью $O({{\tau }^{2}})$ при помощи экстраполяции:

$\begin{gathered} {{{\mathbf{x}}}^{{1\left( 0 \right)}}} = {{{\mathbf{x}}}^{0}} + \tau {{\Gamma }^{0}}{{{\mathbf{u}}}^{0}} + \frac{{{{\tau }^{2}}}}{2}{{\left( {\frac{{{{d}^{2}}{\mathbf{x}}}}{{d{{t}^{2}}}}} \right)}^{0}} + O({{\tau }^{3}}), \\ {{\Gamma }^{{1( - )}}} = {{\Gamma }^{0}} + \tau {{\left( {\frac{{d\Gamma }}{{dt}}} \right)}^{0}} + O({{\tau }^{2}}). \\ \end{gathered} $

Далее, используя ${{{\mathbf{x}}}^{{1\left( 0 \right)}}},$ находим поля ${{{\mathbf{E}}}^{{1(0)}}} = {\mathbf{E}}\left( {{{{\mathbf{x}}}^{{1\left( 0 \right)}}},{{t}^{1}}} \right),$ ${{{\mathbf{B}}}^{{1\left( 0 \right)}}} = {\mathbf{B}}\left( {{{{\mathbf{x}}}^{{1\left( 0 \right)}}},{{t}^{1}}} \right)$ с точностью $O({{\tau }^{3}}).$

I. 3. Находим ${{{\mathbf{u}}}^{{1(0)}}}$с точностью $O({{\tau }^{3}})$следующим образом. Из соотношения (74) и второго уравнения в (70) можно получить следующее векторное уравнение относительно ${{{\mathbf{u}}}^{1}}\,:$

(75)
$\left. \begin{gathered} {{{\mathbf{u}}}^{1}} + \frac{\tau }{2}e{\text{'}}{{\Gamma }^{1}}\left[ {{{{\mathbf{B}}}^{1}}{\mathbf{ \times }}{{{\mathbf{u}}}^{1}}} \right] = {\mathbf{W}} + O({{\tau }^{3}}),\,\,\,\,{\text{г д е }} \hfill \\ {\mathbf{W}} = {{{\mathbf{W}}}^{0}} + {{{\mathbf{W}}}^{1}},\,\,\,{{{\mathbf{W}}}^{0}} = {{{\mathbf{u}}}^{0}} + \frac{\tau }{2}{{\left( {\frac{{d{\mathbf{u}}}}{{dt}}} \right)}^{0}},\,\,\,{{{\mathbf{W}}}^{1}} = \frac{\tau }{2}e{\text{'}}{{{\mathbf{E}}}^{1}}. \hfill \\ \end{gathered} \right\}$

Это уравнение можно представить в форме системы линейных уравнений с матрицей $(\widehat I + \hat {\omega })\,:$

(76)
$\begin{gathered} \left( {\widehat I + \hat {\omega }} \right){\mathbf{u}} = {\mathbf{W}}\;,\,\,\,\, \\ {\text{г д е }}\,\,\,\,\hat {I} - \,\,{е д и н и ч н а я м а т р и ц а и \;}\,\,\hat {\omega }{\mathbf{u}} = \left[ {{\mathbf{\omega }} \times {\mathbf{u}}} \right]{\text{.}} \\ \end{gathered} $

Обратная матрица ${{(\widehat I + \hat {\omega })}^{{ - 1}}}$ имеет вид

$\begin{gathered} {{\left( {\widehat I + \hat {\omega }} \right)}^{{ - 1}}} = \frac{1}{{1 + {{{\left| {\mathbf{\omega }} \right|}}^{2}}}}\left( {\widehat I - \hat {\omega } + {\mathbf{\omega }} \otimes {\mathbf{\omega }}} \right),\,\,\quad{\text{г д е }}\quad \\ {\mathbf{\omega }} \otimes {\mathbf{\omega }} - \,\,{\text{д и а д а ,}}\,\,{\text{о б р а з о в а н н а я в е к т о р о м }}\,\,{\mathbf{\omega }}. \\ \end{gathered} $

Поэтому решение векторного уравнения (76) дается формулой

(77)
$\begin{gathered} {\mathbf{u}} = {{(\widehat I + \hat {\omega })}^{{ - 1}}}{\mathbf{W}} = \\ = \frac{1}{{1 + {{{\left| {\mathbf{\omega }} \right|}}^{2}}}}\left( {{\mathbf{W}} - \left[ {{\mathbf{\omega }} \times {\mathbf{W}}} \right] + \left( {{\mathbf{W}};{\mathbf{\omega }}} \right){\mathbf{\omega }}} \right). \\ \end{gathered} $

Применяя эту формулу к уравнению (76) при ${\mathbf{\omega }} = \frac{\tau }{2}{{\Gamma }^{1}}{{{\mathbf{B}}}^{1}},$ получим выражение для ${{{\mathbf{u}}}^{1}}$ через остальные слагаемые:

(78)
$\begin{gathered} {{{\mathbf{u}}}^{1}} = \\ = \frac{{{\mathbf{W}} - \frac{\tau }{2}e{\text{'}}{{\Gamma }^{1}}\left[ {{{{\mathbf{B}}}^{1}} \times {\mathbf{W}}} \right] + {{{\left( {\frac{\tau }{2}{{\Gamma }^{1}}} \right)}}^{2}}\left( {{\mathbf{W}};{{{\mathbf{B}}}^{1}}} \right){{{\mathbf{B}}}^{1}}}}{{\left( {1 + {{{\left( {\frac{\tau }{2}{{\Gamma }^{1}}\left| {{{{\mathbf{B}}}^{1}}} \right|} \right)}}^{2}}} \right)}} + O({{\tau }^{3}}). \\ \end{gathered} $

Подставляя в эту формулу уже найденные ${{\Gamma }^{{1( - )}}},$ ${{{\mathbf{B}}}^{{1\left( 0 \right)}}},$ ${{{\mathbf{E}}}^{{1(0)}}},$ получим следующие формулы:

(79)
$\left. \begin{gathered} {\mathbf{W}} = {{{\mathbf{W}}}^{0}} + {{{\mathbf{W}}}^{{1(0)}}},\,\,\,\,{{{\mathbf{W}}}^{{1(0)}}} = \frac{\tau }{2}e{\text{'}}{{{\mathbf{E}}}^{{1(0)}}}, \\ {{{\mathbf{u}}}^{{1(0)}}} = \frac{{{\mathbf{W}} - \frac{\tau }{2}e{\text{'}}{{\Gamma }^{{1( - )}}}\left[ {{{{\mathbf{B}}}^{{1(0)}}} \times {\mathbf{W}}} \right] + {{{\left( {\frac{\tau }{2}{{\Gamma }^{{1( - )}}}} \right)}}^{2}}\left( {{\mathbf{W}};{{{\mathbf{B}}}^{{1(0)}}}} \right){{{\mathbf{B}}}^{{1(0)}}}}}{{\left( {1 + {{{\left( {\frac{\tau }{2}{{\Gamma }^{{1( - )}}}\left| {{{{\mathbf{B}}}^{{1(0)}}}} \right|} \right)}}^{2}}} \right)}} + O({{\tau }^{3}}). \\ \end{gathered} \right\}$

I. 4. Используя полученные ${{{\mathbf{u}}}^{{1(0)}}}$ и ${{{\mathbf{B}}}^{{1\left( 0 \right)}}},$ ${{{\mathbf{E}}}^{{1(0)}}},$ с точностью $O({{\tau }^{3}})$ находим ${{\Gamma }^{{1(0)}}}$ по формуле (69), а также ${{\left( {\frac{{d{\mathbf{u}}}}{{dt}}} \right)}^{{1(0)}}}$ по формуле (70) и ${{\left( {\frac{{{{d}^{2}}{\mathbf{x}}}}{{d{{t}^{2}}}}} \right)}^{{1(0)}}}$ по формуле (72). Присваиваем ${{\delta }_{0}} = 1.$

II. Обычная итерация $\kappa \to \kappa + 1.$ Пусть уже известны ${{{\mathbf{x}}}^{{1\left( \kappa \right)}}},$ ${{{\mathbf{u}}}^{{1(\kappa )}}},$ ${{\Gamma }^{{1(\kappa )}}},$ ${{{\mathbf{B}}}^{{1\left( \kappa \right)}}},$ ${{{\mathbf{E}}}^{{1(\kappa )}}},$ ${{\left( {\frac{{d{\mathbf{u}}}}{{dt}}} \right)}^{{1(\kappa )}}},$ ${{\left( {\frac{{{{d}^{2}}{\mathbf{x}}}}{{d{{t}^{2}}}}} \right)}^{{1(\kappa )}}}$ с точностью $O({{\tau }^{3}}).$

II. 1. Находим ${{{\mathbf{x}}}^{{1\left( {\kappa + 1} \right)}}}$ с точностью $O({{\tau }^{4}})$ из соотношения (73):

$\begin{gathered} {{{\mathbf{x}}}^{{1\left( {\kappa + 1} \right)}}} = {{{\mathbf{x}}}^{0}} + \frac{\tau }{2}\left( {{{\Gamma }^{0}}{{{\mathbf{u}}}^{0}} + {{\Gamma }^{{1\left( \kappa \right)}}}{{{\mathbf{u}}}^{{1\left( \kappa \right)}}}} \right) - \\ - \,\,\frac{{{{\tau }^{2}}}}{{12}}\left( {{{{\left( {\frac{{{{d}^{2}}{\mathbf{x}}}}{{d{{t}^{2}}}}} \right)}}^{{1\left( \kappa \right)}}} - {{{\left( {\frac{{{{d}^{2}}{\mathbf{x}}}}{{d{{t}^{2}}}}} \right)}}^{0}}} \right) + O\left( {{{\tau }^{4}}} \right). \\ \end{gathered} $

Далее находим поля ${{{\mathbf{E}}}^{{1(\kappa + 1)}}} = {\mathbf{E}}\left( {{{{\mathbf{x}}}^{{1\left( {\kappa + 1} \right)}}},{{t}^{1}}} \right),$ ${{{\mathbf{B}}}^{{1\left( {\kappa + 1} \right)}}} = {\mathbf{B}}\left( {{{{\mathbf{x}}}^{{1\left( {\kappa + 1} \right)}}},{{t}^{1}}} \right)$ с точностью $O({{\tau }^{4}}).$

II. 2. Находим ${{{\mathbf{u}}}^{{1\left( {\kappa + 1} \right)}}}$ с точностью $O({{\tau }^{3}})$ по формуле (74) полностью аналогично (79):

$\begin{gathered} {\mathbf{W}} = {{{\mathbf{W}}}^{0}} + {{{\mathbf{W}}}^{{1(\kappa + 1)}}},\,\,\,\,{{{\mathbf{W}}}^{{1(\kappa + 1)}}} = \frac{\tau }{2}e{\text{'}}{{{\mathbf{E}}}^{{1(\kappa + 1)}}}, \\ {{{\mathbf{u}}}^{{1(\kappa + 1)}}} = \frac{{{\mathbf{W}} - \frac{\tau }{2}e{\text{'}}{{\Gamma }^{{1(\kappa )}}}\left[ {{{{\mathbf{B}}}^{{1(\kappa + 1)}}}{\mathbf{ \times W}}} \right] + {{{\left( {\frac{\tau }{2}{{\Gamma }^{{1(\kappa )}}}} \right)}}^{2}}\left( {{{{\mathbf{B}}}^{{1(\kappa + 1)}}};{\mathbf{W}}} \right){{{\mathbf{B}}}^{{1(\kappa + 1)}}}}}{{\left( {1 + {{{\left( {\frac{\tau }{2}{{\Gamma }^{{1(\kappa )}}}\left| {{{{\mathbf{B}}}^{{1(\kappa + 1)}}}} \right|} \right)}}^{2}}} \right)}} + O({{\tau }^{3}}). \\ \end{gathered} $

II. 3. Используя полученные ${{{\mathbf{u}}}^{{1(\kappa + 1)}}}$ и ${{{\mathbf{B}}}^{{1\left( {\kappa + 1} \right)}}},$ ${{{\mathbf{E}}}^{{1(\kappa + 1)}}},$ с точностью $O({{\tau }^{3}})$ находим ${{\Gamma }^{{1(\kappa + 1)}}}$ по формуле (69), а также ${{\left( {\frac{{d{\mathbf{u}}}}{{dt}}} \right)}^{{1(\kappa + 1)}}}$ по формуле (70) и ${{\left( {\frac{{{{d}^{2}}{\mathbf{x}}}}{{d{{t}^{2}}}}} \right)}^{{1(\kappa + 1)}}}$ по формуле (72).

II. 4. Находим параметр относительной точности ${{\delta }_{{{\kern 1pt} \kappa + 1}}} = \left\| {{{{\mathbf{x}}}^{{1(\kappa + 1)}}} - {{{\mathbf{x}}}^{{1(\kappa )}}}} \right\|$ и параметр сходимости ${{\gamma }_{{\kappa + 1}}} = {{\delta }_{{\kappa + 1}}} - {{\delta }_{\kappa }}.$ Итерации выполняются до достижения заданной относительной точности ${{\delta }_{{\kappa + 1}}} \leqslant \delta .$ Если условие сходимости итерационного процесса ${{\gamma }_{{\kappa + 1}}} = {{\delta }_{{\kappa + 1}}} - {{\delta }_{\kappa }} < 0$ на некоторой итерации нарушается, или количество итераций превышает заданное максимально допустимое число ${{\kappa }_{0}},$ то текущий шаг по времени τ уменьшается в два раза: $\tau \to {\tau \mathord{\left/ {\vphantom {\tau 2}} \right. \kern-0em} 2},$ и итерационный процесс начинается сначала. В конце выполняем присваивание ${{{\mathbf{x}}}^{1}} = {{{\mathbf{x}}}^{{1\left( {\kappa + 1} \right)}}},$ ${{{\mathbf{u}}}^{1}} = {{{\mathbf{u}}}^{{1(\kappa + 1)}}},$ ${{\Gamma }^{1}} = {{\Gamma }^{{1(\kappa + 1)}}},$ ${{{\mathbf{E}}}^{1}} = {{{\mathbf{E}}}^{{1\left( {\kappa + 1} \right)}}},$ ${{{\mathbf{B}}}^{1}} = {{{\mathbf{B}}}^{{1\left( {\kappa + 1} \right)}}},$ ${{\left( {\frac{{d{\mathbf{u}}}}{{dt}}} \right)}^{1}} = {{\left( {\frac{{d{\mathbf{u}}}}{{dt}}} \right)}^{{1(\kappa + 1)}}},$ ${{\left( {\frac{{{{d}^{2}}{\mathbf{x}}}}{{d{{t}^{2}}}}} \right)}^{1}} = {{\left( {\frac{{{{d}^{2}}{\mathbf{x}}}}{{d{{t}^{2}}}}} \right)}^{{1(\kappa + 1)}}}.$

Отметим, что изложенный алгоритм расчета решения задачи Коши (49), (50) в нерелятивистском случае широко используется нами с 2000-го года (см., например, Бородачев и др., 2003; Мингалев и др., 2006; 2007; 2009; 2012; Malova и др., 2015), был всесторонне протестирован и на практике доказал свою правильность и эффективность. Поэтому необходимо протестировать только новую неявную схему для релятивистского случая. Для этого тестирования лучше всего подходит точное решение задачи Коши (47), (50) в случае постоянных параллельных друг другу магнитного и электрического полей, когда заряд испытывает ускорение электрическим полем. Отметим, что это решение обсуждается в монографии (Ландау, Лифшиц, 2012), но там рассматривается только частный случай, когда начальная скорость заряда ортогональна магнитному полю, а также выбрана неудачная параметризация, которая дает для решения только неявное соотношение, и не позволяет получить решение в явной форме. Авторам не удалось найти в доступных источниках обсуждаемое точное решение. Поэтому мы получили его сами в явной векторной форме и для произвольных начальных условий.

Рассмотрим задачу Коши для безразмерной системы (69), (70) в случае постоянных и параллельных полей ${\mathbf{B}} \equiv {\text{const,}}$ ${\mathbf{E}} \equiv {{E}_{\parallel }}{\mathbf{b}} \equiv {\text{const}}$ (где ${\mathbf{b}} = {{\mathbf{B}} \mathord{\left/ {\vphantom {{\mathbf{B}} {\left| {\mathbf{B}} \right|}}} \right. \kern-0em} {\left| {\mathbf{B}} \right|}}$). Ее можно представить в виде

(80)
$\begin{gathered} \frac{{d{\mathbf{x}}\left( t \right)}}{{dt}} = \frac{{{\mathbf{u}}\left( t \right)}}{{\sqrt {1 + {{{\left| {{\mathbf{u}}\left( t \right)} \right|}}^{2}}} }},\,\,\,\,\frac{{d{\mathbf{u}}\left( t \right)}}{{dt}} = \\ = e{\text{'}}\left( {{{E}_{\parallel }}{\mathbf{b}} + \frac{{\left[ {{\mathbf{u}}\left( t \right) \times {\mathbf{B}}} \right]}}{{\sqrt {{\kern 1pt} 1 + {{{\left| {{\mathbf{u}}\left( t \right)} \right|}}^{2}}} }}} \right),\,\,\,{\mathbf{x}}\left( {{{t}^{0}}} \right) = {{{\mathbf{x}}}^{0}},\,\,\,{\mathbf{u}}\left( {{{t}^{0}}} \right) = {{{\mathbf{u}}}^{0}}. \\ \end{gathered} $

Обозначим

(81)
$\begin{gathered} \Gamma \left( {\mathbf{u}} \right) = \sqrt {1 + {{{\left| {\mathbf{u}} \right|}}^{2}}} = \frac{1}{{\sqrt {1 - {{{\left| {\mathbf{v}} \right|}}^{2}}} }},\,\,\,{{u}_{\parallel }} = \left( {{\mathbf{u}};{\mathbf{b}}} \right), \\ {{{\mathbf{u}}}_{\parallel }} = \left( {{\mathbf{u}};{\mathbf{b}}} \right){\mathbf{b}},\,\,\,\,{{{\mathbf{u}}}_{ \bot }} = {\mathbf{u}} - {{{\mathbf{u}}}_{\parallel }},\,\,\,\,\tau = t - {{t}^{0}}.{\kern 1pt} \\ \end{gathered} $

Полученное нами точное решение задачи (80) удобно представить в следующей форме. Для фактора Лоренца получается формула

(82)
$\begin{gathered} \Gamma \left( {\mathbf{u}} \right) = \gamma \left( \tau \right) = \sqrt {1 + {{{\left| {{\mathbf{u}}_{ \bot }^{0}} \right|}}^{2}} + {{{\left( {u_{\parallel }^{0} + e{\text{'}}{{E}_{\parallel }}\tau } \right)}}^{2}}} = \\ = \sqrt {{{{\left( {{{\gamma }^{0}}} \right)}}^{2}} + e{\text{'}}{{E}_{\parallel }}\tau \left( {2u_{\parallel }^{0} + e{\text{'}}{{E}_{\parallel }}\tau } \right)} ,\,\,\,\,{\text{г д е }}\,\,\,\,{{\gamma }^{0}} = \sqrt {{\kern 1pt} 1 + {{{\left| {{{{\mathbf{u}}}^{0}}} \right|}}^{2}}} . \\ \end{gathered} $

Для фазы циклотронного вращения получаются формулы

(83)
$\begin{gathered} \frac{{d\Phi \left( \tau \right)}}{{d\tau }} = \frac{{\left| {\mathbf{B}} \right|}}{{\gamma \left( \tau \right)}}, \\ \Phi \left( \tau \right) = e{\text{'}}\frac{{\left| {\mathbf{B}} \right|}}{{{{E}_{\parallel }}}}{\text{ln}}\left( {\frac{{u_{\parallel }^{0} + e{\text{'}}{{E}_{\parallel }}\tau + \gamma \left( \tau \right)}}{{u_{\parallel }^{0} + {{\gamma }^{0}}}}} \right). \\ \end{gathered} $

Для безразмерного импульса получается формула, аналогичная формуле для скорости в (67):

(84)
$\begin{gathered} {\mathbf{u}}\left( t \right) = {\mathbf{u}}_{\parallel }^{0} + e{\text{'}}{{E}_{\parallel }}\tau {\mathbf{b}} + \\ + \,\,{\mathbf{u}}_{ \bot }^{0}{\kern 1pt} {\text{cos}}\left( {\Phi \left( \tau \right)} \right) - e{\text{'}}\left[ {{\mathbf{b}} \times {\mathbf{u}}_{ \bot }^{0}} \right]{\kern 1pt} {\text{sin}}\left( {\Phi \left( \tau \right)} \right). \\ \end{gathered} $

Для координаты получается формула

(85)
$\begin{gathered} {\mathbf{x}}\left( t \right) = {{{\mathbf{x}}}^{0}} + \frac{{e{\text{'}}}}{{{{E}_{\parallel }}}}\left( {\gamma \left( \tau \right) - {{\gamma }^{0}}} \right){\mathbf{b}} + \\ + \,\,\frac{{{\mathbf{u}}_{ \bot }^{0}}}{{\left| {\mathbf{B}} \right|}}{\text{sin}}\left( {\Phi \left( \tau \right)} \right) - e{\text{'}}\frac{{\left[ {{\mathbf{b}} \times {\mathbf{u}}_{ \bot }^{0}} \right]}}{{\left| {\mathbf{B}} \right|}}{\kern 1pt} \left( {1 - {\text{cos}}\left( {\Phi \left( \tau \right)} \right)} \right). \\ \end{gathered} $

Были проведены расчеты решения задачи Коши (80) по изложенной выше неявной схеме 2 -го порядка (73), (74) в широком диапазоне релятивистских начальных условий и значений параметра ${{{{E}_{\parallel }}} \mathord{\left/ {\vphantom {{{{E}_{\parallel }}} {\left| {\mathbf{B}} \right|}}} \right. \kern-0em} {\left| {\mathbf{B}} \right|}}$ на время $T = 40{{\Theta }_{0}},$ где ${{\Theta }_{0}} = {{2\pi } \mathord{\left/ {\vphantom {{2\pi } {\frac{{d\Phi \left( 0 \right)}}{{d\tau }}}}} \right. \kern-0em} {\frac{{d\Phi \left( 0 \right)}}{{d\tau }}}} = \frac{{2\pi {{\gamma }^{0}}}}{{\left| {\mathbf{B}} \right|}}.$ Эти расчеты показали, что при достаточно малом шаге по времени численное решение с высокой точностью соответствует точному решению (82)–(85).

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

  1. Бородачев Л.В., Мингалев И.В., Мингалев О.В. Дрейфовый алгоритм расчета движения заряда в дарвиновской модели плазмы // ЖВМ и МФ. 2003. Т. 43. № 3. С. 467–480.

  2. Ландау Л.Д., Лифшиц Е.М. Теория поля. Теоретическая физика, том II. Изд. 8-е стереотипное. М.: Физматлит, 2012. 536 с.

  3. Лифшиц Е.М., Питаевский Л.П. Физическая кинетика. Теоретическая физика. Т. X. М.: Наука, 1979. 528 с.

  4. Мингалев О.В., Мингалев И.В., Мингалев В.С. Двумерное численное моделирование динамики мелкомасштабной неоднородности в околоземной плазме // Космич. исслед. 2006. Т. 44. № 5. С. 416−427.

  5. Мингалев О.В., Мингалев И.В., Малова Х.В., Зеленый Л.М. Численное моделирование плазменного равновесия в одномерном токовом слое с ненулевой нормальной компонентой магнитного поля // Физика плазмы. 2007. Т. 33. № 11. С. 1028−1041.

  6. Мингалев О.В., Мингалев И.В., Малова Х.В., Зеленый Л.М., Артемьев А.В. Несимметричные конфигурации тонкого токового слоя с постоянной нормальной компонентой магнитного поля // Физика плазмы. 2009. Т. 35. № 1. С. 85−96.

  7. Мингалев О.В., Мингалев И.В., Мельник М.Н., Артемьев А.В., Малова Х.В., Попов В.Ю., Шен Чао, Зеленый Л.М. Кинетические модели токовых слоев с широм магнитного поля // Физика плазмы. 2012. Т. 38. № 4. С. 329–344.

  8. Bemporad A. Spectroscopic detection of turbulence in post-CME current sheets // Astrophys. J. 2008. V. 689. Iss. 1. P. 572–584. doi 10.1086/592377

  9. Bian N.H., Kontar E.P. Stochastic acceleration by multi-island contraction during turbulent magnetic reconnection // Phys. Rev. Lett. 2013. V. 110. Iss. 15. id. 151101. doi 10.1103/PhysRevLett.110.151101

  10. Cartwright M.L., Moldwin M.B. Comparison of small-scale flux rope magnetic properties to large-scale magnetic clouds: Evidence for reconnection across the HCS? // J. Geophys. Res.: Space Phys. 2008. V. 113. Iss. A9. Cite ID A09105. doi 10.1029/2008JA013389

  11. Dai L., Wygant J.R., Cattel C.A., Thaller S., Kersten K., Breneman A., Tang X. Cluster observations of fast magnetosonic waves in the heliospheric current sheet // Geophys. Res. Lett. 2014. V. 41. P. 1398–1405. doi 10.1002/2014GL059223

  12. Drake J.F., Swisdak M., Che H., Shay M.A. Electron acceleration from contracting magnetic islands during reconnection // Nature. 2006. V. 433. Iss. 7111. P. 553–556. doi 10.1038/nature05116

  13. Eastwood J.P., Balogh A., Dunlop M.W., Smith C.W. Cluster observations of the heliospheric current sheet and an associated magnetic flux rope and comparisons with ACE // J. Geophys. Res.: Space Phys. 2002. V. 107. Iss. A11. P. SSH 9-1–SSH 9-9. doi 10.1029/2001JA009158

  14. Greco A., Perri S., Servidio S., Yordanova E., Veltri P. The complex structure of magnetic field discontinuities in the turbulent solar wind // Astrophys. J. Lett. 2016. V. 823. L39. Iss. 2. CiteID L39. doi 10.3847/2041-8205/823/2/L39

  15. Khabarova O., Zank G.P., Li G., le Roux J.A., Webb G.M., Dosch A., Malandraki O.E. Small-scale magnetic islands in the solar wind and their role in particle acceleration. I. Dynamics of magnetic islands near the heliospheric current sheet // Astrophys. J. 2015a. V. 808. Iss. 2. Cite ID 181. doi 10.1088/0004-637X/808/2/181

  16. Khabarova O., Zank G.P., Li G., le Roux J.A., Webb G.M., Malandraki O.E., Zharkova V.V. Dynamical small-scale magnetic islands as a source of local acceleration of particles in the solar wind // J. Phys. Conf. Ser. 2015b. V. 642. Cite ID 012033. doi 10.1088/1742-6596/642/1/012033

  17. Khabarova O., Zank G.P., Li G., Malandraki O.E., le Roux J.A., Webb G.M., Dosch A. Small-scale magnetic islands in the solar wind and their role in particle acceleration. II. Particle energization inside magnetically confined cavities // Astrophys. J. 2016. V. 827. Iss. 2. Cite ID 122. doi 10.3847/0004-637X/827/2/122

  18. Khabarova O.V., Zank G.P., Malandraki O.E., Li G., le Roux J.A., Webb G.M. Observational evidence for local particle acceleration associated with magnetically confined magnetic islands in the heliosphere—a review// Sun and Geosphere. 2017. V. 12. Iss. 1. P. 23–30.

  19. Malova H.V., Mingalev O.V., Grigorenko E.E., Mingalev I.V., Melnik M.N., Popov V.Yu., Delcourt D.C., Petrukovich A.A., Shen C., Rong D., Zelenyi L.M. Formation of self-organized shear structures in thin current sheets // J. Geophys. Res. Space Phys. 2015. V. 120. doi 10.1002/2014JA020974

  20. Markidis S., Henri P., Lapenta G., Divin A., Goldman M., Newman D., Laure E. Kinetic simulations of plasmoid chain dynamics // Phys. Plasmas. 2013. V. 20. Iss. 8. Cite ID 082105. http://dx.doi.org/ doi 10.1063/1.4817286

  21. Matthaeus W.H., Ambrosiano J.J., Goldstein M.L. Particle acceleration by turbulent magnetohydrodynamic reconnection // Phys. Rev. Lett. 1984. V. 53. P. 1449–1452. doi 10.1103/PhysRevLett.53.1449

  22. Merkin V.G., Lyon J.G., McGregor S.L., Pahud D.M. Disruption of a heliospheric current sheet fold // Geophys. Res. Lett. 2011. V. 38. id. L14107. doi 10.1029/2011GL047822

  23. Musielak Z.E., Suess S.T. Magnetohydrodynamic bending waves in a current sheet // Sol. Phys. 1988. V. 330. P. 456–465. doi 10.1086/166483

  24. Oka M., Phan T.-D., Krucker S., Fujimoto M., Shinohara I. Electron acceleration by multi-island coalescence // Astrophys. J. 2010. V. 714. Iss. 1. P. 915–926. doi 10.1088/0004-637X/714/1/915

  25. le Roux J.A., Zank G.P., Webb G.M., Khabarova O.A. Energetic ion acceleration by small-scale solar wind flux ropes // Astrophys. J. 2015. V. 801. Iss. 2. id. 112. doi 10.1088/0004-637X/801/2/112

  26. le Roux J.A., Zank G.P., Webb G.M., Khabarova O.V. Kinetic transport theory for particle acceleration and transport in regions of multiple contracting and reconnecting inertial-scale flux ropes // Astrophys. J. 2016. V. 827. Iss. 1. id. 47. doi 10.3847/0004-637X/827/1/47

  27. Ruderman M.S. Nonlinear surface wave propagation on heliospheric current sheet // COSPAR Colloq. Ser. 1990. V. 1. P. 249–252. doi 10.1016/B978-0-08-040780-7.50046-4

  28. Ruderman M.S. On the analogy between a system of harmonic oscillators and resonant Alfvén waves in plasmas // Phys. Plasmas. 1998. V. 5. Iss. 6. P. 2463–2465. doi 10.1063/1.872927

  29. Wang S., Lee L.C., Wei C.Q., Akasofu S.-I. A mechanism for the formation of plasmoids and kink waves in the heliospheric current sheet // Sol. Phys. 1988. V. 117. P. 157–169.

  30. Yamauchi M., Lui A.T.Y. Modified magnetohydrodynamic waves in a current sheet in space // Phys. Plasmas. 1997. V. 4. Iss. 12. P. 4382–4387. doi 10.1063/1.872600

  31. Zank G.P., le Roux J.A., Webb G.M., Dosch A., Khabarova O. Particle acceleration via reconnection processes in the supersonic solar wind // Astrophys. J. 2014. V. 797. Iss. 1. id. doi 10.1088/0004-637X/797/1/28

  32. Zank G.P., Hunana P., Mostafavi P., le Roux J.A., Li Gang, Webb G.M., Khabarova O.V. Particle acceleration by combined diffusive shock acceleration and downstream multiple magnetic island acceleration // J. Phys. Conf. Ser. 2015a. V. 642. Iss. 1. id 012031. doi 10.1088/1742-6596/642/1/012031, http://iopscience.iop.org/1742-6596/642/1/012031.

  33. Zank G.P., Hunana P., Mostafavi P., le Roux J.A., Li Gang, Webb G.M., Khabarova O., Cummings A., Stone E., Decker R. Diffusive shock acceleration and reconnection acceleration processes // Astrophys. J. 2015b. V. 814. Iss. 2. id 137. doi 10.1088/0004-637X/814/2/137

  34. Zharkova V., Khabarova O. Particle acceleration in the reconnecting heliospheric current sheet: solar wind data versus 3D PIC simulations // Astrophys. J. 2012. V. 752. Iss. 1. ID. 35. doi 10.1088/0004-637X/752/1/35

  35. Zharkova V., Khabarova O. Additional acceleration of solar-wind particles in current sheets of the heliosphere // Ann. Geophys. 2015. V. 33. Iss. 457. P. 457–470. doi 10.5194/angeo-33-457-2015

  36. Zheng J., Hu Q. Observational evidence for self-generation of small-scale magnetic flux ropes from intermittent solar wind turbulence // Astrophys. J. Lett. 2018. V. 852. Iss. 2. id. L23. doi 10.3847/2041-8213/aaa3d7

  37. Zhou X., Büchner J., Barta M., Gan W., Liu S. Electron acceleration by cascading reconnection in the solar corona. I. Magnetic gradient and curvature drift effects // Astrophys. J. 2015. V. 815. Iss. 1. id. 6. doi 10.1088/0004-637X/815/1/6

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