Журнал вычислительной математики и математической физики, 2023, T. 63, № 11, стр. 1894-1910

2.1. Исходные уравнения в частных производных и метод их решения

Приведем одномерные уравнения ЭМГД в форме законов сохранения (см. [12, 9]):

${d \mathord{\left/ {\vphantom {d {dt}}} \right. \kern-0em} {dt}} = {\partial \mathord{\left/ {\vphantom {\partial {\partial t}}} \right. \kern-0em} {\partial t}} + {\mathbf{v}} \cdot \nabla $ – полная производная, $n$ – плотность частиц (объемная концентрация) ионов или электронов, ${\mathbf{B}} = ({{B}_{x}},{{B}_{y}},{{B}_{z}})$ – напряженность магнитного поля, и ${\mathbf{v}} = (u,v,w)$ – скорость ионов. Используются безразмерные величины, определяемые через физические так: $x = {{\hat {x}} \mathord{\left/ {\vphantom {{\hat {x}} L}} \right. \kern-0em} L}$, $t = \hat {t}{{\omega }_{0}}$, ${\mathbf{v}} = {{{\mathbf{\hat {v}}}} \mathord{\left/ {\vphantom {{{\mathbf{\hat {v}}}} {{{V}_{A}}}}} \right. \kern-0em} {{{V}_{A}}}}$, ${\mathbf{B}} = {{{\mathbf{\hat {B}}}} \mathord{\left/ {\vphantom {{{\mathbf{\hat {B}}}} {\left| {{{{\mathbf{B}}}_{0}}} \right|}}} \right. \kern-0em} {\left| {{{{\mathbf{B}}}_{0}}} \right|}}$, ${\mathbf{B}} = {{{\mathbf{\hat {B}}}} \mathord{\left/ {\vphantom {{{\mathbf{\hat {B}}}} {\left| {{{{\mathbf{B}}}_{0}}} \right|}}} \right. \kern-0em} {\left| {{{{\mathbf{B}}}_{0}}} \right|}}$, $n = {{\hat {n}} \mathord{\left/ {\vphantom {{\hat {n}} {{{n}_{0}}}}} \right. \kern-0em} {{{n}_{0}}}}$; $L$, ${{\omega }_{0}} = {{{{V}_{A}}} \mathord{\left/ {\vphantom {{{{V}_{A}}} L}} \right. \kern-0em} L}$, ${{n}_{0}}$, ${{\left| {\mathbf{B}} \right|}_{0}}$, ${{V}_{A}} = \left| {{{{\mathbf{B}}}_{0}}} \right|{{[4\pi {{n}_{0}}({{m}_{e}} + {{m}_{i}})]}^{{ - 1/2}}}$ – характерная длина, частота, плотность невозмущенной плазмы, модуль вектора невозмущенного магнитного поля, альвеновская скорость. Здесь ${{m}_{i}}$ и ${{m}_{e}}$ – массы ионов и электронов, соответственно, ${{\nu }_{m}}$ – коэффициент магнитной вязкости, связанной с ионно-электронным трением, ${{b}^{2}}$ – коэффициент сжимаемости электронного газа. Параметры дисперсии ${{R}_{i}}$ и ${{R}_{e}}$ даются формулами: ${{R}_{i}} = {{{{\omega }_{{ic}}}} \mathord{\left/ {\vphantom {{{{\omega }_{{ic}}}} {{{\omega }_{0}}}}} \right. \kern-0em} {{{\omega }_{0}}}}$ и ${{R}_{e}} = {{{{\omega }_{{ec}}}} \mathord{\left/ {\vphantom {{{{\omega }_{{ec}}}} {{{\omega }_{0}}}}} \right. \kern-0em} {{{\omega }_{0}}}}$, где ${{\omega }_{{ic}}} = {{e\left| {{{{\mathbf{B}}}_{0}}} \right|} \mathord{\left/ {\vphantom {{e\left| {{{{\mathbf{B}}}_{0}}} \right|} {({{m}_{i}}c)}}} \right. \kern-0em} {({{m}_{i}}c)}}$ и ${{\omega }_{{ec}}} = {{e\left| {{{{\mathbf{B}}}_{0}}} \right|} \mathord{\left/ {\vphantom {{e\left| {{{{\mathbf{B}}}_{0}}} \right|} {({{m}_{e}}c)}}} \right. \kern-0em} {({{m}_{e}}c)}}$ – ионная и электронная циклотронные частоты соответственно, $e$ – заряд электрона, $c$ – скорость света. Для удобства расчета мы можем взять $L$ таким, что ${{\omega }_{0}} = \sqrt {{{\omega }_{{ic}}}{{\omega }_{{ec}}}} $, тогда ${{R}_{i}} = R_{e}^{{ - 1}} = \sqrt {{{{{m}_{e}}} \mathord{\left/ {\vphantom {{{{m}_{e}}} {{{m}_{i}}}}} \right. \kern-0em} {{{m}_{i}}}}} $. В расчетах далее подразумевается водородная плазма, $R_{e}^{{ - 1}} = 0.02341352$. Для правильной постановки граничных условий при расчетах методом конечных разностей система преобразуется к виду, в котором временные производные в каждом уравнении есть только от одного неизвестного (см. [9]):

