Механика и электродинамика пристеночной плазмы

Содержание
Введение....................................................................5
Глава 1. Механика и электродинамика пристеночной плазмы в молекулярном режиме...........................................................................................10
1.1. Физическая и математическая модель задачи.............................. 10
1.1.1. Кинетическое уравнение Больцмана..................................... 10
1.1.2. Обобщенное уравнение Больцмана....................................... 11
1.1.3. Система уравнений Максвелла.......................................... 15
1.1.4. Математическая модель задачи......................................... 15
1.1.5. Начальные и граничные условия........................................ 19
1.2. Вычислительная модель задачи........................................... 21
1.2.1. Методы численного решения уравнения Власова..........................21
1.2.1.1. Метод крупных частиц Ю.М. Давыдова.................................21
1.2.1.2. Метод характеристик................................................23
1.2.2. Методы численного решения уравнения Пуассона.........................26
1.2.3. Описание алгоритма и компьютерной программы решения системы уравнений Власова-Пуассона................................................................................... 29
1.2.3.1. Основные особенности работы программного блока.....................29
1.2.3.2. Описание алгоритма расчетной программы.............................36
1.2.4. Методические расчеты и сравнение с экспериментом.....................43
13. Результаты вычислительных экспериментов для тела цилиндрической геометрии.....................................................................................................47
1.3.1. Функции распределения заряженных частиц..............................47
1.3.2. Эволюция интегрального тока на тело, зависимость установившегося значения тока от основных параметров расчета..................................................................60
1.3.3. Распределение плотности ионного тока на по обводу цилиндра........... 63
1.3.4. Распределение концентраций ионов, электронов и потенциала самосогласованного электрического поля в пристеночной области.................................................64
1.3.5. Поле средних скоростей ионов в пристеночной области.................. 68
1.3.6. Влияние магнитного поля на структуру возмущенной зоны................ 71
1.4. Результаты вычислительных экспериментов для тела плоской геометрии.. 79
1.5. Выводы из главы 1......................................................82
Глава 2. Механика и электродинамика пристеночной плазмы в режиме сплошной среды........................................................................................ 84
2.1. Состояние проблемы и задачи исследования...............................84
2.2. Физические и математические модели задачи механики и электродинамики столкновительной пристеночной плазмы..........................................................92
2.2.1. Слабоиоиизованная ламинарная плазма..................................93
2.2.2. Слабоиоиизованная турбулентная плазма................................95
2.3. Система начальных и граничных условий..................................96
2.3.1. Начальные условия....................................................96
2.3.2. Граничные условия....................................................97
2.4. Вычислительная модель задачи механики и электродинамики пристеночной плотной плазмы................................................................................. 103
2.4.1. Метод крупных частиц применительно к задачам электродинамики плотной плазмы...................................................................................................... 103
2.4.2. Методы решения уравнений Максвелла...................................109
2
2.5. Методические исследования и тестовые задачи........................... 110
2.6. Результаты математического моделирования обтекания тел
СЛабоИОННЗОВаННОЙ СТОЛКНОВИТСЛЫШЙ П.таЗМОЙ................................. 113
2.6.1. Цилиндрическое тело в поперечном потоке ламинарной столкновителыюй плазмы без магнитного поля................................................. 113
2.6.1.1. Профиль скорости нейтральной компоненты........................... 115
2.6.1.2. Поле скоростей электронов и ионов по обводу цилиндра.............. 116
2.6.1.3. Поле концентраций заряженных частиц............................... 121
2.6.1.4. Изолинии потенциала и распределение напряженности электрического поля. 122
2.6.1.5. Распределение плотности тока по обводу цилиндра................... 124
2.6.2. Цилиндрическое тело в ламинарном потоке плазмы в магнитном поле..... 126
2.6.3. Цилиндрическое тело в турбулентном потоке столкновителыюй плазмы 133
2.6.4. Плоский электрод в потоке слабоионизованной столкновительной плазмы..140
2.6.5. Обтекание цилиндрическою тела потоком слабоионизованной столкновительной плазмы при умеренных числах Рейнольдса.................... 143
2.7. Выводы из главы 2..................................................... 147
Глава 3. Механика и электродинамика пристеночной плазмы в переходном режиме...........................................................149
3.1. Состояние проблемы и задачи исследования.............................. 149
3.2. Математические модели задач механики и электродинамики пристеночной плазмы в переходном режиме..................................................153
3.2.1. Математическая модель задачи с учетом столкновений типа «ион-нейтрал» 153
3.2.2. Математическая модель задачи с учетом столкновений типа «электрон-нейтрал»................................................................... 154
3.2.3. Математическая модель задачи с учетом столкновений типа «ион-ион» и «ион-электрон»............................................................. 155
3.3. Вычислительные модели задач механики и электродинамики пристеночной плазмы в переходном режиме.................................... 157
3.3.1. Вычислительная модель с учетом столкновений типа «ион-нейтрал».......157
3.3.2. Вычислительная модель с учетом столкновений типа «электрон-
нейтрал»................................................................. 161
3.3.3. Вычислительная модель с учетом столкновений типа «ион-ион» и «ион-элсктрон»................................................................ 163
3.4. Результаты вычислительных и физических экспериментов по механике и электродинамике пристеночной плазмы в переходном режиме.................... 168
3.4.1. Влияние ион-атомных столкновений на процессы переноса в пристеночной плазме..................................................................... 168
3.4.2. Влияние «электрон-атомных» столкновений на процессы переноса в пристеночной плазме........................................................ 173
3.4.3. Влияние «ион-ионных» и «ион-электронных» столкновений на процессы переноса в пристеночной плазме............................................. 175
3.5. Выводы из главы 3..................................................... 177
Глава 4. Технические приложения по механике н электродинамике пристеночной плазмы................................‘........................177
4.1. Теория и методика электрического зонда (зоидовая диагност ика
плазмы).....................................................................177
4.1.1. Математические и численные модели зондовых задач.................... 180
4.1.2. Зонд в покоящейся слабоионизованной плазме при условии тонкого слоя объемного заряда............................................................183
3
4.1.3. Зонд в покоящейся плотной слабоионизированной плазме в случае произвольного слоя объемного заряда......................................... 186
4.1.4. Цилиндрический и сферический зонды в ламинарном потоке слабоионизованной столкновительноЙ плазмы без магнитного поля при условии
Ле — I, если слой объемного заряда тонкий столкновительный...................189
4.1.5. Цилиндрический зонд в ламинарном потоке слабоионизованной столкновительноЙ плазмы без магнитного поля в общем случае...................191
4.1.6. Цилиндрический зонд в потоке столкновительноЙ плазмы с магнитным полем 193
4.1.7. Цилиндрический зонд в поперечном потоке бесстолкновительной плазмы 194
4.1.8. Плоский пристеночный зонд в ламинарном потоке слабоионизованной столкновительноЙ плазмы без магнитного поля..................................197
4.1.9. Плоские изолированные и пристеночные зонды в потоке бесстолкновительной плазмы.......................................................................201
4.1.10. Двойные зонды. Их взаимное влияние.................................. 204
4.2. Нестационарный электрический зонд.......................................207
4.2.1. Алгоритм определения температуры ионов................................207
4.2.2. Измерительная аппаратура для работы с нестационарным зондом...........213
4.2.3. Зависимость времени релаксации и интегрального тока от
плазмообразующего вещества...................................................217
4.3. Электромагнитный рельсовый ускоритель с якорем в виде плазменного канала.......................................................................221
4.3.1. Физическая постановка задачи..........................................221
4.3.2. Математическая модель пристеночного слоя плазменного контакта.........222
4.3.3. Вычислительная модель задачи..........................................226
4.3.4. Результаты математического моделирования..............................227
4.4. Математическое моделирование возмущенной зоны и «собственной атмосферы» вблизи КЛА........................................................228
4.5. Электродинамические методы воздействия на ионизованный пограничный
слой........................................................................240
4.5.1 Инжекция в плазму пограничного слоя потока
отрицательных ионов..........................................................240
4.5.2. Воздействие на пограничный слой импульсом поперечного магнитного поля... 242
4.6. О возможности электродинамического управления вектором тяги плазменного движителя
(ГІД)....................................................................... 244
4.7. Выводы из главы 4...................................................... 245
Выводы из диссертации........................................................248
Литература..................................................................252
Приложение.................................................................. 261
4
Введение.
Плазмой называется ионизованный газ, в котором выполняется условие квазинейтральности. Это условие может соблюдаться лишь в объемах, характерные размеры которых велики в сравнении с некоторым масштабом длины М|. и на временах, больших но сравнению с масштабом времени М{. Как показано в классических монографиях по физике плазмы [1}*[3], [44], в качестве Мь выбирается радиус дсбаевского экранирования Го, а в качестве М( - период плазменных (или ленгмюровских) колебаний. В пристеночной области условие квазинейтральности может не соблюдаться. Причин для этого много: влияние потенциала тела, теплового и направленного движение заряженных частиц, различных видов эмиссии электронов или инжекции других частиц с поверхности, химических реакций и т.д. Между тем роль пристеночных слоев, как правило, оказывается определяющей при создании различных приборов и устройств, использующих плазму в качестве рабочего тела (плазменные движители, технологические плазмотроны, плазмохимические реакторы, МГД-генераторы, плазменная электроника), или работающие в плазме (космические летательные аппараты (КЛА) и космические станции (КС), гиперзвуковые летательыне аппараты (ГЛА), системы типа «Буран», «Шатл», «Колумбия» и др.). Особое место занимают методы зондовой диагностики, основанные на детальном исследовании процессов переноса заряда в пристеночной области.
Многочисленные работы в области механики и электродинамики пристеночной плазмы, обзор которых будет дан ниже по главам, не полностью охватывают важные для практики режимы течения, часто исследуют асимптотические и предельные случаи, не носят систематического характера. В диссертации исследованы с единых позиций все три режима течения: молекулярный, переходный и континуальный. Обнаружены и физически обоснованы новые нелинейные эффекты, возникающие при взаимодействии заряженных тел с плазмой. Особое внимание уделено приложениям результатов исследований в области авиационно-космической техники и зондовой диагностики. Основной метод исследования - математическое моделирование. Приведены обширные вычислительные эксперименты в различных режимах течения с учетом разнообразных физических процессов. Среди использованных вычислительных методов особое место занимает алгоритм метода крупных частиц Давыдова Ю.М. и его модификации применительно к задачам механики и электродинамики пристеночной плазмы. Результаты работы опубликованы в трех монографиях, десяти статьях в журналах, рекомендованных ВАК, доложены на международных и всесоюзных конференциях (более 25 опубликованных
докладов и тезисов докладов). Ряд результатов используется в заинтересованных организациях и в учебном процессе МАИ.
Диссертация состоит из введения, четырех глав, списка использованных источников и приложения.
Первая глава посвящена механике и электродинамике пристеночной плазмы в молекулярном режиме. В этой главе дан обзор работ отечественных и зарубежных авторов, внесших вклад в развитие данного направления, сформулирована физическая постановка задачи в достаточно общем виде, когда имеется направленная скорость движения плазмы, возможно патичие внешнего магнитного поля, возможна также эмиссия электронов или инжекция ионов с поверхности тела. Характерный размер обтекаемого плазмой тела, его потенциал относительно плазмы, величина направленной скорости потока, отношение температур ионов и электронов, величина индукции внешнего магнитного поля - варьировались. С целью сокращения размерности задачи в фазовом пространстве рассматривалось тело цилиндрической геометрии, расположенное поперек потока плазмы, а внешнее магнитное поле направлялось вдоль оси цилиндра. В такой постановке задача позволяла выявить вес физические особенности взаимодействия заряженного тела с разреженным потоком плазмы при наличии магнитного поля в и то же время она оказывалась нестационарной и 4х-мсрной в фазовом пространстве. Кроме того, после проведения необходимого количества вычислительных экспериментов появилась возможность развивать теорию цилиндрического зонда, которая в описанной постановке отсутствовала. Рассматривалось также тело плоской геометрии, расположенное вдоль потока плазмы. В заключение главы приведены результаты математического
моделирования, включая структуру возмущенной зоны вблизи заряженного цилиндра (поля скоростей, концентраций, потенциалов) в лобовой, боковой и теневой областях, влияние величины направленной скорости, внешнего магнитного ноля и других
параметров задачи на структуру функций распределения частиц и их моментов.
Вторая глава по своей структуре напоминает первую, но она относится к случаю столкновитсльной плазмы. Глава содержит обзор результатов работ предшествующих авторов, а также физическую, математическую и численную модели задачи.
Рассмагривается цилиндрическое заряженное тело, расположенное поперек потока
столкновитсльной плазмы. Магнитное поле может быть направлено вдоль оси цилиндра. В такой постановке задача оказывается достаточно обшей, и в то же время удобной для математического моделирования. Если длина цилиндра много больше его радиуса, то задача становится нестационарной, двумерной в геометрическом пространстве. Как и в первой главе, исследовалось также и тело плоской геометрии, ориентированное вдоль
6
потока. Использовалась стандартная система дифференциальных уравнений [4-5-6], учитывались все три физических процесса переноса - конвекция, диффузия и подвижность, однако система граничных условий была уточнена. Вычислительная модель основывалась на методе крупных частиц с уточнением на каждом временном слое значения электрического поля путем решения уравнения Пуассона. Использование достаточно общей математической модели и надежных вычислительных алгоритмов позволило получить целый ряд новых нелинейных эффектов, которые не были обнаружены в работах предшествующих авторов вследствие использования ими более грубых моделей. В частности, показано, что в определенных интервалах изменение характерных параметров задачи плотность тока но обводу цилиндра может иметь локальный максимум на боковой поверхности, а также еще один максимум в теневой области. Плотность электронного тока может проходить через максимум в зависимости от параметра Холла (в ограниченном диапазоне изменение чисел Маха и Рейнольдса). Как в молекулярном, так и в континуальном режиме показана возможность искривления возмущенной области за телом под действием осевого магнитного поля.
В третьей главе диссертации рассматривается переходный режим, когда число Кнудсена порядка единицы. Основная цель исследования - показать влияние столкновений на структуру функции распределения заряженных частиц и поведение моментов функции распределения, в первую очередь концентраций заряженных частиц и плотностей их токов. Ввиду сложности математических моделей в переходном режиме рассмотрение ограничивалось случаем покоящейся плазмы в отсутствие магнитного поля. Количество публикаций в случае переходного режима течение значительно меньше, чем в предельных режимах Кп»1 и Кп«1, поэтому результаты, полученные в этом режиме в строгой постановке, представляют существенный интерес. Изложение материала ведется в тон же последовательности, что и в гл. 1,2. Формулируется физическая, математическая и численная модели задачи и затем приводятся результаты вычислительных экспериментов. Отдельно исследуется роль столкновений типа «ион-нейтрал», «электрон-нейтрал» и «ион-ион», «электрон-ион». Если первые два типа столкновений предполагают решение кинетического уравнения Больцмана совместно с уравнением Пуассона, то столкновения между заряженными частицами описываются системой уравнений Фокксра-Планка и уравнения Пуассона.
В результате решения задачи получены функции распределения заряженных частиц с учетом столкновений, распределения концентраций и самосогласованных потенциалов, а также плотности токов на поверхность цилиндра. Показано, что столкновение типа «ион-нейтрал», «электрон-нейтрал» приводят к снижению плотности
7
тока по сравнению с молекулярным режимом, но она оказывается больше, чем в режиме сплошной среды.
Четвертая глава целиком посвящена многочисленным техническим приложениям результатов, полученных в предшествующих главах диссертации.
В параграфе 4.1 рассмотрены вопросы зондовой диагностики плазменных потоков в различных режимах течения. Методики привязаны к цилиндрическому зонду, расположенному поперек потока (рис. 1) и к плоскому зонду (пристеночному или выносному), расположенному вдоль потока. Приведены вольт-амперные характеристики зондов, полученные по уточненным математическим и численным моделям, сформулированным в гл. 1-3. Следует отметить, что в литературе подобного набора характеристик найти не удалось, хотя имеются характеристики для предельных асимптотических случаев, мало пригодные для практики. Попытки применить их для обработки экспериментально полученных зондовых характеристик приводят к ошибкам 200%-^400% и более [7], [8], [9]. Характеристики, приведенные в п. 4.1, позволяют избежать подобных ошибок.
В п.' 4.2 изложены идеи использования нестационарного зонда, когда на зонд подается ступенчатый импульс потенциала с фронтом нарастания Ц><< т, где т -характерное время релаксации в плазме. Как отклик на изменение потенциала происходит эволюция зондового тока, которая определяется концентрациями и температурами тяжелой компоненты плазмы - ионов. Приравнивая значение тчнся., полученное в численном эксперименте С Тэксп.у полученном в натурном эксперименте, получаем еще одно соотношение между параметрами плазмы, из которого можно определить ионную температуру, которая по классической зондовой теории не определяется. Приводятся обширные данные по значениям тЧИсл., полученным в численных экспериментах по моделям гл. 1-3.
В параграфе 4.3 приводятся расчеты плазменного якоря для электромагнитного рельсового ускорителя. Использованы модели главы 2 с учетом термоэмиссии электронов с плоского катода. Получены распределения нарамегров плазмы и электрических полей в пристеночной области плазменного контакта.
В параграфе 4.4. исследуется возмущенная область вблизи спутника цилиндрической формы, который движется в ионосферной плазме таким образом, что вектор скорости перпендикулярен оси цилиндра. Подробно исследованы лобовая, боковая и теневая области вблизи спутника. Получены поля концентраций, скоростей, электрических полей. Исследован случай, когда с поверхности спупшка инжектируется поток заряженных частиц, создавая тем самым особые условия вблизи него. Возмущенная
область плазмы вблизи спутника образует «собственную атмосферу», параметры которой важны для обеспечения жизнедеятельности спутника и проведения на нем физических экспериментов.
Параграф 4.5 посвящен методам электродинамического воздействия на ионизованный пограничный слой, возникающий вблизи ГЛА, в частности, показана возможность создания радиопрозрачпого канала. Рассмотрено два способа: инжекиия в пограничный слон потока отрицательных ионов и создание импульсного поперечного магнитного ноля. Методами математического моделирования, изложенного в гл. 2, показана возможность понижения концентрации электронов на порядок и более с целью его радиопросветления.
В парафафе 4.6. рассмотрена возможность управления вектором тяги плазменного движителя. Отклонение плазменного потока поперечным магнитным полем рассчитано модельно но отклонению в тех же условиях возмущенной зоны («следа»), возникающей за цилиндром, обтекаемым поперечным потоком плазмы. Магнитное поле располагается параллельно оси цилиндра. Рассчитана зависимость угла отклонения «следа» за телом от индукции магнитного поля.
Существует еще целый ряд технических приложений, которые могут быть рассчитаны разработанными в диссертации методами математического моделирования:
• вопросы диспергирования и перемешивания различных компонент электромагнитными силами;
• некоторые вопросы, касающиеся нанатехнологий;
• расчет приборов радиоэлектроники;
• расчет элек тромагнитной совместимости плазменных струй, истекающих из ПД и каналов радиосвязи и др.
Практически не затронут вопрос расчета самих плазменных движителей (ПД) различного типа.
В заключение считаю своим приятным долтм выразить благодарность чл. кор. РАН Пирумову У.Г., всем сотрудникам, преподавателям, аспирантам и студентам кафедры 807 и 806, с которыми проводилась совместная работа и обсуждения по отдельным разделам диссертации.
Особую благодарность хочу выразить Ю.М. Давыдову за плодотворное обсуждение и консультации но методам крупных частиц.
Выражаю благодарность зав. каф. 807 академику РАН Ганиеву Р.Ф. за постоянную поддержку и стимулирование работы над диссертацией.
9
Глава 1.
Механика и электродинамика пристеночной плазмы в
молекулярном режиме.
1.1. Физическая и математическая модель задачи.
Молекулярный режим обтекания внесенных в плазму тел характеризуется условием: число Кнудсена Кп = Х/гр —> оо, где А. - средняя длина свободного пробега частиц, гр - характерный размер задачи. Такой режим встречается при движении летательных аппаратов в ионосфере Земли и за ее пределами; в струях, истекающих из движителей малой тяги; в некоторых типах газового разряда (тлеющий разряд, дуга низкого давления, СВЧ-разряд и др.).
1.1.1. Кинетическое уравнение Больцмана.
Основное уравнение кинетической теории для одночастичной функции распределения f(r,v,t) было получено JI. Больцманом в 1872 году [10]:
Df/Dt = Jst(f), (1.1)
где Jst(f) - интеграл столкновений, D/Dt - субстанциональная производная:
D/Dt = б/dt + \-б/(дг) + F Щд\') (1.2)
r,v - радиус-вектор и скорость частицы, F - результирующая сила, отнесенная к единице массы. Уравнение (1.1) определяет процессы переноса в однокомпонентном газе, если учитываются только парные столкновения. Интеграл столкновений может быть записан в виде:
■Iя = J(f»V - f«f|>)g«ßbdbdvdv() (1.3)
Здесь gup - относительная скорость сталкивающихся частиц сх и ß, b - прицельный параметр, ф - азимутальный угол. Интеграл столкновений (1.3) должен удовлетворять законам сохранения. Для случая упругих столкновений в однокомпонентном газе выполняются соотношения:
J J$Vidv - 0 (i= 1,2,3) (1.4)
Здесь ф| - инварианты столкновений (ф1 = m; фг = rnv; фз = mv2/'2). Нели умножить почленно уравнение (1.1) на инварианты столкновений и проинтегрировать по всем скоростям частиц, то с учетом (1.4) из уравнения Больцмана можно получить
дифференциальные уравнения гидродинамики. Решив уравнение (1.1), можно вычислить моменты функции распределения, например, концентрацию частиц
n = Jfdv, (1.5)
и температуру:
3kT/2 = 0,5mjf-(v - v())2dv (1.6)
Здесь Vo — гидродинамическая скорость течения.
Аналитическое решение уравнения (1.1) удается получить лишь в частных случаях, обзор которых можно найти в работах [11-15]. Одним из них является распределение Максвелла, пригодное в локально термодинамически равновесном газе в отсутствие внешних сил. Для этого случая выполняется равенство
0 = Л (1.7)
из которого следует распределение Максвелла
Iм = п(т/(2якТ))зй exp[-mv2/(2k Г)], (1.8)
где v — тепловая скорость молекул. В случае слабоионизованного газа, когда столкновениями между заряженными частицами можно пренебречь (vc « ÖVca, где vc -частота столкновения между заряженными частицами; vca - частота столкновений между заряженными и нейтральными частицами; Ö - относительное количество энергии, теряемое заряженными частицами за одно столкновение), классическое уравнение Больцмана имеет вид:
Dfe/Dt = öfa/öt + v-Öfa/(ör) + Fa-ö fj(dv) = J„a, a = i,c, (1.9)
где Fa = qaE/ma - сила, действующая на заряженную частицу в электрическом поле в расчете на единицу массы, qtt - заряд частицы, индекс i относится к ионам, е — к электронам, а — к нейтральным частицам. Если плазма настолько разрежена, что столкновениями между заряженными частицами можно пренебречь, то интеграл столкновений Jaa = 0, и уравнение (1.9) упрощается:
3fa/at + v-öfa/(Sr) + Fa-af«/(0v) = 0, a = i,e, (1.10)
Это уравнение известно в литературе, как уравнение Власова.
1.1.2. Обобщенное уравнение Больцмана.
Кинетическое уравнение Больцмана играет огромную роль в физике, особенно в гидродинамике, в теории процессов переноса, в космологии. Однако с момента своего создания и до настоящего времени больцмановская физическая кинетика, основанная на
уравнении (1.1), подвергается критике ввиду имеющих место внутренних противоречий. Отметим некоторые из них, основываясь на анализе, проведенным Б.В. Алексеевым [16-18].
1. Уравнение Больцмана справедливо для масштаба времени т>. (т>. - время между столкновениями частиц) и тг (тг — гидродинамическое время течения), но оно не работает на временах тв(т„- время взаимодействия сталкивающихся частиц). Поэтому ряд физических явлений, в которых существенен масштаб времени тв> выпадает из рассмотрения. Это, в частности, касается теоретических основ современных нанотехнологий.
2. Обычно одночастичная функция распределения, входящая в (1.1), нормируется на число частиц в единице объема. Если частицы рассматривать как материальные точки, то такая нормировка возможна. Но при вычислении интеграла столкновений, например, с использованием модели твердых сфер, необходимо знать сечения столкновений, т.е. рассматривать частицы конечного размера. Попадание центра масс частицы в контрольный объем не означает, что вся частица находится в этом объеме. В произвольный момент времени всегда найдутся частицы, которые располагаются частично внутри, а частично снаружи контрольной поверхности, что ведет к флуктуациям массы, а, следовательно, и других гидродинамических величин. Возможность таких флуктуаций уравнение (1.1) не предусматривает.
3. Как показано в пункте 1.1., из кинетического уравнения Больцмана вытекают гидродинамические уравнения: уравнение неразрывности, движения и энергии. Последнее при некоторых предположения может быть сведено к уравнению теплопроводности параболического типа, связывающему первую производную от температуры по времени и вторую производную но координатам. Его решение показывает, что импульс температуры, созданный в источнике, может дать скачок температуры в любой достаточно удаленной точке за бесконечно малый интервал времени. Это означает, что происходит бесконечно быстрое распространение теплоты, что абсурдно с молекулярно-кинетической точки зрения. Физики «исправили» это противоречие, введя в уравнение теплопроводности вторую производную по времени и сведя тем самым его к уравнению гиперболического типа. «Исправленное» уравнение приводит к конечным скоростям тепловых возмущений. Такое «исправление» неизбежно ведет к новой гидродинамической теории, которая в свою очередь должна быть следствием кинетического уравнения, отличающегося от уравнения Больцмана (1.1).
Приведенных примеров достаточно, чтобы обосновать необходимость вывода более общего уравнения, чем уравнеЕшс (1.1), которое было бы лишено указанных
12
противоречий. Среди большого числа попыток создания обобщенного уравнения Больцмана [19-24] наиболее последовательное и логически обоснованное с нашей точки зрения принадлежит Б.В. Алексееву [16-18]. Ниже будут изложены основные идей и выводы, вытекающие из работ Б.В. Алексеева (Впервые работа Б.В. Алексеева была доложена на лекциях по физической кинетике, прочитанных для научных сотрудников Софийского университета 1987 году, а в 1988 году опубликована в тезисах доклада на седьмом совещании по механике реагирующих сред (Красноярск)). Прежде всего при выводе обобщенного кинетического уравнения были учтены все три характерных масштаба (от масштабов времени перейдем к масштабам длины, умножив масштаб времени на масштаб скорости):
• гв — / — масштаб взаимодействия (/ - длина Ландау, на которой кинетическая энергия теплового движения кТ равна потенциальной энергии взаимодействия зарядов:
кТ = е2/(4лЕоО->/ = е2/(4я£0кТ));
• г* = X - средняя длина свободного пробега;
• гр - Б - гидродинамический масштаб, например, диаметр трубы, в которой течет газ.
Обычно
гв« X « Ь (1.11)
Как и в случае классического уравнения Больцмана (1.1), вывод обобщенного уравнения начинается с уравнения Лиувилля для Ы-частичиой функции распределения
аъ/а+2>г(аг»/дг,) + ^-(з^/зу,) =о (1.12)
1=1 1=1
Затем применяется метод Боголюбова-Борна-Грина-Кирквуда-Ивона (ББГКИ) [25]. Безразмерные уравнения для Б-часгичной функции распределения fs записываются с использованием пространственного масштаба г„. Далее применяется метод многих масштабов [26]. В данном случае имеется зри характерных масштаба: г„ = /, = А., гг = Ь.
В методе многих масштабов представляется в виде асимптотического ряда:
4 = (1.13)
1=0 А .
в котором ^ зависит от всех трех групп переменных. Крышка сверху означает
безразмерную величину (в дальнейшем будет опускаться), е* - малый параметр (е = пг„3 -число частиц в объеме взаимодействия).
13
Далее осуществляется стандартная процедура перехода к одночастичной функции распределения, в результате чего получается обобщенное уравнение Больцмана для слабоионизованного газа:
I- D/Dt(xaa- Dfa/Dt) = J„a, а = i,e (1.14)
В (1.14) таа - среднее время между столкновениями нейтральных (а) и заряженных (а) частиц. Интеграл столкновений ){ia - остается в больцмановской форме (1.3).
Рассмотрим случай максвелловского закона взаимодействии между заряженными и нейтральными частицами, когда сила межмолекулярного взаимодействия обратно пропорциональна пятой степени расстояния:
Fa. = ХтеЛ5- (1-15)
В этом случае таа может быть вычислено [17]:
Таа = |0,422-2я(х*(ша + ша)/( ma-m«))10-«!*]’1- (116)
Нели ввести
А = (8/3)(я)|/2Г(5/2)0,422(та + m,)'(x/(m.M))l/2 ,
М = ma/(ma + ma),
то частота столкновений v заряженных и нейтральных частиц может быть представлена в виде
v = n*(ma + ша>А, (1.18)
1фИ ЭТОМ Таа = V*1.
Обобщенное уравнение Больцмана в данном случае записывается следующим образом:
di'jdi + va.5fa/(ara) + Fa-dfa/(öva) -TaaC^fa/St2 + 2 Fa-^/CßVaa) + V„) : FaF„) = J„.
Здесь двоеточие означает двойное скалярное произведение тензоров.
В общем случае при использовании уравнения (1.14) следует иметь ввиду, что
D/Dt = а/аt+у„-а/(Эга) + Fe£ a/(3v0); F„ = f« + f«., (i .20)
где Ka = qa-E/ma - внешняя сила, Faa- сила межмолекулярного взаимодействия.
14
1.1.3. Система уравнений Максвелла.
Внешнее электрическое ноле, входящее в уравнения (1.9), (1.19) (а также силы, возникающие в магнитном поле, если оно существует), определяются решением системы дифференциальных уравнений Максвелла (27J:
div В = 0 rot Е = -dB/dt
div D = р rot Н = jnp + dDldt (1.21)
В = ppoH, D = ес0Е, jI1? = аЕ В случае использования обобщенного уравнения Больцмана система уравнений Максвелла должна быть модифицирована. Ее традиционная формулировка (1.21) не содержит уравнения неразрывности, однако получено с его использованием:
dp/dt + div jnp = 0, (1.22)
где р - заряд единицы объема, j - плотность тока проводимости. Величины р и j получены без учета флуктуаций плотности, возникающих вследствие конечных размеров
элементарных зарядов. Если учесть флуктуации, то в уравнениях (1.21) следует писать: р = рср - р*; j = jCp — j*, где р* и j* вычисляются в рамках обобщенной больцмановской теории [18].
Как показано в [28], в случае зондовых задач зондовые токи достаточно малы, и их магнитными полями можно пренебречь. В этих условиях процессы переноса заряда определяются электрическим полем, зависящим от потенциала зонда и пристеночного слоя объемного заряда. Поэтому система уравнений Максвелла сводится к уравнению 1 lyaccoim:
Дф = е(пе - n,)/8o, Е = -Уф, (1.23)
где ф - потенциал самосогласованного электрического ноля, п1>е - концентрации
заряженных частиц, е= 1,6-10'19 Кл, со = 8,85-1 О*12 Ф/м.
1.1.4. Математическая модель задачи.
В общем случае задача механики и электродинамики пристеночной плазмы в молекулярном режиме зависит от шести фазовых переменных и времени. С целью сокращения необходимых ресурсов ЭВМ, по в то же время для сохранения достаточной общности задачи были выбраны следующие постановки:
1. Внесенное в разреженную плазму тело имеет форму цилиндра радиуса гр, длины Ьр » гр и находится под потенциалом <рр.
Рис. 1.1. Расположение тела в потоке.
Направленная скорость у» указана на рисунке. Внешнее магнитное поле индукции В может быть направлено вдоль оси цилиндра. В этих условиях с учетом симметрии задачи и условия Ц » Гр В цилиндрической системе координат (Г,9,7,Уг,У0,У7) функции распределения заряженных частиц зависят только от четырех фазовых переменных и времени:
/а = и г,0,уг,уо,0 (1.24)
Кинетическое уравнение в классической больцмановской постановке в указанных координатах записывается так [5]:
дг
+ V
дг аг
г
от
г ае
— 4-
4 + Да.Е,
г т..
ду.
+
Ча р
0
ЧШа Г
ЗГ.
= 0
о\(
(1.25)
Оно дополняется уравнением Пуассона для самосогласованного электрического поля:
(1.26)
а2<р 1 аш 1 а2ф 1 ^
г + -——Хч«“« ;Е = .УФ,
о ео
аг2 г аг г2 а 2
где
пЛг>е>0 =
2 к Г
Ч та
+00 +00
\ /Га(Г>Мг^0>1)^гС1У0.
(1.27)
11лотность тока частиц сорта а на поверхность цилиндра и интегральный ток на цилиндр единичной длины имеют следующий вид:
16
ІаМ) =
2кТ.
і
N 2 0 +оо
Ш
а у
-00 -00
1а(0 = ГР /Ісх (1,0>ІО
(1.29)
Система начальных и граничных условий будет обсуждена ниже в п. 1.1.5.
Если в окрестности цилиндра существует продольное однородное магнитное ноле, то уравнение (1.25) выглядит гак:
Ёк+У ЁЪ. + Ь.ЁЬ_ +
Зі ' дг г 30
Я^Ео_^_Я^Вуг
гп г т „
+
Яа г; , Яа
Ш
Ег +
В V,
т
дГ
дч.
+
д1
(1.25')
а _
д\,
= О
2. Рассматривается пластина достаточно большого размера, нормаль к которой перпендикулярна вектору направленной скорости плазмы. На пластине выделим участок в виде удлиненного прямоугольника, электрически изолированный от пластины. Направим ось У по нормали к пластине, ось X - вдоль вектора скорости, ось Z - вдоль удлиненной стороны прямоугольника (рис. 1.2).
Рис. 1.2. Плоский пристеночный электрод
1-активная поверхность электрода;
2-обтекаемая плазмой пластина.
Размер короткой стороны прямоугольника 2гр, его потенциал <рр. Выделенный прямоугольник можно рассматривать, как часть боковой поверхности спутника, а в зондовой теории - как пристеночный зонд ленточного типа. Если имеется магнитное поле, то оно должно быть постоянным и направленным вдоль оси Z. В данной постановке
17
задача оказывается четырехмерной в фазовом пространстве и представляется достаточно общей. Запишем систему уравнений Власова-Пуассона для элемента конструкции ленточного типа в декартовой системе координат (Х,У,2):
ас
ді
<2
— + V.
Ё!*
дх
+ V.
Я. і ч“
ду ша I х йу/ ч х
д\
= 0
у
о ф д (р 1-е-.
(1.30)
(1.31)
Ґ /-ч» «т. Л 2 +со +0°
Па =
2кТ
V а / -оо-со
| ]Х(х,у,ух,уу,і)сІУх<ІУу.
(1.32)
Л* =
2кТ
і
\Т
Ч /
+оо О
Ча / /їа(Хр.Ур^х,Уу,і)уу<1УхСІУ;
(1.33)
Система (І.ЗО)-(І.ЗЗ) составляет математическую модель данной задачи, начальные и граничные условия будут рассмотрены ниже.
Системы уравнений, приведенные выше, приводились к безразмерному виду. При этом использовалась следующая система масштабов:
масштаб потенциала М<р = кТУе, (1.34)
масштаб длины Мь = ((еокТ,/(с2п„))|/2, (1.35)
масштаб концентрации м„ — Прэ, (1.36)
масштаб скорости Муа = (2кТа/ш«)1/2, а = і,с, (1.37)
масштаб функции распределения МГа = Мп/(Муа)3, (1.38)
масштаб времени М( = Мь/Муь (1.39)
масштаб напряженности электрического поля Ме - Мф/Мь (1.40)
масштаб плотности тока = еМпМУІ, (1.41)
масштаб интегрального тока Мі = Мі(Мь)2, (1.42)
масштаб магнитной индукции М» = 2Ме/Муі. (1.43)
Учитывая приведенные выше масштабы введем следующие обозначения:
• го = г^/Ме - безразмерный радиус цилиндра;
• фо = Фр/Мф - безразмерный потенциал тела;
• vo = УоJ М* - безразмерная направленная скорость плазмы;
• Во = В/Мв - безразмерная величина магнитной индукции.
(1.43')
18
Приведем для примера безразмерную систему уравнений (1.25'),(1.26)-(1.29), полученную при использовании масштабов (1.34) - (1.43):
д{.
ді
дг
г дв
V?
ч/57^-+л/57|^Ёг + ^Ву0 Г 2єа ца
дї.
ду.
— 4-
А А
V V,
а 1 Н-а
зг
ау,
^ = о
(1.44)
а^ф і аф і а2ф
^2 1 ; а; 1 ;2 -*А2 ^ і »
ар г а? г" ае-
(1.45)
+О0 +<©
пв(г,М) = | |Га(г,6,Уг,У0,ї)сІУгсІу,
-00 -со
(1.46)
О +х
ІІ«(ї»0)И>/8 І |га(г0,6,уг>у0,{)угсіугсіу(
(1.47)
(1.48)
є« = Ца= 6а = г« = 1 при а = і,
тс шс Тс т.
£а=ТГ,^и=—,Йа=-Г — Т. ГИ: Т; ГП
, га = -1 при а = с.
І І 1 с
В дальнейшем символ безразмерной величины «Л» будем опускать.
1.1.5. Начальные и граничные условия.
Для решения системы дифференциальных уравнений математических моделей, сформулированных в п. 1.1.4., необходимо задать начальные и граничные условия. Для уравнений (1.25), (1.30) необходимо задать начальную функцию распределения для частиц каждого сорта. Кроме того, должно быть задано граничное условие для Га на границах расчетной области. Для решения уравнения Пуассона (1.26), (1.31) необходимы два граничных условия для потенциала: значение <р при г = гр и его значение на внешней границе расчетной области, которое, как правило, считается нулевым. В ряде случаев [29]
19
в качестве начальных распределений частиц выбирают равновесное распределение Максвелла
где п10 - концентрация частиц в нсвозмущенной плазме, Та - температура компоненты а, V - скорость частицы, \т - вектор скорости набегающего потока.
Как показали методические исследования, начальное распределение существенно влияет на процесс установления решения, стационарное решение от начальной функции распределения зависит мало. Поскольку точное значение начального распределения, как правило, неизвестно, задачу приходится решать дважды. Сначала берется в качестве начального максвелловское распределение. Полученное на этом этапе стационарное решение берется за начальное распределение на втором этане.
На внешней границе возмущенной зоны используется либо функция распределения Максвелла, либо иное распределение, вычисленное в отдельной задаче. Выпишем выражение (1.49) для случая обтекания цилиндрического тела:
Граничные условия на поверхности обтекаемого плазмой тела зависят от физических процессов, которые могут происходить на ней. Для металлической поверхности в простейшем случае ставится условие идеального поглощения, согласно которому электрон, попавший на стенку, пропадает, а ион отдает свой заряд и возвращается в плазму в виде нейтрального атома. При этом отсутствуют потоки заряженных частиц от стенки в плазму.
В реальных условиях вслсдствии фото-, вторичной или термоэмиссии возможны электронные потоки с поверхности тела, функции распределения которых достаточно хорошо изучены [5,29,30].
Г« = Поо(та/(2лкТ«))3'2 схр[-та(у - <)2/(2кТа)],
(1.49)
1а(0,Г,0,Уг,У0) = (1Ью/я)(ша/(2кТа))3/2.
ехр[-т«{(уг + у^собО)2 + (^о - ужБІпО)2}/(2кТа)],
(1.50)
и участка пластины ленточного типа
(1.51)
20
1.2. Вычислительная модель задачи
1.2.1. Методы численного решении уравнения Власова.
1.2.1.1. Метод крупных частиц Ю.М. Давыдова.
Для решения уравнения Власова удобно использовать явную схему метода крупных частиц Ю.М. Давыдова [31-351, дополненную алгоритмом вычисления на ка*ждом временном слое самосогласованного электрического поля. Метод первоначально был разработан применительно к задачам газовой динамики для расчета сжимаемых течений сплошной среды [31] и в дальнейшем распространен на самый широкий круг задач, в который вошли и задачи расчета пристеночных слоев, возникающих вблизи заряженных тел, внесенных в плазму [28], [4], [5].
В соответствии с методом моделируемая среда заменяется системой частиц, совпадающих в данный момент времени с ячейками эйлеровой сетки. Расчет каждого временного шага разбивается на три этапа. На первом этапе рассматривается изменение за время А1 импульса и энергии крупной частицы, при этом её граница не смещается относительно начального положения. Па втором этапе моделируется движение частиц через границы эйлеровых ячеек и происходит перераспределение частиц по пространству. На третьем этапе рассчитывается перераспределение массы и заряда но пространству.
Выполнение граничных условий задачи обеспечивается введением рядов фиктивных ячеек, шаг по времени выбирается из условия Куранта (как показано в п.2, условие Куранта оказывается слишком жестким для задач электродинамики пристеночных слоев плазмы).
Рассмотрим подробно алгоритм метода крупных частиц для решения уравнения Власова на примере задачи обтекания плазмой неподвижной заряженной сферы [36]. Система уравнений в безразмерном виде имеет следующий вид:
Фазовое пространство в данный момент времени разбивается на ячейки с по.мошыо эйлеровой сетки:
Т- + 0-У)
(1.52)
/1,=0 = к'гП ехр(-у2).
(1.53)
где г е [г0, г„], и= [О, уе [-1,1], И, чу шаги эйлеровой сетки с числом узлов Л^г.г.г
Для каждой ячейки сетки ф... вводится величина
1|К
01*=$/''*** 7. о»у*
где л* - номер шага но времени, Мгскх1у ~ элемент объема в фазовом пространстве, ]=&кгг1у2. Таким образом, бде - нормированное количество частиц в ячейке в
момент времени г\
Ввиду неразрывности плотности среды [36], уравнение для 0=#АгАу£{, соответствующее уравнению Власова (1.52), имеет вид
Эб Эб уЕдб „ 2Ч(V Е\д<2 „ (V Е\-
+ + + ) -“ + — Ь^ = 27 ~ + — б- (1.54)
Э* Эг 2 да ч/* 2а/ ду чг IV)
В отличии от/величина О уже не постоянна вдоль траекторий в криволинейном пространстве из-за появления в (1.54) правой части типа источника.
Стационарное решение задачи, если оно существует, получается в результате установления, поэтому весь процесс вычислений сводится к многократному повторению вычислительного алгоритма.
Расчет каждого временного шага в соответствии с общим алгоритмом производится в три этапа.
1. Эйлеров этап.
Пренебрегаем всеми эффектами, связанными с перемещением крупных частиц. Потоков массы, заряда, импульса через границы ячеек нет. В соответствии с заданным распределением 0^к находится поле концентраций лД а затем с помощью численного решения уравнения Пуассона (п. 1.2.2.) определяется электрическое иоле Е,\
2. Лагранжсв этап.
На этом этапе вычисляются эффекты переноса, связанные с перемещением центров крупных частиц под действием электрического поля. Перемещение ячейки (у,А) на текущем шаге по времени находится с помощью уравнений движения
й=?’
с начальными условиями
г(г') = П, &(/’)= »1, У(‘‘) = УК
При этом используется напряженность электрического ноля, вычисленная на эйлеровом этапе.
Величина О меняется в соответствии с уравнением
22
тогда лагранжена ячейка будет содержать в момент I = ? + А( количество частиц
йук 1 + й(* * + ДО* Таким образом, получаем смещенное положение крупных частиц в
момент Г5 + Ал
3. Заключительный этап.
Величина (2£1 в ячейках исходной эйлеровой сетки определяется суммарным
вкладом всех ячеек, наложившихся на данную ячейку в момент Ґ + А/. По известным распределениям 0£1 вычисляются значения макропараметров в момент Ґ + А Л Например,
Схема метода крупных частиц для уравнения Власова консервативна, что позволяет повысить точность вычислений. Кроме того, в [36] показано, что численная диффузия идет в направлении сноса и слабо зависит от размерности фазового пространства.
Точность метода существенно зависит от размеров области, на которой располагается эйлерова сетка. Размер области по г уточняется в процессе эволюции из-за появления дрейфа носителей заряда. При этом можно использовать алгоритм плавающей сетки, разработанный в [36], что позволяет существенно снизить размер расчетной сетки в просгранстве скоростей.
1.2.1.2. Метод характеристик.
Уравнение Власова принадлежит к числу квазилинейных дифференциальных уравнений в частных производных первого порядка с нулевой правой частью. Для решения уравнений такого вида широко применяется численный метод характеристик. Алгоритм метода основан на том, что функция распределения в рассматриваемом случае постоянна вдоль характеристик, которые одновременно являются траекториями движения частиц плазмы в фазовом пространстве координат-скоростей. Этот факт позволяет на любом текущем времени расчета последовательно перебрать все узлы расчетной сетки и для каждого узла проинтегрировать систему уравнений характеристик от текущего момента времени до нуля. Получив решение системы для нулевого момента времени и
23
подставив его в начальное условие для уравнения Власова мы получаем искомое значение функции распределения для рассматриваемого узла расчетной сетки. При этом возможны следующие случаи:
1. В процессе интегрирования решение системы уравнений характеристик попало за внешнюю границу расчетной области гж. В этом случае в соответствии с тем, что за пределами расчетной области плазма считается невозмущенной и там отсутствует электрическое поле, значения скоростей, соответствующие данной характеристике, уже не могут изменится при дальнейшем интегрировании системы, поэтому процесс интегрирования можно продолжать только для тех уравнений характеристик, которые связаны с изменением координат.
2. Решение системы уравнений характеристик попало в область расположения обтекаемого плазмой тела. В этом случае существование рассматриваемой частицы плазмы не имеет физического смысла, т.к. тело по условию задачи не является источником заряженных частиц. В соответствии с этим прекращаем процесс интегрирования системы уравнений характеристик и функцию распределения соответствующей частицы приравниваем к нулю.
3. Решение системы уравнений характеристик попало за пределы расчетной области в пространстве скоростей, фиксированной для начального момента времени расчета. Если при этом начальным условием для уравнения Власова является сеточная функция распределения, то прекращаем процесс интегрирования системы уравнений характеристик и функцию распределения соответствующей частицы приравниваем к нулю (носитель функции распределения по скоростям изначально должен выбираться гак, ч тобы за его пределами данная функция была достаточно мала, и можно считать сё равной нулю без ущерба для точности расчета). Если начальным условием является аналитическое выражение (например, распределение Максвелла), то процесс интегрирования можно продолжать, поскольку в этом случае функцию распределения можно вычислить на интервале скоростей от -со до -со.
Теперь рассмотрим более подробно этапы метода характеристик на примере задачи обтекания цилиндра разреженной плазмой. Система уравнений в цилиндрических координатах имеет вид (1.25)-(1.29) с начальным максвелловским распределением (1.50). Уравнению Власова (1.25) соответствует система уравнений характеристик:
24
(1.55)
Условия на текущий момент времени
г’(г)-г. в*(о=е, уг\о=у„ од-у,.
(1.5 6)
В качестве 01раниченного носителя функции распределения рассмотрим следующий компакт в фазовом пространстве:
Ч»|-К.СОвв]. ОД-^ + ЦзтО. Р^+К.ЗШО]}.
Здесь го - радиус цилиндра, г<с - условная внешняя граница возмущенной зоны, Утах -граница обрезания «хвоста» максвелловского распределения. По переменной 0 используется симметрия возмущенной зоны. В пространстве скоростей носитель сдвинут так, чтобы центр тяжести начальной функции распределения оказался в центре квадрата со стороной 2 Утах.
Шаг по времени выбирается максимально большим при условии приемлемой точности расчета. Для выбора оптимального шага по времени проводятся методические расчеты.
Вводятся сеточные функции (/а)сИ.(£г)у»(£о)^- Индексы .у,/,у,к,1 соответствуют I, г, О, Уг, У(I. Значения (^г)^»(^в)у вычисляются после решения уравнения Пуассона для самосогласованного электрического ноля с правой частью, определяемой значениями функции распределения па предыдущем шахе по времени. При расчете нет необходимости хранить значения функций распределения на предыдущих шагах по времени, т.к. они не используются при интегрировании уравнений характеристик.
Решаем «обратную» задачу Коши для системы (1.55) с условиями (1.56) на отрезке времени [/5,0]. В итоге значение функции распределения в узле (г,-, 0,, У*, Уф) будет определятся по формуле
Для выбора оптимального метода интегрирования системы уравнений характерисгик проводятся методические расчеты. Метод выбирается из соображений достаточной точности расчета при установлении за приемлемое машинное время. В работе использовался метод Рунге-Кутта-Фсльберга пятого порядка [38]. В процессе
*ирр/в = {(г,е.^,Ув):г€[г0,г-1, вб[0,я). Уг еНи-У.СО50,
25