¿J МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ, ЧИСЛЕННЫЕ МЕТОДЫ И КОМПЛЕКСЫ ПРОГРАММ (в строительстве)
УДК 537.525
Сафиуллин Р.К. - доктор физико-математических наук, профессор E-mail: rksaf@mail.ru Салаватуллин A.A. - аспирант E-mail: inurrick@gmail.com
Казанский государственный архитектурно-строительный университет
Адрес организации: 420043, Россия, г. Казань, ул. Зеленая, д. 1
Компьютерное моделирование и расчеты характеристик тлеющего разряда в потоке газа
Аннотация
Тлеющий разряд в потоке газов (ТРП) широко используется для создания активных сред молекулярных лазеров. В мощных проточных ССЬ-лазерах применяется поперечный (относительно направления потока газа) тлеющий разряд. Для повышения устойчивости тлеющего разряда в этом случае применяют сплошной анод и секционированный (в форме набора штырей или ножей) катод. Представляет значительный интерес информация о пространственном распределении заряженных частиц (свободных электронов, положительных и отрицательных ионов) и электрического поля внутри таких разрядных камер. Для решения системы из четырех эллиптических уравнений (для концентраций электронов, положительных и отрицательных ионов и потенциала электрического поля) в работе использован итерационный метод переменных направлений. В результате получены двумерные распределения указанных величин внутри разрядных камер такого типа.
Ключевые слова: тлеющий разряд, плотность числа частиц, электрический потенциал, пространственное распределение.
Математическая модель зарядовой кинетики в газоразрядной плазме в электроотрицательном газе
Теоретическая модель ТРП должна включать уравнения неразрывности для заряженных частиц, уравнение Пуассона для потенциала электрического поля, закон Ома в дифференциальной форме, уравнение для функции распределения электронов по энергии (ФРЭЭ) и уравнения газодинамики [1, 2].
Принятая в работе система уравнений имеет вид:
<ИуГе = ((%- ао)пе + ас,пп - Р^пепр+ д, (1)
сИУГр = (%пе - релепр - Рцпрп„ + д, (2)
ШуГ„ = аапе - ас,п„ - (5цпрп„ , (3)
А<р = -(пр-пе- п„)е/е0, (4)
Ге = пеи - Пц/ИцЕ - Д Упе, Гр = при + пр/ИрЕ - Д РНр, Г„ = л „и - п^Е - Д, Е - У(р, (5) у = е(Гр-Ге- Г'г) = е^рЦр+ПеЦе+ПпЦ^Е + е(пр-пе-п1г)и +е(Ое Упе + Д, Уп„- I),, Упр). (6)
Здесь Ге, Гр, Г„ - векторы плотности потоков электронов, положительных и отрицательных ионов, соответственно; а, , аа , а,/ - частоты ионизации молекул газа электронным ударом, прилипания электронов к электроотрицательным молекулам и отлипания электронов от отрицательных ионов; Д„ Д, - коэффициенты электрон-ионной и ион-ионной рекомбинации; <:/ - интенсивность внешнего источника ионизации; пс, пр, п„ - концентрации электронов, положительных и отрицательных ионов; у - вектор плотности тока; Д, /),„ Д, /4. ¡ир, //„ - коэффициенты диффузии и подвижности соответствующих заряженных частиц; и - вектор скорости потока газа в данной точке разрядной камеры (РК); Е - вектор напряженности электрического поля; £,, -электрическая постоянная; V и А - операторы градиента и Лапласа, соответственно; е -абсолютная величина заряда электрона.
Иавесшш КГАСУ, 2014, № 2 (28)
Математическое мо9епароваяие, численные мето9ы а комплексы программ (в строшпельстве)
В условиях квазинейтральности плазмы тлеющего разряда первый член правой части выражения (6) является доминирующим, поэтому, с учетом того, что //е » ¡ир, ¡и„. получаем у ~ епср1сЕ. Из уравнений (1-3) следует соотношение ¿йту = 0. Иногда для расчета электрического потенциала его применяют к приближенному выражению у ~ епе/иеЕ [3, 4]. Однако, как будет видно из дальнейшего, более точно распределение электрического потенциала внутри РК можно рассчитать по уравнению Пуассона (4).
Нами была рассмотрена РК с поперечным ТРИ, через которую газ прокачивался в направлении оси X, а электроды располагались, как показано на рис. 1. Использовался сплошной анод и секционированный катод. Катодная плата представляла собой набор узких пластинчатых (ножевых) катодов, каждый из которых подключался к источнику питания через балластное сопротивление. Как известно, секционирование катода применяется для повышения устойчивости тлеющего разряда и увеличения энерговклада в разряд [5, 6].
Поток raía
ЧК
сатоды
Ж
\\Ч
X анод
b
I
Рис. 1. Схема РК. Межкатодное расстояние а = 4 см, Ъ = 3 см, с = 8 см, L = 32 см, длина катодов / = 3,6 см
При переходе в приближение амбиполярной диффузии, которое часто используется для описания процессов в газоразрядной плазме [5-7], уравнения (1-3) для рассматриваемой РК приобретают следующий вид:
//.Г/7, Гх + iivdne/dy = Дм 0njdx2 + (9njdy~) +(а,-аа)пе+ас1п„-f3einenp+ q, (7)
и.Гпп Гх + nvdn/dy = Dpa ((fnp/dx2 + cfn/dy2) +а,пе- f3einenp- Рцпрп„ + q, (8)
ихдп„/дх + Uydn/dy = D„a (tfn^dx2 + tfn/dr) + aane- adn„- Рцпрп№. (9)
Электрическое поле входит в уравнения (7-9) неявным образом через величины a¡, аа, Д, и //е. Наиболее сильно от параметра E/N зависит величина а, . Здесь введены следующие обозначения: Dea, Dpa, D„a - коэффициенты амбиполярной диффузии для электронов, положительных и отрицательных ионов, соответственно. Общие выражения для них приведены в [7, 8].
Аналогичным образом могут быть записаны уравнения для энергии (энтальпии или температуры) газа, а также уравнения для колебательных температур в случае С02 - лазеров или для населенностей колебательных уровней молекул СО и N2в случае СО-лазеров.
Следует отметить, что подобная РК рассматривалась в работах [3, 4], в которых
считалось, что плазма всюду квазинейтральна (пр
п„). т.е. рассматривался только
положительный столба разряда. Более того, в этих работах использовалось точное условие нейтральности плазмы, так как одно из уравнений зарядовой кинетики, а именно, уравнение для концентрации положительных ионов пр из рассмотрения исключалось. В нашей же модели решаются все уравнения зарядовой кинетики, что позволяет оценить степень отклонения от квазинейтральности как в приэлектродных областях, так и в положительном столбе разряда.
В данной работе мы ограничились рассмотрением системы уравнений (4), (7-9). Течение газа считалось ламинарным. Профиль скорости газа на входе и внутри РК задавался пуазейлевским или полагалось для простоты их = const., uv = 0.
Метод решения системы уравнений зарядовой кинетики и электрического поля
Система эллиптических уравнений (4, 7-9) решалась итерационным методом переменных направлений [9] в плоскости ХУ. Расчетная область представляла собой
Иавесшш КГАСУ, 2014, № 2 (28)
Математическое мо9епароваяие, численные мето9ы а комплексы программ (в строшпельстве)
прямоугольник со стороной d = 8 см в направлении оси X и стороной 6 = 3 см, параллельной оси Y. Она содержала один или несколько пластинчатых катодов толщиной 1 мм, расположенных на верхней стороне прямоугольника, и часть сплошного анода (нижняя сторона прямоугольника). Эта область покрывалась равномерной разностной сеткой 100 х 100. Граничные условия для задачи выбирались такими же, как в работах [3, 4]. Помимо этого граничные условия варьировались в широких пределах с целью выяснить зависимость от них характера распределения концентраций заряженных частиц и электрического потенциала. Итерации применялись к каждому из уравнений (4,7-9) в отдельности и к системе в целом, пока не достигалась относительная погрешность вычислений в пределах 1% . Следует отметить, что для сравнения с результатами работы [4] распределение потенциала рассчитывалось нами не только согласно уравнению (4), но также из уравнения div(neW (р) = 0 , примененного в [4].
Результаты расчетов распределения электрического потенциала внутри РК
Ниже приведены некоторые результаты ряда расчетов. Они получены при задании граничных условий Дирихле для искомых величин на всей границе рассматриваемой прямоугольной области, причем везде на границе использовалось условие квазинейтральности плазмы. Это означает, что рассматривался положительный столб разряда, который занимает практически весь объем камеры (толщина прикатодного и прианодного слоев в таких РК не превышает 1 мм [5]). На входе (х = 0) и выходе из расчетной области (х = d = 8 см) было принято: пс{у) = 106-108 см°, пр{у) = 1,1 пе(у), п„{у) = 0,1 пе(у), (р(у) = - 1400 у/Ъ (от 0 на аноде
10 3 10 3
до - 1400 В при у = Ь). На аноде (у = 0) полагалось пе(х) = 10 см°, пр{х)= 1,5 10 см°, п„{х) = 0,5 101" см"3, <р(х) = 0. На катодной плате (у = Ъ, рис. 2-4) было положено:
пе(х) = пе0( 1 - //1 ún[n(2xk/d - 1)] | ), пр(х) = лр0( 1 - //1 ún[n(2xk/d - 1)] | ), (10)
п„(х) = пп0(1 - //1 sm\n(2xkcl - 1)] | ), <р(х) = - 1400 + 500 | sm\n(2xk с/ - 1)] | , (11) где к полагалось равным 1, 2 и 3, а параметр // - варьировался в диапазоне от 0,1 до 0,9. Такой вид граничных условий для концентраций заряженных частиц был выбран нами в соответствии с результатами экспериментов, описанными в [10], где было показано, что пе(х) вдоль такой РК имеет ярко выраженный колебательный характер.
При заданных граничных условиях программа позволяет рассчитывать двумерные распределения потенциала и плотностей заряженных частиц внутри РК. Наблюдается качественное и полуколичественное согласие с экспериментом, выполненным для такой РК в работе [10].
Проведенные расчеты показывают сильную зависимость полученных распределений от характера граничных условий. Расчеты показывают также, что в отличие от результатов одномерного анализа [6] подобных РК, поле Е максимально вблизи катода и убывает в направлении к аноду. Такая же картина была получена ранее в экспериментах [10]. Это свидетельствует о том, что ввиду сильной неоднородности электрического поля в подобных РК одномерный подход здесь является, по-видимому, неадекватным.
Расчеты проводились для скоростей потока в диапазоне (30-250) м/с. Напряжение между анодом и катодом задавалось, как и в эксперименте, равным 1400 В при статическом давлении газа р в диапазоне (4-10) кПа и температуре Т = 293 К. Оказалось, что характер решения сильно зависит не только от вида граничных условий, но и от значений коэффициентов амбиполярной диффузии Dea, Dpa, D„a. Согласно проведенным расчетам, степень отклонения от квазинейтральности в пределах положительного столба достигала 3-4%, а в приэлектродных областях сильно зависела от характера граничных условий и могла быть на порядок больше. В качестве примера на рис. 2-4 приведены рассчитанные распределения электрического потенциала внутри РК. Они соответствуют граничным условиям для искомых величин (10, 11). Видно, что напряженность электрического поля вблизи катода в несколько раз больше, чем в окрестности анода.
Известил КГАСУ, 2014, № 2 (28)
Математическое мо9епироваяие, численные мето9ы и комплексы программ (в строительстве)
Катоды Катоды
Катоды Катод
10-1.26Е+03 9-1.12Е+03 8 -9 79Е+02 7 -8 39Е*02 6 -6.99Е*02 5-5 60Е+02 4 Ч20Е+02 3 -2.80Е+02 2-1.40Е+02 1 О.ООЕ'ОО fi
Анод
Рис. 3. Электрический потенциал внутри РК (к = 2.11 = 0,8)
10-1 2БЕ-03 9-1.12Е+03 8 -9 80Е+02 7 -8 40Е+02 6-7.00Е+02 5-5 60Е+02 4-4.20Е+02 3-2 80Е+02 2-1 40Е-02 1 0.00Е+00 fi
Анод
Рис. 2. Электрический потенциал внутри РК (k= l,ri = 0,8)
10 -1.25Е+03 9 -1.11Е+03 8 -9.79Е+02 7 -8.39Е-»02 6 -6 99Е+02 5-5.59Е+02 4-419Е*02 3 -2.79Е+02 2 -1.39Е*02 1 О.ООЕ+ОО
fi
Анод
Рис. 4. Электрический потенциал внутри РК (к = 3, ц = 0,8)
Ю-1.26Е<03 9-1.12Е-»03 8-9.80Е<02 7-8.40Е*02 6 -7.00Е+02 5-5.60Е+02 4 -4.20Е+02 3 -2.80Е<02
1 О ООЕ+ОО
fi
Анод
Рис. 5. Электрический потенциал внутри РК
Для сравнения на рис. 5 приведено рассчитанное распределение потенциала, соответствующее другим граничным условиям на боковых границах и на катодной плате. Здесь было принято, что при х = 0 и х = d потенциал изменяется как ф(у) = - 900 у/b (от 0 на аноде до - 900 В при у = Ь). На катодной плате для концентраций заряженных частиц
9 3 10 3 9 3
было принято: пе(х) = 2-10 см"3, пр{х)= 10 см°, п„(х) = 8-10 см°, <р(х) = изменялось линейно от - 900 В при х = 0их=8смдо- 1400 В на пластинчатом катоде.
Таким образом, методом переменных направлений численно решена двумерная система уравнений зарядовой кинетики и уравнения для электрического потенциала внутри РК с поперечным ТРП. При этом значения коэффициентов ионизации и прилипания электронов к молекулам вычислялись в результате расчета ФРЭЭ методом, описанным в работах [11, 12]. Достигнуто качественное и полуколичественное согласие с имеющимся экспериментом. Численный метод оказался устойчивым при варьировании коэффициентов уравнения (7-9) в пределах: 0 < а, < 1,5- 10 с" , аа < 1,2-104с"1, q < 5-Ю20 мдс~\ Данный метод может быть рекомендован для решения более полной системы уравнений, включащей уравнения для газовой и колебательных температур или населенностей колебательных уровней молекул.
Список библиографических ссылок
1. Сафиуллин Р.К. Распределение концентраций заряженных частиц и потенциала электрического поля в тлеющем разряде в потоке газа. // Известия вузов (Проблемы энергетики), 2002, № 1-2. - С. 69-77.
2. Сафиуллин Р.К. К расчету камер с поперечным тлеющим разрядом в потоке
Навестим КГАСУ, 2014, № 2 (28)
Математическое мо9епироваяие, численные мето9ы и комплексы программ (в строительстве)
электроотрицательного газа. // Известия вузов (Проблемы энергетики), 2003, № 3-4. -С. 140-145.
3. Басыров Р.Ш., Гайсин Ф.М., Миннигулов A.M., Тимеркаев Б.А. Пространственное раиределение параметров тлеющего разряда в потоке электроотрицательного газа. // ТВТ, 1994, 32, № 2. - С. 334-338.
4. Басыров Р.Ш., Тимеркаев Б.А. Модель тлеющего разряда в поперечном потоке электроотрицательного газа. //ТВТ, 28, № 1. - С. 30-34.
5. Голубев B.C., Пашкин С.В. Тлеющий разряд повышенного давления. - М.: Наука, 1990.- 336 с.
6. Велихов Е.П., Ковалев А.С., Рахимов А.Т. Физические явления в газоразрядной плазме. - М.: Наука, 1987. - 160 с.
7. Rogoff G.L. Ambipolar Diffusion Coefficients for Attaching Gases. // J. Phys.D. Appl. Phys., 1985, 18, № 8. - P. 1533-1545.
8. Сафиуллин P.К. Амбиполярная диффузия в газоразрядной плазме. // Известия вузов (Проблемы энергетики), 2000, № 5-6. - С. 110-112.
9. Самарский А.А. Теория разностных схем. - М.: Наука, 1977. - 656 с.
10. Миннигулов A.M. Характеристики электронного газа и распределение потенциала электрического поля в тлеющем разряде в потоке воздуха. Дисс. на соискание ученой степени кандидата физико-математических наук. - Казань, КАИ, 1980. - 186 с.
11. Арасланов Ш.Ф., Сафиуллин Р.К. Энергетическое распределение электронов в плазме тлеющего разряда. // Известия вузов (Проблемы энергетики), 1999, № 7-8. - С. 61-68.
12. Сафиуллин Р.К. Расчет констант скоростей ионизации и диссоциативного прилипания электронов к молекулам в газоразрядной плазме. // Известия вузов (Проблемы энергетики), 2001, № 7-8. - С. 55-63.
Safiullin R.K. - doctor of physical and mathematical sciences, professor E-mail: rksaf@mail.ru Salavatullin A.A. - post-graduate student E-mail: inurrick@gmail.com
Kazan State University of Architecture and Engineering
The organization address: 420043, Russia, Kazan, Zelenaya St., 1
Computer simulation and calculations of the characteristics of a glow discharge in a gas flow
Resume
Glow discharge in a gas flow is widely used for the creation of active mediums of powerful molecular gas lasers. In the paper transverse glow discharge in electronegative gas flow was investigated numerically by alternating direction method. Discharge chambers with continuous anode and partitioned cathodes were considered. As a result, two-dimensional distributions of charged particles densities (for electrons, positive and negative ions) and of electric potential inside such chambers were obtained.
A system of elliptic equations for charged particles (electrons, positive and negative ions) and for electric potential was solved by iterative alternating direction method in two-dimensional subspace XY. Calculations were carried out for a rectangle which contained 1 mm width lamellate cathodes and continuous anode. The difference grid (100x100) was applied. The boundary conditions for densities of electrons, positive and negative ions and for electric field potential were widely varied with a view to evaluate their influence on the solution.
Under given boundary conditions the program elaborated enables to calculate two-dimensional distributions of electric potential and of charged particles densities inside discharge chambers. The calculations showed the strong dependence of these distributions on boundary conditions. Keywords: glow discharge, number density, electric potential, spatial distribution.
HaBecmuq EffACY, 2014, № 2 (28)
MameMamtwecKoe MoQentipoBaiiiie, <iticnefifibie MemoQbi a KOMnneKCbi npospaMM (b cmpoumenbcmBe)
Reference list
1. Safiullin R.K. Spatial distributions of number densities of charged particles and of electric potential in glow discharge in a gas flow // Izvestya vusov. Problemy Energetiky, 2002, № 1-2.-P. 69-77.
2. Safiullin R.K. Numerical investigation of chambers with transverse glow discharge in a flow of attaching gas // Izvestya vusov, Problemy of Energetiky, 2003, № 3-4. - P. 140-145.
3. Basyrov R.Sh., Gaisin F.M., Minnigulov A.M., Timerkaev B.A. Space distribution of glow discharge parameters in attaching gas flow // Physics of high temperatures (Russia), 1994, Vol. 32, № 2. - P. 334-338.
4. Basyrov R.Sh., Timerkaev B.A. A model of glow discharge under transverse flow of attaching gas // Physics of high temperatures (Russia), 1994, Vol. 28, № 1. - P. 30-34.
5. Golubev V.S., Pashkin S.V. Glow discharge at increased pressure. - M.: Publishers Nauka, 1990. - 336 p.
6. Velikhov E.P., Kovalev A.S., Rakhimov A.T. Physical phenomena in gas discharge plasma. - M.: Publishers Nauka, 1987. - 160 p.
7. Rogoff G.L. Ambipolar diffusion coefficients for attaching gases. // J. Phys. D. Appl. Phys., 1985, 18, № 8. - P. 1533-1545.
8. Safiullin R.K. Ambipolar diffusion in gas discharge plasma // Izvestya vusov, Problemy Energetiky, 2000, № 5-6. - P. 110-112.
9. Samarsky A.A. The theory of difference schemes. - M.: Publishers Nauka, 1977. - 656 p.
10. Minnigulov A.M. Properties of electron gas and electric potential distribution in glow discharge in air flow, kandidat thesis. - Kazan, KAI, 1980. - 186 p.
11. Araslanov Sh., F., Safiullin R.K. Energetic distribution of electrons in glow discharge plasma // Izvestya vusov. Problemy Energetiky, 1999, № 7-8. - P. 61-68.
12. Safiullin R.K. Calculation of rates of ionization and rates of electron dissociative attachment to molecules in gas discharge plasma // Izvestya vusov. Problemy energetiky, 2001, № 7-8.-P. 55-63.