Для расчета данных уравнений была разработана численная схема с центральными разностями по пространству и аппроксимацией временных производных с помощью метода Рунге–Кутты четвертого порядка [9]. В работе [15] для уравнений вида ${{a}_{t}} + c{{a}^{{(m)}}} = 0$, $m$ – порядок производной по $x$, показано, что схема на основе метода Рунге–Кутты для временных производных и аппроксимации по пространству центральными разностями условно устойчива для методов 1–4 порядков, но оптимально использовать методы 3 и 4 порядков, поскольку у них корректная схемная вязкость, что позволяет вести расчет неограниченное время при всех значениях $m$ и использовать естественное условие устойчивости $\tau < C{{h}^{m}}$, $\tau $ и $h$ – шаги по времени и пространству. Причем для метода третьего порядка схемная диссипация имеет третий порядок малости, а для метода четвертого порядка – пятый. Поэтому этот метод был рекомендован для широкого класса уравнений. Здесь использовался один из вариантов метода четвертого порядка. Поясним, как это делалось в общем виде.

Пусть численная схема описывается системой уравнений ${{{\mathbf{f}}}_{t}} = {\mathbf{F}}(t,{\mathbf{f}})$, где ${\mathbf{f}}$ – вектор, возникший в результате дискретизации по пространству. Здесь система автономная, аргумент $t$ включен для общности и понимания смысла шагов. Согласно используемому методу Рунге–Кутты аппроксимацию по времени можно представить формулами:

где $k$ – временной индекс. Однако в данном случае некоторые из уравнений имеют вид

Пространственно-временную аппроксимацию делаем так

где $l$ – пространственный индекс. Под звездочкой для $f$ понимаются те же индексы 1, 2, 3, 4, а для $n$, ${\mathbf{f}}$ и $t$ – те же значения аргументов функции ${\mathbf{F}}$, что и в выражениях (2.3). На каждом шаге метода Рунге–Кутты расчет последовательно выполняется для всех уравнений (2.2), поэтому величины ${{n}_{*}}$ на каждом шаге вычисляются раньше, чем ${{f}_{*}}$ в (2.4). Аналогично поступаем c уравнениями вида

Это приводит к неявным алгебраическим системам уравнений, которые решаются методом прогонки. Аналогично можно включать неявную аппроксимацию и для других членов, например вязкости, для улучшения условия устойчивости. Такой подход ранее успешно использовался как в классическом варианте (труба с газом), так и при наличии пространственно-временных производных (труба с несжимаемой жидкостью) [15]. Система (2.2) получена подстановкой уравнений простого неконсервативного вида, приведенных в [7], в консервативные уравнения (2.1). При программировании эта подстановка заново была сделана с заменой производных разностными операторами, что существенно упростило написание программы и исключило ошибки.

Схема в основном явная, неявные аппроксимации применялись только для пространственно-временных производных. Численный эксперимент показал, что для устойчивости и сходимости в недиссипативном случае следует придерживаться условия $\tau < Ch$.

Программа ориентирована на решение задач распространения уединенных волн и задачу о распаде произвольного разрыва. На границах, удаленных от исследуемой области, в программе ставятся условия фиксации заданных значений всех неизвестных, одних и тех же для трех приграничных ячеек. При анализе аппроксимации это эквивалентно заданию краевых значений, заданию нулевых значений первых производных по $x$ для $u$ и $n$, а также нулевых значений первых и вторых производных для остальных неизвестных. Это нужно для используемой численной схемы. Избыточность числа граничных условий [9] не создает проблем при решении таких задач. Программа оборудована средствами, вырезать фрагмент, добавить фрагмент. Это позволяет устранять рост возмущений вблизи границы, последствия отражения волн от границ, а также при распаде разрыва позволяет проводить расчет отдельно взятой структуры.

2.2. Дисперсионное соотношение и прогнозируемые типы волн

Приведем дисперсионное соотношение линеаризованного варианта системы (2.1) в недиссипативном случае при $n = 1$, $\left| {\mathbf{B}} \right| = 1$, а также скорости длинных быстрых и медленных магнитозвуковых длинных волн, а также альвеновских волн [7]:

Здесь $\theta $ – угол между вектором напряженности магнитного поля и осью $Ox$, $\omega $ – циклическая частота, $k$ – волновое число. Есть быстрая и медленная магнитозвуковые и альвеновская ветви дисперсионного соотношения $\omega = \omega (k)$. Скорости распространения быстрых и медленных магнитозвуковых, а также альвеновских волн для данных уравнений конечны при любых значениях длины волны ($k \in R$). При $k \to 0$ скорости распространения волн совпадают с соответствующими скоростями в МГД, при $k \to \infty $ скорости быстрых магнитозвуковых волн стремятся к $u \pm b$, а скорости медленных и альвеновских стремятся к $u$.

