Содержание к диссертации
Введение
ГЛАВА I. ФИЗИЧЕСКАЯ И МАТЕМАТИЧЕСКАЯ ПОСТАНОВКА ЗАДАЧИ 16
1.1 Обоснование выбора магнитогазодинамической модели с учетом диффузии магнитного поля 16
1.2 Анализ пригодных для данного класса задач вычислительных моделей и обоснование выбора высокоточных монотонных разностных схем 22
ГЛАВА II. ПОСТРОЕНИЕ ВЫСОКОТОЧНЫХ МОНОТОННЫХ РАЗНОСТНЫХ СХЕМ ДЛЯ ГИПЕРБОЛИЧЕСКИХ СИСТЕМ УРАВНЕНИЙ 24
2.1 Монотонные разностные схемы высокого порядка аппроксимации для одномерных систем уравнений гиперболического типа 24
2.1.1 Монотонные схемы второго-третьего порядка аппроксимации для простейшего уравнения переноса 25
2.1.2 Обобщение монотонных схем второго-третьего порядка аппроксимации на случай одномерной гиперболической системы уравнений 33
2.2 Обобщение полученных высокоточных монотонных схем на случай трех пространственных переменных 37
ГЛАВА III. ОБОБЩЕНИЕ И ТЕСТИРОВАНИЕ ПОЛУЧЕННЫХ ВЫСОКОТОЧНЫХ ЧИСЛЕННЫХ МЕТОДОВ НА СЛУЧАЙ УРАВНЕНИЙ ГАЗОВОЙ ДИНАМИКИ 40
3.1 Построение монотонных разностных схем высокого порядка аппроксимации для одномерных систем уравнений газовой динамики 40
3.2 Применение энтропийной коррекции при численном решении уравнений газовой динамики 50
3.3 Построение монотонных разностных схем высокого порядка аппроксимации для многомерных систем уравнений газовой динамики 58
ГЛАВА IV. ОБОБЩЕНИЕ И ТЕСТИРОВАНИЕ ПОЛУЧЕННЫХ ВЫСОКОТОЧНЫХ ЧИСЛЕННЫХ МЕТОДОВ НА СЛУЧАЙ УРАВНЕНИЙ МАГНИТНОЙ ГИДРОДИНАМИКИ 66
4.1 Построение монотонных разностных схем высокого порядка аппроксимации для одномерных систем уравнений магнитной гидродинамики 66
4.2 Построение монотонных разностных схем высокого порядка аппроксимации для многомерных систем уравнений магнитной гидродинамики 77
ГЛАВА V. РЕЗУЛЬТАТЫ ЧИСЛЕННОГО ИССЛЕДОВАНИЯ ПОВЕДЕНИЯ ПЛАЗМЕННОГО ОБЛАКА В ВЕРХНЕЙ ИОНОСФЕРЕ 80
5.1 Алгоритм выделения компонент односкоростного плазменного течения 80
5.2 Результаты расчёта динамики плазмы умеренного начального энергосодержания и их анализ 82
5.3 Результаты расчёта динамики плазмы при больших энерговыделениях 93
ЗАКЛЮЧЕНИЕ 106
ЛИТЕРАТУРА 107
- Обоснование выбора магнитогазодинамической модели с учетом диффузии магнитного поля
- Монотонные разностные схемы высокого порядка аппроксимации для одномерных систем уравнений гиперболического типа
- Построение монотонных разностных схем высокого порядка аппроксимации для одномерных систем уравнений газовой динамики
- Построение монотонных разностных схем высокого порядка аппроксимации для одномерных систем уравнений магнитной гидродинамики
- Алгоритм выделения компонент односкоростного плазменного течения
Введение к работе
Актуальность темы
В последние десятилетия проведено большое количество экспериментов с выбросом в верхнюю атмосферу и магнитосферу плазменных и слабоионизованных сгустков и струй сравнительно высокой массы. К наиболее интересным из подобного рода геофизических экспериментов относятся эксперименты с использованием мощных возмущений. Научное и прикладное значение подобных экспериментов зависит, прежде всего, от того насколько удаётся понять и численно промоделировать их результаты. Как известно, эволюция крупномасштабных возмущений в верхней атмосфере и магнитосфере сопровождается весьма сложными и многообразными физическими процессами, по некоторым из которых пока нет достаточно полного понимания. Кроме того, эволюция образующейся возмущённой области, вообще говоря, имеет пространственный трёхмерный характер, и лишь на отдельных временных стадиях своего развития может приближенно моделироваться геометрией меньшей размерности. Эти два обстоятельства существенно осложняют процесс численного моделирования эволюции плазменных сгустков в ионосфере и магнитосфере, имеющих высокую начальную энергию и массу (эрг, г). Проделать это на основе имеющихся комплексов программ не представляется возможным, так как они обычно применимы лишь для слабоградиентных течений идеальных сплошных сред, и не учитывают определяющих внутренних процессов, связанных с неравновесной ионизацией и диффузией магнитного поля.
Цель работы
Целью работы является разработка высокоточных численных методов, их программная реализация для моделирования сильных возмущений в ионосфере Земли и последующей эволюции образующихся плазменных сгустков.
В данной работе были поставлены следующие задачи:
создание высокоточных численных методов и реализация их в виде программного комплекса для трехмерных уравнений магнитогазодинамики с учетом диффузии магнитного поля;
всестороннее тестирование программного комплекса и его апробация на прикладных задачах;
разработка удобного интерфейса для представления и визуализации результатов расчетов.
Помимо этого, естественной задачей данной работы является уточнение и развитие известных ранее соответствующих математических моделей и численных методов с целью повышения их эффективности.
Научная новизна
Для трехмерных уравнений магнитной газодинамики (МГД - модели с учетом диффузии магнитного поля) разработаны новые модификации монотонных консервативных вариантов сеточно-характеристического метода и метода типа Годунова 2-3-го порядка аппроксимации.
Общий алгоритм перехода от известного состояния на текущем временном слое к искомому состоянию в следующий момент времени включал расщепление по пространственным переменным, а, при наличии разрывов большой интенсивности, также по физическим процессам («газодинамический» этап, «магнитный» этап, этап расчета диффузии магнитного поля).
Выполнены численные исследования эволюции сильных возмущений плазмы в околоземном космическом пространстве до времен порядка 100 секунд. Показано существенное влияние неоднородности экспоненциальной атмосферы и геомагнитного поля на структуру течения плазмы.
Практическая ценность
Построен и реализован в виде комплекса программ эффективный монотонный численный алгоритм решения трёхмерных задач плазмодинамики, обладающий высоким порядком аппроксимации. Алгоритм применён для решения крупномасштабных МГД - задач о поведении плазмы в околоземном космическом пространстве при существенно различных условиях – на высотах от 150 до 1000 км.
Определены основные особенности и закономерности в поведении плазменных областей. В частности, показано, что для 150 км магнитное поле начинает заметно влиять на поведение плазменной области уже при c, и с этого же времени характер её поведения становится трёхмерным.
В результате выполненных в данной работе исследований достаточно подробно выяснена физическая картина развития возмущённой области и плазмы взрыва для умеренных и больших значений начального энерговыделения и для широкого диапазона высот.
На защиту выносятся следующие положения
1. Построен и реализован в виде комплекса программ эффективный монотонный численный алгоритм решения трёхмерных задач плазмодинамики, обладающий высоким порядком аппроксимации.
2. Разработанный алгоритм является эффективным для решения крупномасштабных МГД – задач о поведении плазмы в околоземном космическом пространстве при существенно различных начальных условиях – на высотах от 150 до 1000 км, особенно при расчетах на большие интервалы времени.
3. Комплекс программ применен для исследования основных особенностей и закономерностей в поведении плазменных областей в зависимости от:
начальной энергии взрыва;
высоты над поверхностью земли;
угла между точкой взрыва над поверхностью земли и плоскостью экватора.
Публикации
Научные результаты диссертации опубликованы в [1-11]. В совместных работах автору принадлежит разработка новых модификаций монотонных консервативных вариантов сеточно-характеристического метода и метода типа Годунова 2-3-го порядков аппроксимации, в том числе для уравнений трехмерной МГД, а также их алгоритмическая реализация в виде комплекса программ для решения трёхмерных задач плазмодинамики.
Апробация
Результаты, полученные в работе, докладывались на конференциях:
53rd Annual Meeting of the Division of Fluid Dynamics, November 19-21 2000, Washington, DC, USA.
First MIT Conference on Computational Fluid and Solid Mechanics, June 12 - 15 2001, the Massachusetts Institute of Technology, Cambridge, MA, USA.
Международный российско–японский симпозиум – Актуальные проблемы вычислительной механики, 5-10 августа 2002, г. С-Петербург, Россия.
Second MIT Conference on Computational Fluid and Solid Mechanics, June 17-20 2003, the Massachusetts Institute of Technology, Cambridge, MA, USA.
30th EPS Conference on Controlled Fusion and Plasma Physics, July 7–11 2003, St Petersburg, Russia.
18th International Conference on Numerical Simulation of Plasmas, September 7-10 2003, at the Oceanfront Resort and Conference Center in Falmouth, Massachusetts, USA.
31st EPS Conference on Plasma Physics, 28 June - 2 July 2004, Imperial College, London, Great Britain.
Third MIT Conference on Computational Fluid and Solid Mechanics, June 14-17 2005, the Massachusetts Institute of Technology, Cambridge, MA, USA.
Структура и объём диссертации
Диссертация состоит из введения, пяти глав и заключения. Общий объём диссертации 152 страницы. Список использованных источников содержит ссылки на 88 публикаций.
Обоснование выбора магнитогазодинамической модели с учетом диффузии магнитного поля
В [4] приведена постановка задачи для наиболее сложной по содержанию физических процессов стадии торможения плазмы в результате взаимодействия её ионов и атомов с ионами и атомами фона, а также с геомагнитным полем. Однако подробные расчёты выполнены лишь до времени завершения инерционной стадии, когда заканчивается формирование зарядового состава плазмы, определяющего характер её дальнейшего взаимодействия с фоном. Так как в настоящее время рассмотрение общего многоскоростного МГД -приближения является все ещё достаточно сложной задачей, расчёты стадии торможения выполнены лишь для вариантов со сравнительно невысокой начальной энергией плазмы, разлетающейся в геомагнитное поле. В этом случае применимо односкоростное приближение с учётом всего комплекса внутренних процессов, определяющих изменение степени ионизации и электронную и ионную температуры в плазме.
Как показано в [59], при высокой степени ионизации плазмы, когда доля тепловой энергии в её внутренней энергии невелика, использование дивергентной формы уравнения энергии для определения ионной и электронной температуры (а именно Те необходима для расчёта ионизационных процессов и диффузии магнитного поля), вообще говоря, не может дать удовлетворительной точности для их определения. Поэтому для ранней стадии разлета плазмы рассчитывались непосредственно уравнения для электронной и ионной температур. Однако для более поздней стадии развития возмущенной области, когда степень ионизации среды уже сравнительно мала, а трёхмерный характер её поведения становится определяющим, использование дивергентной формы уравнений позволяет, с одной стороны, надеяться на удовлетворительную точность расчёта внутренней энергии , плотности р и скорости У, а с другой - использовать разностные схемы, построенные для сквозного расчета разрывных решений с использованием именно дивергентной формы уравнений. Возникающая при этом неточность в определении температуры (в случае существенно разных долей кинетической и других видов энергии) обычно соответствует общему приближённому характеру постановки физической части задачи, которая определялась, прежде всего, ограниченностью вычислительных возможностей решения трёхмерных задач плазмодинамики. Так как основная цель данной работы состояла в создании трёхмерного алгоритма расчёта параметров крупномасштабного плазменного течения в верхней атмосфере и магнитосфере, то, учитывая ограниченность вычислительных возможностей, был сделан ряд упрощающих предположений по физике процессов, не влияющих существенно на общие закономерности поведения параметров возмущённой области:
- для плазмы и фона рассматривалась сплошная односкоростная плазменная среда;
- при t=0 плазма и фон разделены не только заданием различных термодинамических характеристик, но и степенью ионизации а = пе/п 5 где пе и п - концентрации электронов и тяжелых частиц;
- течение плазмы и фона предполагается замороженным, т.е. процессы ионизации и рекомбинации не учитываются и в каждой частице ее = const.
Монотонные разностные схемы высокого порядка аппроксимации для одномерных систем уравнений гиперболического типа
При численном решении задач, содержащих в области интегрирования разрывы искомых функций и другие особенности, важную роль играют монотонные разностные схемы (по другой терминологии - схемы с положительной аппроксимацией или мажорантные схемы), которые позволяют получать монотонные профили в разрывных решениях. Для уравнений и систем гиперболического типа понятие о таких схемах впервые введено Фридрихсом [63]. К этому классу принадлежат известные методы Куранта - Изаксона - Риса [64], С.К. Годунова [65], П. Лакса [66], сеточно-характеристический [67] и ряд других известных численных методов. Как было доказано в [65] для простейшего одномерного линейного уравнения переноса и для явных двухслойных линейных разностных схем на регулярных сетках монотонные по С.К. Годунову линейные схемы (2.2), для которых (и наоборот), могут иметь только первый порядок точности.
Эта теорема, как показано в [68,69], справедлива и для линейных схем с положительной по Фридриху аппроксимацией [63] в самом общем случае в том числе для многослойных и неявных схем. Это, конечно, несколько ограничивает область применения таких схем, однако, переход к нелинейным разностным схемам (в том числе для линейного уравнения (2.1)): аУм =а(т,Н,и)Ум, позволяет строить монотонные схемы с более высоким, чем первый, порядком аппроксимации. К этому классу монотонных схем относятся гибридные схемы ([70], [71] и др.), TVD схемы ([72] и др.), в которых, для схем вида (2.2), (2.4), монотонность определяется из условий.
В [68,69,73-76] для уравнений и систем гиперболического типа была построена общая теория линейных разностных схем с положительной аппроксимацией и близких к ним схем более высокого, чем первый порядка точности, основанная на использовании метода неопределенных коэффициентов и введении линейных пространств этих коэффициентов
В данной главе с использованием результатов [62] для простейшего уравнения переноса (2.1) и систем уравнений гиперболического типа рассматриваются монотонные схемы второго-третьего порядка аппроксимации, опирающиеся на характеристическое свойство точного решения (2.1), для которых, в случае выполнения условия Куранта 0 а = Лт/И \ выполняется характеристическое условие монотонности:
Построение монотонных разностных схем высокого порядка аппроксимации для одномерных систем уравнений газовой динамики
Обсудим теперь более подробно, почему так важно в случае присутствия разрывов большой интенсивности в решении, предиктор в (3.7) рассчитывать по схеме Годунова [65]. Как известно, если интенсивность ударной волны стремиться к бесконечности, то её ширина в рамках уравнений Навье-Стокса с линейной вязкостью стремиться к нулю, в то время как, для квадратичной вязкости ширина ударной волны не зависит от её интенсивности. В работе [82] было показано, что разностные схемы, использующие точное решение задачи Римана при расчете уравнений Эйлера (3.1), имеют на сильных ударных волнах схемную вязкость близкую к квадратичной. В то время как разностные схемы, построенные с использованием формул (3.7)-(3.8), являются схемами с линейной математической вязкостью. Близость схемной вязкости на сильных ударных волнах к квадратичной позволяет вести расчет ударных волн любой интенсивности с временным шагом, определяемым из условия Куранта. В случае схем с линейной вязкостью, начиная с некоторой интенсивности ударной волны, расчет с таким временным шагом становится невозможным из-за сильных численных колебаний. Шаг по времени приходится значительно уменьшать, что приводит не только к увеличению объёма вычислений, но и к уменьшению точности расчета, особенно в случае схем первого порядка точности. Поэтому, исходя из этих соображений, некоторые вычислители предлагают для расчета сильных ударных волн вводить в разностную схему искусственную квадратичную вязкость [83]. Это конечно помогает увеличить шаг численного интегрирования, но при этом заметно понижает точность расчета.
Запуск итерационного процесса при решении задачи Римана в методе Годунова [65] можно избежать в большинстве расчетных точек при помощи предварительной оценки значений решения. Как известно, решение задачи Римана можно свести к нахождению значения давления на поверхности тангенциального разрыва.
Построение монотонных разностных схем высокого порядка аппроксимации для одномерных систем уравнений магнитной гидродинамики
Для уравнений МГД с сильными разрывами, в принципе, можно построить аналог схемы С.К.Годунова [65] и использовать его в качестве предиктора в рассматриваемых монотонных схемах высокого порядка аппроксимации (2.35) или в качестве самостоятельной схемы, однако, технически это достаточно сложная задача, т.к. необходимо предусмотреть более сотни возможных комбинаций из ударных волн, контактных разрыв, волн разрежения, альфвеновских волн, быстрых и медленных магнитозвуковых волн и ряда вырожденных ситуаций (в отличие от уравнений газовой динамики, в которых таких ситуаций существенно меньше - различные комбинации из ударных волн, контактных разрывов и волн разрежения). Поэтому не случайно, что в литературе по МГД отсутствует описание подобных реализаций.
Вместе с тем, возможен промежуточный вариант, заключающийся в использовании в качестве предиктора (в общем алгоритме с использованием монотонных схем высокого порядка аппроксимации) схемы [65] для газодинамической части уравнений и построении аналога этой схемы (или использовании линейного приближения этой схемы, например, консервативного варианта сеточно-характеристического метода [68]) для остальной части уравнений (дополняющей до полных МГД уравнений их части). В данной работе на каждом шаге по времени по каждому из пространственных направлений использовалось именно такое расщепление «по физическим процессам» полных уравнений МГД на газодинамическую (1-й этап шага по времени) и дополняющую части. Поскольку магнитное поле Земли является относительно слабым, для «негазодинамической» части уравнений (2-й этап шага по времени) в качестве предиктора использовался консервативный вариант сеточно-характеристического метода [67].
Алгоритм выделения компонент односкоростного плазменного течения
В сделанных расчетах масштабы внешнего МГД-возмущения значительно превышали масштабы внутренней плазменной области, вследствие этого требовалось охватить большие пространственные масштабы ( 3000-г4000 км), достаточно подробно описывая вместе с тем, внутреннюю плазменную область, т.е. требовалась достаточно подробная расчётная сетка. Так как работоспособность алгоритма могла проверяться только на численном эксперименте с вариацией (E,h), где Е - энергия взрыва, оставшаяся в плазме после завершения инерционной стадии разлёта, а каждая версия программы требовала специальной работы по распараллеливанию, то исследование поздней стадии развития течения (/ 5 секунд) для ряда начальных условий производилась по двухмерному варианту программы. Сравнение с трёхмерными расчётами для умеренного энерговыделения даёт в среднем различие 25-г40 %, для очень больших энерговыделений - различие до 2-х раз.
Так как оптические, термодинамические и ионизационно-химические характеристики плазмы могут заметно отличаться от характеристик воздуха, то было целесообразно также усовершенствовать алгоритм в части выделения отдельных компонент односкоростного плазменного течения. Каждую компоненту среды удобно описывать относительной массовой концентрацией для которой уравнение неразрывности с учётом того, что —+div(pu) = О, преобразуется к виду Где St - скорость изменения /-ой компоненты в ионизационно-химических и диффузионных процессах. В представленных расчётах она полагалась равной нулю. Для численного решения этих линейных по концентрациям компонентов уравнений переноса использовались версии тех же, что и для МГД уравнений, численных методов и алгоритмов. Следует отметить, что сквозной расчет уравнений (5.3) на большие времена требует использования разностных схем высокого порядка аппроксимации (известная проблема сильного «размазывания» разрывов схемами низкого порядка аппроксимации из-за больших коэффициентов схемной вязкости в таких схемах).
В представленных ниже расчетах начальные данные для концентрации плазмы полагались единицей в области начального возмущения и нулем в невозмущенной атмосфере, а для атмосферного воздуха, наоборот: нулем в области начального возмущения и единицей в невозмущенной атмосфере. Это позволило следить с достаточной точностью за положением границы продуктов взрыва плазмы. В основных расчетах использовалась монотонная схема 2-3 порядка аппроксимации (2.35).
Возникающее под действием разлетающейся плазмы МГД-возмущение распространяется на тысячи километров, поэтому существенное влияние на его динамику оказывает не только ближняя от центра энерговыделения область ( 100-f300 км), ионизованная жёстким излучением, но и ионизация в дальней зоне, соответствующая естественным условиям. Поэтому в расчётах в качестве внешних условий задавались естественные значения электронной температуры Тф{К) и степени ионизации а (/г), через которые определялся коэффициент диффузии магнитного поля. Отметим, что на /г 50 км степень ионизации в естественной атмосфере близка к нулю и на эти высоты приходят лишь возмущения геомагнитного поля, среда же остаётся покоящейся. Для правильного расчёта этих возмущений необходимо знать проводимость грунта, и требуется специальная доработка самого расчётного алгоритма. В данной работе вопросы, связанные с поведением возмущённого геомагнитного поля вблизи поверхности Земли не рассматривались.
