6
МАТЕМАТИЧЕСКОЕ И КОМПЬЮТЕРНОЕ
МОДЕЛИРОВАНИЕ MODELING AND SIMULATION
УДК 004.942, 533.9
МОДЕЛИРОВАНИЕ ИМПУЛЬСНОГО ПРОБОЯ В ГЕЛИИ С ИСПОЛЬЗОВАНИЕМ АДАПТИВНЫХ МЕТОДОВ С.И. Елисеев3, В.И. Демидовь, c, А.С. Чирцов1 c, В.И. Колобовй, А.А. Кудрявцев1 c, Е.А. Богдановь
a СПбГУ, 198504, Санкт-Петербург, Россия, step.eliseev@yandex.ru ь Университет Западной Вирджинии, WV 26506, Morgantown, West Virginia, USA c Университет ИТМО, 197101, Санкт-Петербург, Россия d CFD Research Corporation, AL 35806, Huntsville, Alabama, USA
Аннотация. Исследуется возможность использования адаптивных методов для моделирования процессов, происходящих при электрических пробоях в газах. Численное моделирование выполнялось для разряда в гелии при атмосферном давлении между игольчатыми электродами. Физическая модель сопровождающих пробой процессов строилась на базе самосогласованной системы уравнений непрерывности для потоков заряженных частиц (электронов и положительных ионов) и уравнения Пуассона для электрического потенциала. Резкая неоднородность плазмы в области стримеров потребовала использования при моделировании адаптивных алгоритмов построения вычислительных сеток. Дано описание метода адаптивного построения сетки и обоснование его эффективности для решения резко нестационарных задач пробоя газа при атмосферном давлении. Численное моделирование процессов осуществлялось в расширенном варианте свободно распространяемого пакета Gerris. Первоначально ориентированный на решение нелинейных задач динамики жидкостей пакет оказался пригодным для моделирования процессов в нестационарной плазме, описание которых строится на базе уравнений непрерывности. Использование адаптивных сеток позволило получить адекватную численную модель развития пробоя в системе игольчатых электродов. Динамика пробоя проиллюстрирована контурными графиками распределения концентраций электронов и напряженности электрического поля, полученными в ходе решения. Показан и проанализирован механизм пробоя с образованием положительного и анодонаправленного отрицательного стримеров. Показано соответствие между адаптивным построением вычислительной сетки и образующимися в ходе решения градиентами параметров плазмы. Результаты работы могут быть взяты за основу при проведении полномасштабных численных экспериментов по пробою газового промежутка.
Ключевые слова: моделирование плазмы, пробой газа, адаптивные методы, стримеры, импульсный разряд.
SIMULATION OF PULSED BREAKDOWN IN HELIUM BY ADAPTIVE METHODS S.I. Eliseev3, V.I. Demidovb,c, A.S. Chirtsova,c, A.A. Kudryavtseva,c, V.I. Kolobovd, E.A. Bogdanovb
a Saint Petersburg State University, 198504, Saint Petersburg, Russia, step.eliseev@yandex.ru b West Virginia University, WV 26506, Morgantown, West Virginia, USA c ITMO University, 197101, Saint Petersburg, Russia d CFD Research Corporation, AL 35806, Huntsville, Alabama, USA
Abstract. The paper deals with the processes occurring during electrical breakdown in gases as well as numerical simulation of these processes using adaptive mesh refinement methods. Discharge between needle electrodes in helium at atmospheric pressure is selected for the test simulation. Physical model of the accompanying breakdown processes is based on self-consistent system of continuity equations for streams of charged particles (electrons and positive ions) and Poisson equation for electric potential. Sharp plasma heterogeneity in the area of streamers requires the usage of adaptive algorithms for constructing of computational grids for modeling. The method for grid adaptive construction together with justification of its effectiveness for significantly unsteady gas breakdown simulation at atmospheric pressure is described. Upgraded version of Gerris package is used for numerical simulation of electrical gas breakdown. Software package, originally focused on solution of nonlinear problems in fluid dynamics, appears to be suitable for processes modeling in non-stationary plasma described by continuity equations. The usage of adaptive grids makes it possible to get an adequate numerical model for the breakdown development in the system of needle electrodes. Breakdown dynamics is illustrated by contour plots of electron densities and electric field intensity obtained in the course of solving. Breakdown mechanism of positive and negative (orientated to anode) streamers formation is demonstrated and analyzed. Correspondence between adaptive building of computational grid and generated plasma gradients is shown. Obtained results can be used as a basis for full-scale numerical experiments on electric breakdown in gases.
Keywords: plasma modeling, gas breakdown, adaptive methods, streamers, pulsed discharge.
Методы численного моделирования играют важную роль при описании низкотемпературной плазмы и явлений переноса в газах, так как сложность и разнообразие происходящих в них процессов не позволяют получить их всестороннего описания обычными экспериментальными и теоретическими методами [1-4]. Так, например, в случае тлеющего разряда в воздушной смеси при пониженных давлениях
Введение
различия между результатами численных расчетов простейших электронных характеристик и данными эксперимента могут различаться более чем в 2-3 раза [5]. При подходе к решению задач физики плазмы зачастую критическим является вопрос выбора способов расчета. Так, в классическом тлеющем разряде низкого давления изменение физических величин происходит на различных пространственных масштабах. Характеристики большей части плазмы меняются на расстояниях порядка Ь длины разрядного промежутка. Исключение составляют слои пространственного заряда, примерное местоположение которых в рассматриваемой системе нетрудно указать заранее. При численном моделировании тлеющих разрядов зачастую используется метод конечных элементов. Шаг разрешения вычислительной сетки должен быть малым в областях основных неоднородностей (слоев и т.д.) и большим в области плазмы. Наличие предварительных данных о структуре разряда позволяет сильно сэкономить время расчета.
Разряды атмосферного давления характеризуются стримерами, нитями, искрами и другими динамически развивающимися пространственными структурами, характеризующимися резкими пространственными градиентами и быстро протекающими плазмохимическими процессами. Обычно стримеры и искры развиваются в сильно неоднородных электрических полях около заостренных электродов. Они распространяются в виде узких каналов с плотностью плазмы в них, варьирующейся в диапазоне концентраций частиц 1019-1021 м . Электрическое поле из-за высокой проводимости плазмы принимает большие значения в области головки стримера, где сильная ионизация электронными столкновениями приводит к его дальнейшему развитию. По своей природе стример является примером сильно неоднородной, неравновесной, нестабильной плазмы, требующей многоуровневого анализа.
Основой описания явления пробоя газового промежутка является классическая теория Таунсенда, берущая за основу лавинный механизм образования электронов в результате прямой ионизации. Когда напряженность электрического поля такова, что позволяет электрону на длине свободного пробега набрать энергию, достаточную для последующей ионизации, происходит пробой, выражающийся в резком экспоненциальном росте концентрации электронов и плотности тока. Следствием этого является закон Пашена, выражающий зависимость пробойного напряжения от произведения давления на длину межэлектродного промежутка и согласующийся с соответственными экспериментально полученными кривыми Пашена. В случае, когда источником вторичных электронов является электронная эмиссия ионами с поверхности катода, реализуется таунсендовский или диффузионный пробой, критерием которого является выражение
у(еаа -1) = 1,
где у - коэффициент вторичной электронной эмиссии; а - ионизационный коэффициент Таунсенда; ё- межэлектродное расстояние.
Когда число электронов в лавине достаточно велико, в ней происходит разделение заряда, экранирующее внешнее поле, и образуется стример, движущийся со скоростью, превышающей скорость дрейфа электронов. Источником вторичных электронов являются фотоионизация, выталкивание электронов на фронте головной части лавины и убегание быстрых электронов. Критерий лавинно-стримерного перехода весьма условен и основывается на экранировке внешнего поля. Ситуация еще больше усложняется в случае импульсных разрядов ввиду их разнообразия. Возможность применения законов подобия для исследования импульсных разрядов обсуждается в [6]. Особое значение при этом приобретает пространственное распределение заряженных частиц и потенциала в различные моменты времени. Стримерный пробой реализуется при высоких давлениях, когда концентрации заряженных частиц велики. Однако в случае малых межэлектродных промежутков лавинно-стримерный переход может не осуществиться. Как известно, напряжение пробоя межэлектродного промежутка описывается кривыми Пашена, дающими зависимость его значения от параметра рё - произведения давления на расстояние между электродами (рис. 1).
В районе минимума кривой Пашена осуществляется стабильное горение тлеющего разряда. При атмосферных давлениях в гелии это будет соответствовать расстояниям между электродами порядка 400 мкм. А так как механизм пробоя тлеющего разряда является лавинным по своей природе, то в этих условиях представляется возможным выявить и исследовать лавинно-стримерный переход. Из-за малых межэлектродных расстояний экспериментальное исследование динамики пробоя в таких системах затруднено и может дать лишь ограниченную информацию.
Целью настоящей работы является исследование динамики пробоя в гелии при атмосферном давлении с использованием численных методов. Пробой рассматривался в системе игольчатых электродов с малым межэлектродным расстоянием. В связи с тем, что геометрия исследуемого разрядного промежутка и сама газоразрядная среда в области существования движущегося в пространстве стримера оказываются сильно неоднородными, корректное численное моделирование изучаемых плазменных процессов требует развития специальных методов, обеспечивающих необходимую точность расчетов, с одной стороны, и разумную экономию вычислительных ресурсов - с другой.
106
105
03
(U
s
и
и 104
&
Л
я
103
102
КГ1 10° 101 102 103
Произведение давления на межэлектродное расстояние, торрсм
—Не - N6 - Аг Н2 — N2
Рис. 1. Вид кривых Пашена для различных газов
Описание задачи моделирования пробоя с использованием гидродинамической модели плазмы
В ходе работы было проведено численное моделирование импульсного пробоя в гелии при атмосферном давлении в системе электродов «игла-игла» с использованием гидродинамической модели, широко использующейся при описании процессов в плазме стационарных разрядов [7-13]. Она состоит из уравнений баланса электронов и ионов и уравнения Пуассона на электростатический потенциал: Дф = 4пе(пе — щ)
где I - источник заряженных частиц за счет прямой ионизации; Я - гибель заряженных частиц в результате объемной рекомбинации; пе - концентрация электронов; П; - концентрация ионов; ф - потенциал; Е, - напряженность электрического поля, равная градиенту потенциала — Vф; - подвижность элек-
тронов , определенная как
Че mevm
(qe - заряд электрона, те - масса электрона, vm - частота столкновений
электронов с атомами газа, равная 6,1109(сторр)-1Р, Р - давление газа в торрах); £>е - коэффициент диффузии электронов, определяемый через соотношение Эйнштейна как Ие^^е, где постоянная Больцмана кБ равна 8,6210-5 эВК-1, Те - температура электронов, взятая равной 3 эВ. Ана-
логично, подвижность ионов цг и коэффициент диффузии ионов £>е определялись как
41
(Mi - масса
иона гелия) и соответственно (Т\ = 0,025 эВ - температура ионов, равная температуре газа).
Задача решалась в осесимметричной постановке. Геометрия задачи представлена на рис. 2: рассматривается система игольчатых электродов, находящихся друг от друга на расстоянии 600 мкм. Пробой в системе осуществлялся при подаче короткого (около 40 нс) импульса напряжения, вид которого представлен на рис. 3. Ввиду этого главную роль в динамике пробоя играли быстропротекающие процессы прямой ионизации гелия, скорость которой задавалась при помощи коэффициента ионизации Таунсенда а [11] в зависимости от локального значения напряженности электрического поля:
( вр '
I = аи8 = и8 А ехр I —:
R
где А = 3 (смторр)-1, В = 34 В/(смторр); Р - давление, равное атмосферному; Е - абсолютная величина напряженности электрического поля [12]. Рекомбинационные процессы в этих условиях вступают в силу только при больших концентрациях заряженных частиц. Скорость этих процессов описывалась в пороговом приближении:
0, пе <1014см"3 И01Опе,пе >1014 см"3 .
В качестве граничных условий на открытых границах использовались условия Неймана на концентрации частиц и потенциал. Правый электрод заземлен, на левый подавался импульс напряжения длительностью 40 нс (рис. 3). Условия на потоки ионов и электронов задавались в соответствии с [13]. Данная модель реализована в плазменном модуле программного пакета ветБ, анализ возможностей которого и сопоставление с другими пакетами моделирования плазмы проведено в [14].
4
Рис. 2. Геометрия и изначальная сетка: 1, 2, 3 - открытые границы; 4 - ось симметрии; 5, 6 - электроды
5000 4500
4000
3500
3000
2500
PQ
<D
S
и
H
¡^ 2000 К 1500 1000 500
2
3
1
0 0,5 1 1,5 2 2,5 3 3,5 4 Время /, с х104
Рис. 3. Профиль импульса напряжения Особенности численного моделирования неоднородной плазмы в условиях электрического пробоя
Использование метода конечных элементов с фиксированной вычислительной сеткой является нецелесообразным, так как простые оценки пространственных и временных масштабов изменения физических величин в стримерах [11, 15] показывают, что при рассмотрении развития стримера в атмосферном воздухе на расстоянии 1 см потребовало бы вычислительной сетки, состоящей из 109 ячеек.
При моделировании явлений газового пробоя представляется целесообразным применение метода конечных элементов с использованием адаптивных алгоритмов построения вычислительной сетки [16]. Такой подход позволяет автоматически генерировать вычислительную сетку с минимальным вмешательством со стороны пользователя и динамически подстраивать ее под электростатический потенциал и потоки заряженных частиц. Первоначально адаптивные методы построения вычислительных сеток разрабатывались для решения резко нестационарных задач гидродинамики и аэродинамики с решением уравнений Эйлера или Навье-Стокса. Гидродинамическое описание потоков заряженных частиц в плазме позволило распространить эти методы на задачи низкотемпературной плазмы и газовых разрядов.
Суть метода заключается в следующем. Вся вычислительная область представляет собой ячейку нулевого уровня (квадрат), которая затем может улучшаться разбиением на 4 ячейки первого уровня и т.д., образуя древовидную структуру с корнями и листьями. Разбиение ячеек производится в зависимости от получающихся в ходе решения градиентов физических величин, которые определяются как разность их значений в соседних ячейках. Затем локальные значения этих градиентов используются для подсчета коэффициента адаптации а. Затем, если его значение превышает 1, то обе соседние ячейки разбиваются на 4. В случае, когда значение а лежит в диапазоне от % до 1, то ячейки остаются неизменными. Если же значение а меньше % для всех ячеек, имеющих общий корень, то ячейки объединяются в одну. Размеры соседних ячеек могут отличаться максимум в два раза, что обеспечивает плавные градиенты, значительно упрощает процедуру интерполяции и ускоряет расчеты. Отметим, что выбор вида зависимости коэффи-
циента адаптации от градиентов искомых величин зависит от пользователя и зачастую является ключевым фактором, влияющим на точность и скорость решения.
На рис. 4 представлен пример построения адаптивной сетки и направление расчета градиентов. Вычисление неизвестных величин происходит исключительно в листьях древовидной структуры уровней.
Уровень 3
Рис. 4. Пример адаптивного построения вычислительной сетки: направление расчета градиентов величин при вычислении коэффициента адаптации (а); древовидная структура уровней и направление движения
по сетке при проведении вычислений (б)
Адаптация вычислительной сетки во время решения происходила в зависимости от параметра а [12], заданного как
а = 20Г ^ ) + "( ) 1 + 3 (к>&0(ие) + ) + ) + Ь^Щ)) ,
^ max(ne) max(иe))
Каждая отдельная ячейка проходила процедуру улучшения тогда, когда градиент а превосходил фиксированное значение Vа >1 в этой ячейке. Загрубление сетки происходило при значениях Va <1/4 в этой ячейке. Отметим, что данный вид зависимости коэффициента адаптации от параметров плазмы был получен авторами [12] в ходе многочисленных пробных расчетов.
Результаты и анализ
В ходе моделирования была получена детальная динамика пробоя, начавшемся при напряжениях порядка 1,5 кВ.
На рис. 5 представлены напряженности электрического поля и вычислительные сетки в разные моменты пробоя. Конфигурация сетки изменялась в соответствии с вышеописанным алгоритмом и обеспечивала требуемую точность вычислений.
На рис. 6 представлено пространственное распределение концентрации электронов в разные моменты пробоя в разрядном объеме и на оси симметрии. Как видно, пробой имеет стримерный характер, динамика которого может быть описана следующим образом. В предпробойной стадии процессы ионизации преобладают в областях у электродов с наибольшим значением напряженности поля. Дрейф электронов от катода к аноду объясняет образование локального максимума концентрации электронов и малые значения у анода. В начале пробоя у анода образуется слой положительного заряда, который постепенно смещается к катоду за счет дрейфа ионов и диффузии электронов, образуя положительный стример, экранирующий внешнее поле. Разбиение в области около катода становится меньше для лучшего расчета образующегося градиента напряженности электрического поля. Экранирование поля стримером, в свою очередь, приводит к увеличению напряженности поля у катода, интенсификации дрейфа электронов и образованию анодонаправленного стримера. Слияние двух стримеров (рис. 5, г) приводит к образованию локального максимума концентрации электронов. Адаптация сетки в левой части модели происходит в соответствии с образующимися градиентами концентраций в слое у катода и на границе анодонаправ-ленного стримера.
Е, кВ/см
I
14
12 10 8
Рис. 5. Контурные графики распределения напряженности электрического поля и изменение вычислительной сетки в процессе пробоя в моменты времени: 7,8 нс (а); 9 нс (б); 9,4 нс (в); 9,8 нс (г)
б
пе, см
I
101
101
101
101
101
109
108
10'
106
105
Рис. 6. Распределение концентрации электронов в моменты пробоя: 7,4 нс (а); 9 нс (б); 9,4 нс (в); 9,8 нс (г)
6
4
2
0
3
а
в
г
Заключение
Проведены расчеты динамики пробоя в гелии при атмосферном давлении с использованием адаптивных алгоритмов. Продемонстрирован в действии процесс адаптации сетки под градиенты получающихся в ходе решения параметров плазмы - напряженности электрического поля и концентрации заряженных частиц. Показано формирование и распространение стримеров в процессе пробоя. Относительно простая система уравнений позволила описать процесс пробоя и формирование стримера. Динамика пробоя представляла собой образование катодо- и анодонаправленных стримеров, движущихся навстречу друг другу, с характерным экранированием поля пространственным зарядом, образующимся в головках стримеров. Следует отметить, что полученная динамика пробоя в системе не учитывала процессы фотоионизации, которые зачастую считаются основными в формировании стримеров.
Результаты моделирования могут являться основой для дальнейшего исследования явлений электрического пробоя путем учета сложных плазмохимических процессов в многокомпонентных газах.
Литература
1. Bogdanov E.A., Kolobov V.I., Kudryavtsev A.A., Tsendin L.D. Scaling laws for oxygen discharge plasmas // Technical Physics. 2002. V. 47. N 8. P. 946-954.
2. Bogdanov E.A., Kudryavtsev A.A., Tsendin L.D., Arslanbekov R.R., Kolobov V.I., Kudryavtsev V.V. Substantiation of the two-temperature kinetic model by comparing calculations within the kinetic and fluid models of the positive column plasma of a DC oxygen discharge // Technical Physics. 2003. V. 48. N 8. P. 983994.
3. Bogdanov E.A., Kudryavtsev A.A., Tsendin L.D., Arslanbekov R.R., Kolobov V.I., Kudryavtsev V.V. Scaling laws for the spatial distributions of the plasma parameters in the positive column of a DC oxygen discharge // Technical Physics. 2003. V. 48. N 9. P. 1151-1158.
4. Gutsev S.A., Kudryavtsev A.A., Zamchiy R.Yu., Demidov V.I., Kolobov V.I. Diagnostics and modeling of a short (without positive column) glow discharge in helium with nonlocal plasma // Proc. 40th European Physical Society Conference on Plasma Physics. Finland, 2013. N 06.502.
5. Чернышева М.В., Марек В.П., Чирцов А.С., Швагер Д.А. Компьютерное моделирование при изучении физических процессов в тлеющем разряде в воздушных смесях при низких давлениях // Научно-технический вестник информационных технологий, механики и оптики. 2014. № 3 (91). С. 140-146.
6. Mesyats G.A. Similarity laws for pulsed gas discharges // Physics-Uspekhi. 2006. V. 49. N 10. P. 10451065.
7. Bogdanov E.A., Chirtsov A.S., Kudryavtsev A.A. Fluxes of charged particles in two-chamber ICP discharge in oxygen // IEEE Transactions on Plasma Science. 2011. V. 39. N 11 part 1. P. 2562-2563.
8. Chirtsov A.S., Kapustin K.D., Kudryavtsev A.A., Bogdanov E.A. Nonlocal behavior of electron fluxes and excitation rates for «local» EEDF in moderate and high pressures DC positive column plasmas // IEEE Transactions on Plasma Science. 2011. V. 39. N 11 part 1. P. 2580-2581.
9. Rafatov I. Bogdanov E.A., Kudryavtsev A.A. On the accuracy and reliability of different fluid models of the direct current glow discharge // Physics of Plasmas. 2012. V. 19. N 3. Art. 033502.
10. Bogdanov E.A., Chirtsov A.S., Kudryavtsev A.A. Fundamental nonambipolarity of electron fluxes in 2D plasmas // Physical Reviews Letters. 2011. V. 106. N 19. Art. 195001.
11. Райзер Ю.П. Физика газового разряда. 3-е изд. Долгопрудный: Интеллект, 2009. 736 с.
12. Kolobov V.I., Aslanbekov R.R. Simulations of low-temperature plasmas with adaptive cartesian mesh // ASP Conference Series. 2012. V. 459. P. 328-333.
13. Hagelaar G.J.M., Kroesen G.M.W., Van Slooten U., Schreuders H. Modeling of the microdischarges in plasma addressed liquid crystal displays // Journal of Applied Physics. 2000. V. 88. N 5. P. 2252-2262.
14. Kolobov V.I., Arslanbekov R.R., E Bogdanov.A., Eliseev S., Kudryavtsev A.A. Comparison of computational tools for simulations of glow and corona discharges // Proc. XXXI International Conference on Phenomena in Ionized Gases (ICPIG). Granada, Spain, 2013. PS2-024.
15. Базелян Э.П., Райзер Ю.П. Искровой разряд. М.: МФТИ, 1997. 320 с.
16. Kolobov V.I., Arslanbekov R.R. Towards adaptive kinetic-fluid simulations of weakly ionized plasmas // Journal of Computational Physics. 2012. V. 231. N 3. P. 839-869.
Елисеев Степан Иванович - инженер-исследователь, СПбГУ, 198504, Санкт-Петербург, Россия,
step.eliseev@yandex.ru
Демидов Владимир Иванович - кандидат физико-математических наук, доцент, Research Professor,
Университет Западной Вирджинии, WV 26506, Morgantown, West Virginia
Чирцов Александр Сергеевич - кандидат физико-математических наук, доцент, профессор,
Университет ИТМО, 197101, Санкт-Петербург, Россия; доцент, СПбГУ, 198504, Санкт-Петербург, Россия, Alex_chirtsov@mail.ru
Кудрявцев Анатолий Анатольевич
Колобов Владимир Иванович
Богданов Евгений Анатольевич
Vladimir I. Demidov Alexander S. Chirtsov
Anatoly A Kudryavtsev
Vladimir I. Kolobov Evgeny A Bogdanov
кандидат физико-математических наук, доцент, старший научный сотрудник, Университет ИТМО, 197101, Санкт-Петербург, Россия; доцент, СПбГУ, 198504, Санкт-Петербург, Россия, akud53@gmail.ru кандидат физико-математических наук, технический сотрудник, CFD Research Corporation, AL 35806, Хантсвиль, Алабама, США, vik@cfdrc.com
кандидат физико-математических наук, доцент, доцент, СПбГУ, 198504, Санкт-Петербург, Россия, eugene72@mail.ru
PhD, Associate professor, Research Professor, West Virginia University, WV26506, Morgantown, West Virginia, USA
PhD, Associate professor, ITMO University, 197101, Saint Petersburg, Russia; Associate professor, Saint Petersburg State University, 198504, Saint Petersburg, Russia, Alex_chirtsov@mail.ru PhD, Associate professor, ITMO University, 197101, Saint Petersburg, Russia; Associate professor, Saint Petersburg State University, 198504, Saint Petersburg, Russia, akud53@gmail.ru
PhD, technical fellow, CFD Research Corporation, AL 35806, Huntsville, Alabama, USA
PhD, Associate professor, Associate professor, Saint Petersburg State University, 198504, Saint Petersburg, Russia, eugene72@mail.ru
Принято к печати 02.07.14 Accepted 02.07.14