У быстрой магнитозвуковой ветви есть точка перегиба при $k > 0$ в случае $\theta < {{\theta }_{c}}$ (см. [7]), где $\theta = {\text{arctg}}\left( {{{{{{\left( {B_{y}^{2} + B_{z}^{2}} \right)}}^{{1/2}}}} \mathord{\left/ {\vphantom {{{{{\left( {B_{y}^{2} + B_{z}^{2}} \right)}}^{{1/2}}}} {{{B}_{x}}}}} \right. \kern-0em} {{{B}_{x}}}}} \right)$, ${{\theta }_{c}}$ – критическое значение. Прямая, касательная к этой ветви в начале координат, пересекает ветвь при $k > 0$. Есть также прямая ${{V}^{{rm}}} = {\omega \mathord{\left/ {\vphantom {\omega k}} \right. \kern-0em} k}$, касательная дисперсионной ветви при $k > 0$. Поэтому при ${{V}^{{rm}}} > V > {{V}^{ + }}$ соответствующая прямая два раза пересекает ветвь при $k > 0$. Далее под понятием пересечение будем иметь в виду пересечение при $k > 0$. Совпадение скоростей волн при наличии рационального отношения длин волн будем называть резонансом. Логично думать, что и в случае нелинейных уравнений будет наблюдаться аналогичный резонанс и его можно предсказать. В разд. 3 будет показано, что встречаются еще резонансы, связанные с нелинейными свойствами уравнений, не определяемые из дисперсионного соотношения. При $\theta > {{\theta }_{c}}$ точки перегиба и пересечения при $V > {{V}^{ + }}$ и $k > 0$ нет. В случае медленных магнитозвуковых волн при $V < {{V}^{ - }}$ имеются два пересечения: с медленной магнитозвуковой и с альвеновской ветвью. В [7] методом исследования нормальных форм в случае $\theta > {{\theta }_{c}}$ доказано существование солитонных решений с малой амплитудой волн. При $\theta < {{\theta }_{c}}$, а также для медленных магнитозвуковых волн согласно теории возможны обобщенные уединенные волны, т.е. комбинации уединенной и периодической волн. Во всех случаях магнитозвуковые волны малой амплитуды могут быть описаны уравнением Кортевега–де Вриза, при этом в случае обобщенной уединенной волны минимальная амплитуда периодической составляющей экспоненциально мала [16]. Для быстрых магнитозвуковых волн доказывалось также существование так называемых 1 : 1 уединенных волн (волны с осцилляциями, затухающими на бесконечности), наличие которых связано с существованием прямой, проходящей через начало координат и касающейся дисперсионной ветви при $k > 0$, у этих волн $V > {{V}^{{rm}}}$. В [10]– [12] были описаны численные методики анализа решений уравнений бегущих волн, корректные для нахождения симметричных обычных и 1 : 1 уединенных волн в случае, когда прямая, соответствующая скорости волны, не пересекает дисперсионные ветви при $k > 0$ и для нахождения обобщенных уединенных волн в случае, когда имеется пересечение.

Для альвеновских волн колебательного типа в зависимости от значений $\theta $ и $b$ возможно пересечение касательной в начале координат с быстрой магнитозвуковой ветвью при $\theta < {{\theta }_{{cam}}}$, соответственно есть значение ${{V}^{{ra}}}$, аналогичное значению ${{V}^{{rm}}}$, отсутствие пересечения при ${{\theta }_{{cam}}} < \theta < {{\theta }_{{caa}}}$, пересечение с альвеновской ветвью при $\theta > {{\theta }_{{caa}}}$ (см. [11]):

Данные формулы применимы в случае, когда аргумент арккосинуса меньше 1. В противном случае значение величины критического угла формально считать нулем.

В случае колебательных альвеновских волн численные расчеты задачи о распаде разрыва уравнений ЭМГД показывают, что нет обычных или обобщенных уединенных волн, в том числе и приближенных, что связано с особенностями нелинейности [9].

Можно при $\theta = 0$ рассматривать еще вращательные альвеновские волны. В [8] для них приведены точные решения, описывающие периодические волны произвольной амплитуды, причем их скорость от этой амплитуды не зависит, что позволяет предположить, что для этих волн нет и уединенной волны огибающей.

2.3. Метод исследования бегущих волн

Приведем уравнения волн, бегущих со скоростью $V$:

Индексом $0$ обозначены величины для состояния равновесия. При получении численных решений (2.5) решалось квадратное уравнение (2.6) для определения $n$, а затем находилось $u$. Знак при корне определялся в результате проверки того, что начальные значения в решении соответствуют заданным. Фактически нужно брать отрицательный знак при $V > b$ и положительный при $V < b$ [7], ниже при $\theta $ вблизи ${\pi \mathord{\left/ {\vphantom {\pi 2}} \right. \kern-0em} 2}$ использовался отрицательный знак, а вблизи $0$ – положительный.

Согласно методике поиска периодических волн, разработанной в [13] и примененной в [11] к уравнениям (2.5) без вязкости, находим состояние равновесия, выбираем среди начальных данных симметризуемую величину, варьируемую с некоторым шагом (внешний цикл программы), например ${{B}_{y}} = {{B}_{{y0}}} + Z$. Варьируем с некоторым шагом другую симметризуемую величину (внутренний цикл), например $v = {{v}_{0}} + Y$. Остальные начальные данные соответствуют состоянию равновесия. Рассчитывая (2.5), выбираем ${{x}_{{*p}}}$ такое, что антисимметризуемая величина, например $w$, обращается в ноль $p$ раз. Выбираем значение $Y$ так, чтобы вторая антисимметризуемая величина ${{B}_{z}}$ приняла значение ноль. Симметризуя/антисимметризуя решение относительно ${{x}_{{*p}}}$, находим периодическое решение на периоде. Симметричные и антисимметричные переменные могут быть определены и противоположным способом, наоборот, но при исследовании предельных солитонных решений уравнений (2.5) предпочтительнее первый вариант. В результате работы программы для каждого $p$ из некоторого заданного диапазона значений создается файл со столбцами данных $(Z,Y)$. Для расчета применялся метод Рунге–Кутты второго порядка.

Делаем поточечные отображения данных файлов на плоскости $(Z,Y)$. Для разных значений $p$ можно использовать разные цвета или разную жирность точки. В результате на этой плоскости образуются ветви. Анализируем ветви различными способами.

Точка, соответствующая одному и тому же решению, может отображаться при разных значениях $p$, поэтому ищем для него минимальное значение $p$, обозначим его через $q$. Для этого вначале отображаются точки с большим значением $p$, а затем – с меньшим, точки с большим значением $p$ будут перекрыты точками с меньшим значением. По величине $q$ можно судить об отношении периодов волн двоякопериодического решения. Резонансная двухволновая ветвь с соотношением длин волн $n$ может соответствовать $q = n$, одноволновые ветви – $q = 1$, хотя это не всегда так. В частности, ниже при исследовании резонансов коротких волн с длинными альвеновскими и магнитозвуковыми волнами максимальное $p$ бралось равным 10 и этого было достаточно, хотя рассматривались резонансы с большими значениями отношений длин волн. Для одноволновых, коротковолновых альвеновских и быстрых магнитозвуковых ветвей $q$ действительно было равным единице.

Строятся графики решений для различных неизвестных и соотношение длин волн определяется визуально. Проводится анализ расположения ветвей линеаризованных решений с тем, чтобы выявить закономерности и продолжить их на нелинейные. В случае однократного пересечения с дисперсионными ветвями прямой, соответствующей скорости $V$, через точку равновесия проходит одноволновая ветвь линейных волн. В случае двукратного пересечения – две одноволновые ветви. Анализируется порядок расположения ветвей. Одноволновые ветви являются первичными. В случае устойчивости участка ветви он пересекается вторичными двухволновыми резонансными ветвями. Вторичные ветви в случае устойчивости могут пересекаться ветвями фазовых колебаний двух волн. Фазовые колебания резонансных решений неинтегрируемых по Лиувиллю систем рассматривались в [17]. Делается также проверка устойчивости ветви путем численного эксперимента. Результатом развития неустойчивости периодических решений вторичной ветви может быть приближенное решение в виде уединенной волны, огибающей решения с фазовыми колебаниями [12, 11], первичной коротковолновой ветви – в виде обобщенной уединенной волны. Однако при большом отношении длин волн провести такой детальный анализ затруднительно, ниже применялся иной подход, основанный на общем анализе расположения точек с выявлением точек равновесия (2.5) и коротковолновых, одноволновых ветвей.

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

Анализ числа пересечений прямой $V = {\omega \mathord{\left/ {\vphantom {\omega k}} \right. \kern-0em} k}$ с дисперсионными ветвями в данном методе не обязателен, но полезен, при двух пересечениях ожидаются двоякопериодические решения, при одном пересечении – однопериодные (двоякопериодические возможны), но при отсутствии пересечений те и другие тоже могут быть.

Подробное описание результатов исследования ветвей периодических решений для быстрых магнитозвуковых волн при наличии двухволнового резонанса дано в [11].

3.1. О существовании обобщенных уединенных медленных магнитозвукоых волн

Согласно теории [7] для медленных магнитозвуковых волн можно допустить наличие обобщенных уединенных волн и структур разрывов типа переходов между периодическим и двоякопериодическим состоянием. Но при расчете задачи о распаде разрыва такие структуры не должны наблюдаться. Для медленных магнитозвуковых волн прогнозируется стохастическое решение с излучением коротких альвеновских волн внутрь волновой зоны [12]. В расчетах уравнений ЭМГД [9] при малых амплитудах скачка наблюдалась структура солитонного типа. При увеличении амплитуды возникало излучение внутрь волновой зоны, но оно расценивалось как схемный эффект, поскольку обеспечить достаточно мелкий пространственный шаг для аппроксимации коротких альвеновских волн не представлялось возможным. Метод, основанный на исследовании (2.5) посредством анализа собственных значений и линейной аппроксимации периодической компоненты, показывал отсутствие обобщенных уединенных волн [9].

Было осуществлено исследование методом анализа периодических решений для начального состояния равновесия ${{n}_{0}} = 1$, ${{u}_{0}} = 0$, ${{w}_{{z0}}} = 0$, ${{B}_{x}} = \cos \theta $, ${{B}_{{y0}}} = \sin \theta $, ${{B}_{{z0}}} = 0$ назовем это состоянием 1; скорость волн $V = {{V}^{ - }}(1 + \mu )$. В отличие от случая быстрых магнитозвуковых волн, где можно выбрать значения $\mu $, при которых длины коротких и длинных волн резонансных волн сопоставимы и можно выделить резонансные ветви, в случае резонанса медленных магнитозвуковых и альвеновских волн отношение периодов велико при всех значениях $\mu $ (порядка 1000). При этом ветви сливаются и точки, отображающие их на графиках, как бы группируются, создавая подобия областей, хотя с математической точки зрения это множества меры ноль, выделяются только коротковолновые одноволновые альвеновские ветви. Но при исследовании чисто альвеновских ветвей для случая, когда касательная в начале координат пересекает альвеновскую ветвь при $k \ne 0$, как и для быстрых магнитозвуковых ветвей, можно подобрать скорость так, чтобы длины волн стали сопоставимы и выделить резонансные ветви. В процессе исследования проводились расчеты с постепенным увеличением значения $\mu > 0$. При таком выборе знака $\mu $ начальная точка равновесия должна быть ассоциирована с подошвой обобщенной медленной магнитозвуковой уединенной волны. Согласно проведенному исследованию “область” двоякопериодических ветвей, содержащая одноволновую ветвь медленных магнитозвуковых волн, практически доходит до точки равновесия при умеренных значениях $\mu $, поэтому можно получить решение, практически неотличимое от чисто солитонного. Детальные расчеты показывают, что эта ветвь все же отделена от точки равновесия и уединенной волны нет. Если допустить, что и при малых значениях $\mu $ медленная магнитозвуковая ветвь и резонансные ветви отделены от точки равновесия, то точного решения в виде обычной или обобщенной уединенной волны, связанной с резонансом медленных магнитозвуковых и альвеновских волн, при малых значениях $\mu $ нет. Если же при малых $\mu $ обобщенная уединенная волна существует, то минимальная амплитуда ее периодической компоненты настолько мала, что ее невозможно определить численно. Заметим, что аналогичные расчеты, см. ниже, проводились и для случая $\theta $, близкого к нулевому значению, для быстрых магнитозвуковых волн и также на основе численного анализа периодических решений однозначно нельзя сказать, существует или нет обобщенная уединенная волна, связанная с саморезонансом быстрых магнитозвуковых волн, хотя для $\theta $, близкого к ${{\theta }_{c}}$, ее существование убедительно подтверждалось [9].

3.2. Взаимодействие длинных альвеновских и длинных медленных магнитозвуковых волн

Особый интерес представляют исследования взаимодействия медленных магнитозвуковых и альвеновских волн в случае, когда величина $\theta $ близка к ${\pi \mathord{\left/ {\vphantom {\pi 2}} \right. \kern-0em} 2}$ (${{V}^{a}} = {{V}^{ - }} = 0$ при $\theta = {\pi \mathord{\left/ {\vphantom {\pi 2}} \right. \kern-0em} 2}$), примеры расположения ветвей показаны на фигурах ниже, $\theta = 1.525$, ${{b}^{2}} = 0.5$. При этом величины ${{V}^{ - }}$ и ${{V}^{a}}$ близки, оказалось, что возможны одновременные взаимодействия длинных медленных магнитозвуковых, длинных альвеновских волн и коротких альвеновских волн, т.е. трехволновой резонанс. Кроме того, с точки зрения физики при $\theta \approx {\pi \mathord{\left/ {\vphantom {\pi 2}} \right. \kern-0em} 2}$ затухание Ландау слабое. В этом случае имеется три точки равновесия системы (2.5). На фигурах ниже точки равновесия пронумерованы, проставлены обозначения $A$ или $M$, указывающие на то, что через них проходит длинноволновая альвеновская или медленная магнитозвуковая ветвь, $(M)$ указывает на то, что длинноволновая медленная магнитозвуковая ветвь подходит близко к этой точке, но не доходит до нее.

На фиг. 1а, $\mu = 0.5$, показано расположение медленных магнитозвуковых и альвеновских ветвей в случае, когда у медленных магнитозвуковых волн есть взаимодействие только с короткими альвеновскими волнами. Из точки 1 выходят ветви коротких альвеновских волн, но вблизи располагаются также резонансные ветви длинных магнитозвуковых и коротких альвеновских волн. Из точки 2 выходят ветви длинных магнитозвуковых волн и коротких альвеновских. Верхняя точка, точка 3, соответствует альвеновским волнам. Из нее выходят короткая и длинная альвеновские ветви. Есть длинноволновое двоякопериодическое солитоноподобное решение с преобладанием длинных альвеновских волн. Согласно МГД эволюционными являются переходы $1 \to 2$, $3 \to 1$. Под эволюционностью здесь понимается выполнение условия того, что число граничных условий на разрыве равно числу уходящих характеристик упрощенной системы, которой здесь является система МГД-уравнений, плюс единица. Это необходимое условие для корректности постановки задачи для МГД с разрывом. Здесь подразумевается, что разрывы в смысле эволюционности являются аналогами ударных волн в газовой динамике: две приходящие характеристики одного типа по разные стороны разрыва, а у характеристик остальных типов – с одной стороны уходящая, а с другой – приходящая, класс 0 согласно классификации по числу свободных параметров [14]. Причем для первого перехода есть приближенная система периодических медленных магнитозвуковых волн (фактически там присутствует еще компонента коротких альвеновских волн малой амплитуды), позволяющая в недиссипативном и слабодиссипативном случае построить приближенную солитонную структуру такого перехода с волновой зоной, описываемой усредненными уравнениями. При малых $\mu $ медленная магнитозвуковая уединенная волна на границе волновой зоны может быть описана уравнением Кортевега–де Вриза.

На фиг. 1б, $\mu = 0.575$, показан случай, когда скорость приблизилась к скорости альвеновских волн. Есть солитоноподобные двоякопериодические решения, ассоциированные с медленными магнитозвуковыми волнами, и гибридные, ассоциированные с альвеновскими и медленными магнитозвуковыми волнами одновременно, фиг. 2, кривые 1 и 2 (для кривой 2 показан полупериод). Для гибридного солитоноподобного решения в отличие от медленного магнитозвукового решения асимптотического приближения не найдено, это существенно нелинейное решение. Помимо системы приближенных периодических волн, для перехода $1 \to 2$ есть и приближенная система двоякопериодических волн для перехода $3 \to 1$ (если говорить точно, то это троякопериодические решения). Среди них есть решения, содержащие одновременно оба типа описанных выше приближенных обобщенных уединенных волн, фиг. 2, кривая 3, в рамке в увеличенном виде детально показаны малоамплитудные альвеновские периодические короткие волны в окрестности $x = 0$ для этого решения. В реальности в масштабе, использованном для графиков на фиг. 2, размер рамки по оси $x$ должен быть равен 1, а по оси ${{B}_{y}}$ равен 0.01. Отметим, что во всех случаях присутствуют эти короткие альвеновские волны, не видимые на графиках 1, 2, 3. Заметим, что хотя приближенное медленное магнитозвуковое солитонное решение выявлено здесь при довольно большом значении $\mu $, при дальнейшем его незначительном увеличении ($\mu = 0.72$) амплитуда коротковолновой альвеновской составляющей становится уже существенной, а амплитуда солитоноподобной медленной магнитозвуковой компоненты снижается [11]. Когда скорость почти совпадает с альвеновской, но все же меньше ее ($\mu = 0.73$), по-прежнему эволюционны переходы $1 \to 2$, $3 \to 1$. Приближенные солитонные решения или классические кинки (структуры класса –1 [14] с дополнительным условием на разрыве) здесь не обнаруживаются, в окрестностях точек 1 и 3 выявляются только двоякопериодические и троякопериодические решения, комбинации из длинной магнитозвуковой и короткой альвеновской волны.

На фиг. 3а показан случай, когда скорость незначительно превышает альвеновскую, ${{B}_{y}} - {{B}_{{y0}}} < 0$, $v > 0$, $\mu = 0.8$. Точка 1 теперь ассоциирована с альвеновскими волнами и из нее выходят коротковолновая и длинноволновая альвеновские ветви. Эволюционны переходы $3 \to 2$ и $1 \to 3$. Также как и в случае $\mu = 0.575$ имеются два типа приближенных солитонных решений и для обоих переходов есть система приближенных периодических и двоякокопериодических волн, позволяющая осуществить этот переход в слабодиссипативном случае. По мере увеличения $\mu $ точки 2 и 3 сближаются, амплитуда медленного магнитозвукового разрыва становится сколь угодно малой, т.е. соответствующее приближенное солитонное решение становится асимптотически точным.

При дальнейшем увеличении скорости (фиг. 3б, $\mu = 0.9$) сохраняется только альвеновская точка равновесия и двухволновые альвеновские решения. Приближенных солитонных решений здесь нет, что соответствует результатам расчетов уравнений ЭМГД в [9].

Отметим, что при выходе коротковолновой альвеновской ветви на границу зоны существования периодических решений (это понятие визуальное, область заполненная точками, с математической точки зрения это, очевидно, множество меры ноль) возникает чисто вращательная альвеновская волна, аналитически описанная в [8]. Обнаружены и двоякопериодические вращательные альвеновские волны.

Некоторые из границ зон существования двоякопериодических решений (в частности граница в виде практически прямой наклонной линии) соответствуют заострению волн коротковолновой альвеновской компоненты. Заострения длинных медленных магнитозвуковых и длинных альвеновских волн при наличии резонанса обнаружено не было. Проводилось исследование расположения альвеновских ветвей и в случае, когда касательная к альвеновской дисперсионной ветви пересекает быструю магнитозвуковую ветвь. При $U < {{V}_{a}}$ были выявлены резонансные двоякопериодические решения, причем граница их существования также, как и в случае пересечения с альвеновской дисперсионной ветвью ассоциирована с заострением коротких волн, но в этом случае быстрых магнитозвуковых. Конец одноволновой быстрой магнитозвуковой ветви также ассоциирован с заострением. Можно предположить, что при наличии резонанса опрокидывание длинных волн компенсируется излучением коротких волн и в этом смысле этот эффект аналогичен действию вязкости. Короткие волны могут опрокинуться.

Уравнения (2.5) анализировались на наличие всех состояний равновесия в общем случае. Когда $\theta $ не близко к ${\pi \mathord{\left/ {\vphantom {\pi 2}} \right. \kern-0em} 2}$ при малых $\mu $ выявляется только два состояния равновесия, связанных с медленными магнитозвуковыми волнами, при увеличении $\mu $ число реальных, т.е. с $n > 0$, состояний равновесия может увеличиваться до четырех, можно допустить существование и более сложных нелинейных взаимодействий волн.

3.3. Применение более простой системы уравнений бегущих волн

В случае применения более простого варианта системы (2.1) с $R_{e}^{{ - 1}} = 0$ система уравнений бегущих волн выводится из (2.5) и имеет второй порядок, как и система бегущих волн для классического и модифицированного уравнения Кортевега–де Вриза, но это система с более сложной нелинейностью:

Значения ${{B}_{y}}$ и ${{B}_{z}}$ из нижних двух уравнений подставляются в верхние два. В случае расчета с магнитной вязкостью требуется еще решить два алгебраических уравнения для нахожения ${{d{v}} \mathord{\left/ {\vphantom {{d{v}} {dx}}} \right. \kern-0em} {dx}}$ и ${{dw} \mathord{\left/ {\vphantom {{dw} {dx}}} \right. \kern-0em} {dx}}$. При использовании газодинамической вязкости получается система третьго порядка, поскольку (2.5) дополняется уравнением для $u$ [9]. Это физическое, но не асимптотическое приближение. На фазовом портрете системы уравнений бегущих волн в недиссипативном случае положения равновесия, отмеченные на фиг. 1 и 3 буквами $A$ и $M$ становятся эллиптическими точками, а $(M)$ – гиперболическими. Появляются солитонные решения. Здесь при $V \ne {{V}^{a}}$ имеются точные медленные магнитозвуковые и гибридные альвено-магнитозвуковые уединенные волны. Напомним,что рассмотрение полной системы показывает, что не всегда рассматриваемые уединенные волны существуют даже в приближенном виде. При добавлении слабой магнитной вязкости можно построить структуры слабодиссипативных разрывов. Некоторые из них показаны на фиг. 4а – для $\mu = 0.575$, $V < {{V}^{a}}$, ${{\nu }_{m}} = 0.01$, б – для $\mu = 0.73$, $V < {{V}^{a}}$, ${{\nu }_{m}} = 0.001$, в – для $\mu = 0.8$, $V > {{V}^{a}}$, ${{\nu }_{m}} = 0.01$. При этом в случае, соответствующем фиг. 4а, можно найти солитонную структуру медленного магнитозвукового разрыва $1 \to 2$ при любых значениях ${{\nu }_{m}}$, и приближенную структуру альвеновского разрыва $3 \to 1$ при ${{\nu }_{m}} \to 0$, показанную на фиг. 4а. Участок альвеновской структуры условно показан сплошной линией, но колебания вблизи $x = 0$, показанные пунктирной линией, сохранены для демонстрации того, что это решение приближенное (фактически расчет проводился из окрестности точки равновесия 3 в направлении, обратном оси $x$). В этом решении колебательная зона двухволновая, а при приближении к точке равновесия 1 наблюдаются приближенные уединенные волны двух разных типов, т.е. возможны две структуры с разными уединенными волнами на границе. Чем меньше значение ${{\nu }_{m}}$, тем точнее волны на границе приближаются к уединенным волнам. При приближении $V$ к значению ${{V}^{a}}$ можно найти еще гибридную альвено-медленномагнитозвуковую структуру $3 \to 2$ при малых значениях ${{\nu }_{m}}$, в последнем решении имеется внутренний свободный параметр, на рис. 4б приведено именно это решение. Оно содержит внутри себя приближенные структуры первых двух типов, описанные выше. Стоит обратить внимание на то, что здесь амплитуды приближенных магнитозвуковых и гибридных уединенных волн уже близки. Здесь в обозначениях разрывов вида ${{n}_{1}} \to {{n}_{2}}$ условно принято, что разрыв движется вправо, первое состояние находится справа, второе – слева, по аналогии с описанием ударных волн в газовой динамике, движущихся вправо, где есть понятия состояние до разрыва и после разрыва.

Разрыв $3 \to 1$ оценивается как альвеновский, поскольку с обеих сторон имеются приходящие альвеновские характеристики. Разрыв $3 \to 2$ оценивается как гибридный, поскольку с обеих сторон приходят и альвеновские, и медленные магнитузвуковые характеристики. Заметим, что разрывы с двумя приходящими характеристиками с двух сторон ранее рассматривались для уравнений, описывающих двумерную эволюцию фронта волны на поверхности жидкости, имеющей в сечении по вертикали форму уединенной волны, а также аналогичные разрывы были обнаружены для системы уравнений, описывающей установившиеся периодические волны на поверхности жидкости [12]. Эти разрывы именовались пересечениями волн, фактически рассматривались двузначные и двухволновые решения по обе стороны разрыва, для них существуют эволюционные структуры разрывов. Указанные структуры описываются уравнением Кадомцева–Петвиашвили и нелинейным уравнением Шрёдингера. Для перехода между двумя периодическими состояниями внутри гибридной структуры просматривается также и некоторая аналогия с трехволновыми конфигурациями для тех же моделей. Действительно, справа одноволновая зона, а слева – двухволновая. Решение в виде гибридной стационарной структуры имеет один дополнительный параметр, т.е. это семейство, но это внутренний параметр решения, а не параметр МГД-уравнений по обе стороны разрыва, поэтому структура неэволюционна.

При $V > {{V}^{a}}$ точки 1 и 3 как бы меняются местами, графики решений получаются как бы центрально отраженными относительно точки 2. На фиг. 4в показана гибридная структура $1 \to 2$, содержащая приближенные альвеновскую структуру $1 \to 3$ и медленную дозвуковую $3 \to 2$. Фазовая траектория этого решения показана на фиг. 5, гиперболическая точка 3 отмечена крестом. При дальнейшем увеличении $V$ точки 2 и 3 сливаются и исчезают, пропадают и структуры разрывов.

Заметим, что в случае модифицированного уравнения Кортевега–де Виза с кубической нелинейностью также имеется три точки равновесия [12], при определенном выборе знаков коэффициентов уравнения, средняя точка эллиптическая, а две крайние – гиперболические, при этом при некотором значении $V$ возникает так называемый кинк – переход между крайними гиперболическими точками, структура разрыва с дополнительным условием, класс –1 (см. [14]). Может также быть средняя гиперболическая точка и крайние эллиптические, тогда есть только уединенные волны. В этом фазовый портрет рассматриваемой системы принципиально отличается, одна из крайних точек гиперболическая – остальные точки эллиптические. Но при $V = {{V}^{a}}$ тоже есть подобие особого решения, в каком-то смысле аналог кинка. При подходе к крайним точкам возникает периодическая волна, похожая на синусоиду, период волны стремится к бесконечности, а максимумы и минимумы соответствуют крайним точкам равновесия.

3.4. Расчет задачи о распаде разрыва

Проводились также расчеты исходных уравнений в частных производных с начальной сглаженной ступенькой с параметрами, соответствующими альвеновской и гибридной структуре. Применялась численная схема, разработанная в [9]. Результаты показаны на фиг. 6, $\mu = 0.8$, ${{\nu }_{m}} = 0.01$, $t = 50\,000$ (пунктиром), $t = 10\,000$ (сплошная). Вместо неэволюционной гибридной структуры возникает нестационарная расширяющаяся со временем структура, в которой имеются как альвеновские, так и медленные магнитозвуковые волны, фиг. 6а. Но некоторая аналогия с фиг. 4в все же есть, в частности заметно наличие двухволновой зоны. При перемене местами параметров для правой и левой сторон возникают медленная магнитозвуковая волна и альвеновская структура, фиг. 6б. Они трудно отделимы друг от друга возможно потому, что при наличии вязкости альвеновские структуры как бы расплываются [9]. Вместо эволюционной приближенной альвеновской структуры возникают две структуры: медленная магнитозвуковая типа простой волны и расширяющаяся со временем альвеновская структура максимальной амплитуды с поворотом поперечного магнитного поля на 180 градусов, фиг. 6в. Фактически здесь альвеновская структура – нестационарный аналог кинка, особого разрыва с дополнительным условием. Интересно, что если поменять параметры для правой и левой стороны местами, то аналогичная кинкоподобная альвеновская структура тоже возникает фиг. 6г, но вместо простой медленной магнитозвуковой волны здесь солитонная структура. Подобным свойством обладают решения о распаде разрыва для модифицированного уравнения Кортевега–де Вриза [12].

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

3.5. Об использовании различных систем уравнений для расчетов МГД и ЭМГД

Расчеты эволюции начальных данных в виде сглаженной ступеньки, соответствующей рассмотренным выше структурам разрывов, проводились для уравнений (2.1) и уравнений с $R_{e}^{{ - 1}} = 0$, как с магнитной вязкостью, так и без нее. Результаты расчетов были схожими, поскольку при использовании системы (2.1) короткие волны не выявлялись из-за недостаточно мелкого пространственного шага сетки, а влияние допущения $R_{e}^{{ - 1}} = 0$ на длинные волны незначительное из-за малости этой величины по сравнению с $R_{i}^{{ - 1}}$. Добавление небольшой вязкости в (2.1) позволяет делать и формально правильные расчеты полных уравнений, поскольку короткие волны должны подавляться вязкостью. В то же время очевидно, что допущение $R_{e}^{{ - 1}} = 0$ существенно влияет в случае больших градиентов и большой амплитуды разрывов, поскольку опрокидывание волн может устраняться за счет излучения коротких волн. Это допущение совсем не применимо к быстрым магнитозвуковым волнам при ${\pi \mathord{\left/ {\vphantom {\pi 2}} \right. \kern-0em} 2} \geqslant \theta \geqslant {{\theta }_{c}}$, поскольку оно приводит для них к прямо противоположному типу дисперсии. Подчеркнем, что расчеты полной системы, системы с $R_{e}^{{ - 1}} = 0$, и классических уравнений МГД проводились одной программой, во всех случаях для устойчивости схемы в недиссипативном случае временной шаг брался пропорциональным пространственному, поэтому принципиальной разницы в расходе машинного времени для этих расчетов нет. Хотя, например, для расчетов с $R_{e}^{{ - 1}} = 0$, экспериментально выявленный коэффициент пропорциональности пришлось брать меньшим, чем в случае $R_{e}^{{ - 1}} \ne 0$. Но применение системы с $R_{e}^{{ - 1}} = 0$ полезно для корректного анализа свойств длинных медленных магнитозвуковых и альвеновских волн. В случае недиссипативных уравнений МГД используемая численная схема сама может устранить опрокидывание и создать схемные волновые зоны, но более естественно включать слабую физическую или искусственную вязкость, чтобы получить классические структуры разрывов [9]. В случае наличия дисперсионных членов в системе уравнений включение вязкости также может потребоваться для устранения опрокидывания разрывов большой интенсивности.

3.6. Взаимодействие длинных альвеновских и длинных быстрых магнитозвуковых волн

В связи с выявленным нелинейным резонансом альвеновских и медленных магнитозвуковых волн аналогичные исследования были проведены и для взаимодействия альвеновских и быстрых магнитозвуковых волн при $\theta $, близком к нулю (${{V}^{a}} = {{V}^{ + }} = 1$ при $\theta = 0$). Затухание Ландау при этом, как и в случае $\theta \approx {\pi \mathord{\left/ {\vphantom {\pi 2}} \right. \kern-0em} 2}$, слабое. Примеры расположения ветвей при $\theta = 0.3$ приведены на фиг. 7, a – для $\mu = - 0.0025$, б – для $\mu = 0.01$. Здесь $V = {{V}^{a}}(1 + \mu )$. Картины расположения ветвей аналогичны случаю взаимодействия альвеновских и медленных магнитозвуковых волн. Фиг. 7а и фиг. 7б, в аналогичны фиг. 3а и фиг. 1б. Соответственно здесь также наблюдается резонанс длинных альвеновских и длинных быстрых магнитозвуковых волн и существуют структуры разрывов, аналогичные выявленным при длинноволновом нелинейном резонансе альвеновских и быстрых магнитозвуковых волн. При уменьшении $\mu $ точки 2 и 3 сливаются и исчезают, расположение ветвей сходно с фиг. 1а. Остаются только длинные альвеновские волны и короткие магнитозвуковые волны. При дальнейшем уменьшении исчезают альвеновские волны, а затем и быстрые магнитозвуковые. При увеличении $\mu $ альвеновская точка 3, как и в случае фиг. 3б, удаляется от точек 1 и 2, взаимодействия между быстрыми магнитозвуковыми и альвеновскими волнами нет.