Физика плазмы, 2023, T. 49, № 4, стр. 381-390

1. ВВЕДЕНИЕ

В задачах по численному моделированию газового разряда необходимо рассчитывать кинетику электронов в самосогласованном электрическом поле. Наиболее полное решение возможно в терминах уравнения Больцмана, которое в общем случае представляет собой интегро-дифференциальное уравнение для функции распределения электронов (ФРЭ) в шестимерном фазовом пространстве $(r,p)$. Для решения уравнения Больцмана в общем случае требуются вычислительные ресурсы, намного превышающие возможности современных персональных компьютеров. По этой причине описание развития газовых разрядов ведется в терминах уравнений для моментов ФРЭ [1]. Так, эволюция концентрации электронов n_e в пространстве и во времени описывается уравнением непрерывности

где ${{\Gamma }_{e}}$ – поток электронов, R_e – источник или сток электронов. Обычно поток электронов аппроксимируется диффузионно-дрейфовым уравнением где E – вектор напряженности электрического поля, n_e – концентрация электронов, μ_e и D_e – коэффициенты подвижности и диффузии электронов. В общем случае коэффициент диффузии электронов является тензором, и в отсутствии магнитного поля имеет два компонента: продольный D_L и поперечный D_T  относительно вектора E.

Источник или сток электронов R_e в уравне-нии (1) определяется реакциями, протекающими в плазме. В случае неона это реакции ударной ионизации атомов электронами и рекомбинации положительных ионов с электронами:

здесь N – концентрация атомов, k_ion – коэффициент ионизации, k_pe – коэффициент рекомбинации, n_p – концентрация положительных атомарных или молекулярных ионов (в неоне молекулярные ионы могут образовываться в процессе ${\text{N}}{{{\text{e}}}^{ + }} + {\text{Ne}} + {\text{Ne}} \to {\text{Ne}}_{2}^{ + } + {\text{Ne}}$). Здесь мы не будем конкретизировать тип процесса рекомбинации (диссоциативная рекомбинация, тройная рекомбинация, фоторекомбинация) поскольку далее этот процесс рассматриваться не будет. Отметим лишь, что константы различных процессов рекомбинации зависят от температуры электронов.

В отсутствии рекомбинации для решения системы уравнений (1)–(2) необходимо знание транспортных коэффициентов μ_e и D_e и скорости ионизации k_ion, которые являются интегралами ФРЭ. Чаще всего при моделировании газовых разрядов используется приближение локального поля, в рамках которого полагается, что коэффициенты μ_e, D_e и k_ion являются функциями локальной приведенной напряженности поля E/N:

Данные зависимости могут быть рассчитаны теоретически методом Монте-Карло (МК), либо с помощью уравнения Больцмана, или измерены экспериментально. Однако применение в уравнениях (1)–(2) экспериментальных зависимостей сталкивается с трудностями. Так, в экспериментах невозможно прямо получить некоторые величины. Вместо частоты ионизации (k_ionN) измеряется первый коэффициент Таунсенда α_ion т.е. количество электронов рождаемых на единице длины, а вместо коэффициентов продольной D_L и поперечной диффузии D_T измеряются отношения D_L/μ_e и D_T/μ_e. Также значения дрейфовой скорости электронов, полученные в разных типах экспериментов: time-of-flight (TOF), steady-state Townsend (SST), pulsed Townsend (PT), могут существенно отличаться друг от друга [2]. Естественно возникает вопрос о том, какие экспериментальные данные можно использовать в уравнениях (1)–(2). Кроме того, в работе [3] показано, что для гелия имеющиеся экспериментальные данные в области больших значений E/N содержат значительную ошибку, обусловленную тем, что в экспериментах в сильных полях не достигается равновесия ансамбля электронов с полем. Поскольку неон часто используется в различных газоразрядных устройствах, а в области больших значений E/N, экспериментальные данные могут содержать заметную ошибку, то возникает необходимость верификации экспериментальных значений транспортных и кинетических коэффициентов электронов в неоне.

Поскольку в [3] было установлено, что в гелии диффузионно-дрейфовое приближение для концентрации электронов становится неприменимым в полях больших ≈150–200 Тд, то целесообразно определить область применимости данного приближения и для неона. Также с целью верификации нашей МК-модели рассчитаны значения дрейфовой скорости, первого коэффициента Таунсенда и отношений D_L/μ_e и D_T/μ_e для экспериментальных условий и выполнено сравнение с имеющимися в литературе данными измерений.

Расчеты транспортных и кинетических коэффициентов электронов в неоне выполнялись и ранее. Например, в работе [4] были выполнены расчеты по программе BOLSIG+ с различными наборами сечений. Эта программа часто используется для расчетов транспортных и кинетических коэффициентов электронов в различных газах, в ней решается кинетическое уравнение Больцмана в рамках двучленного приближения для ФРЭ (приближение Лоренца) [5]. Как известно в сильных электрических полях двучленное приближение нарушается [6], поэтому целесообразно сравнить результаты наших расчетов с расчетами BOLSIG+.

Целью данной работы является расчет методом МК-значений дрейфовой скорости, частоты ионизации, коэффициентов продольной и поперечной диффузии лавины электронов в неоне, которые непосредственно можно использовать в уравнениях (1)–(2).

2. ПАРАМЕТРЫ ЛАВИНЫ ЭЛЕКТРОНОВ

Подробный алгоритм программы расчета транспорта электронов в электрическом поле на основе метода Монте-Карло (МК) изложен в [7]. Здесь мы приведем лишь данные для сечений взаимодействия электронов с атомами неона. Всего учитываются 11 процессов взаимодействия электронов с атомами, которые перечислены в таблице 1. Интегральное и дифференциальное сечения упругого рассеяния в области энергий электрона меньших 50 эВ взяты из [8]; для области 50 эВ–20 кэВ сечения рассчитаны по программе ELSEPA [9]. Полное сечение ионизации в области энергий электрона меньших 40 эВ взято из [10], в диапазоне 40−1000 эВ из [11], в области 1−20 кэВ использованы данные [12]. Интегральные сечения возбуждения взяты из базы данных Biagi v7.1 [13]. Интегральные сечения различных процессов, использованные в нашей модели, приведены на рис. 1. Внизу на этом рисунке показано сравнение полного сечения рассеяния, полученного путем суммирования всех сечений, использованных в нашей модели, с рекомендуемыми значениями из книги [14].

Дифференциальное сечение ионизации было рассчитано по нерелятивистской BED-модели Кима и Рада [15]. В области энергий электрона меньших 50 эВ дважды дифференциальные сечения ионизации по углу вылета и энергии электрона взяты из [16]. В области больших энергий сечения рассчитаны по формулам первого Борновского приближения, так же как это было сделано для гелия [7], путем масштабирования сечений, рассчитанных для атома водорода (по формуле (13) из [17]), при этом использовались следующие параметры: число электронов на оболочке 2p n = 6, средняя кинетическая энергия электронов на оболочке 2p U = 116.02 эВ [15]. Кроме того, в нашей модели полагалось, что в процессе возбуждения электронных уровней атомов, электроны рассеиваются точно также как в процессе упругого рассеяния. Поскольку, полное сечение возбуждения гораздо меньше сечений упругого рассеяния и ионизации (как можно видеть на рис. 1), то данное упрощение не будет влиять на результаты расчетов.

Расчеты параметров лавины электронов в неоне по методу МК были выполнены в той же постановке, что и для гелия [3]. В начальный момент времени t = 0 в точке $r = (0,0,0)$ в однородном электрическом поле $E(r) = - E{{e}_{z}}$ задавался изотропный моноэнергетический источник, содержащий 1000 электронов с энергией ε₀ = 1 эВ. Концентрация атомов неона полагалась равной числу Лошмидта N_L, и атомы считались неподвижными (температура газа равна нулю). Численное моделирование велось до момента времени t_run, когда число электронов достигало значения ≈3 × 10⁶, данное число достаточно, чтобы за время моделирования достигалось равновесие ансамбля электронов с полем, т.е. функция распределения электронов по энергиям становилась стационарной. В табл. 2 приведены рассчитанные зависимости средней энергии электронов $\bar {\varepsilon }$, дрейфовой скорости υ_d, коэффициентов диффузии D_TN_Lи D_LN_L и ионизации k_ion от приведенной напряженности поля в диапазоне от 15 до 1500 Тд. Формулы для расчета данных величин приведены в [3].

На рис. 2 показана функция распределения электронов по энергиям (ФРЭЭ) в широком диапазоне значений E/N. Для сравнения на данном рисунке также приведено распределение Максвелла со значением средней энергии электронов из МК-расчетов. Отметим интересную особенность, что если в случае гелия рассчитанное энергетическое распределение электронов в области E/N ≤ 150 Тд очень близко к распределению Максвелла [3], то в случае неона в области E/N ≤ ≤ 250 Тд распределение Максвелла не достаточно точно описывает расчетные данные в этой области. При E/N = 500 Тд, наблюдается хорошее согласие, но в области больших значений E/N распределение Максвелла снова не согласуется с расчетными данными.

Как говорилось во Введении, в случае гелия диффузионно-дрейфовое приближение для расчета концентрации электронов становилось неприменимым в полях больших ≈150−200 Тд. Это было связано с появлением большого числа быстрых электронов, энергия которых заметно выше среднего значения. При этом эти быстрые электроны опережали основную массу электронов в лавине. В результате угловое распределение электронов становилось существенно анизотропным, и таким образом нарушалось приближение Лоренца, − двучленное разложение функции распределения электронов (ФРЭ) по углам. На рис. 3 для момента времени t_run приведено угловое распределение электронов, а также распределение линейной концентрации электронов вдоль оси симметрии лавины, которое для удобства представления поделено на число электронов в момент времени t_run. Пространственное распределение электронов сравнивается с распределением, полученным из решения диффузионно-дрейфовых уравнений (ДДУ) (1)–(2) [3], с коэффициентами, полученными в МК-расчетах.

Как можно видеть на рис. 3 в неоне диффузионно-дрейфовое приближение применимо вплоть до E/N = 500 Тд, хотя при этом значении угловое распределение электронов уже заметно анизотропно. В области больших полей уже наблюдаются заметные отличия между решением ДДУ- и МК-расчетами: здесь мы видим, что появляется большое число быстрых электронов, которые опережают основную массу электронов, и пространственное распределение становится несимметричным и вытянутым вдоль направления действия электрической силы.

Выполним теперь сравнение полученных нами данных с результатами расчетов по программе BOLSIG+ [5]. Расчеты по этой программе были выполнены с тем же набором сечений, который использовался в нашей МК-модели с той лишь поправкой, что в программе BOLSIG+ вместо полного сечения упругого рассеяния использовалось транспортное сечение. Также отметим, что расчеты по BOLSIG+ выполнены с использованием опции “gradient expansion with PT”, которая позволяет вычислять два типа дрейфовой скорости (подвижности) электронов “flux” и “bulk”. Разница между ними заключается в том, что тип дрейфовой скорости электронов “flux” это равновесная средняя скорость ансамбля электронов, а “bulk” это равновесная скорость центра масс лавины электронов [18]. В отсутствии ионизации значения данных величин будут совпадать. Отметим, что в табл. 2 приведены значения скорости центра масс лавины электронов [3]. На рис. 4 в диапазоне E/N = 15–1000 Тд приведены зависимости средней энергии, подвижности и дрейфовой скорости электронов, а также коэффициента ионизации, рассчитанные методом МК и по программе BOLSIG+. Как можно видеть, в области E/N < 150 Тд значения подвижности “bulk” электронов, полученные по обеим методикам, близки между собой, относительная разница здесь не превышает 10%, но с ростом E/N разница увеличивается, достигая почти 60% при E/N = 1000 Тд. Значения средней энергии электронов, рассчитанные методом МК, заметно превышают значения, полученные с помощью BOLSIG+, начиная с E/N ≈ 500 Тд, и эта разница также увеличивается с ростом E/N, достигая примерно 50% при E/N = = 1000 Тд. Эта связано с тем, что как говорилось выше, при больших значениях E/N нарушается приближение Лоренца. Отметим также, что относительная разница в значениях коэффициента ионизации, вычисленных по обеим методикам, не превышает 25%.

3. СРАВНЕНИЕ С ЭКСПЕРИМЕНТАЛЬНЫМИ ДАННЫМИ

К сожалению, в отличие от гелия, для неона экспериментальные данные по транспортным коэффициентам электронов в области средних и больших значений E/N крайне скудны. Имеются измерения [19] дрейфовой скорости электронов и отношения D_L/μ_e в диапазоне значений приведенной напряженности E/N от 1.4 до 339 Тд, выполненные методом TOF, а также измерение [20] отношения D_T/μ_e в области 1.41 < E/N < 5650 также в рамках эксперимента TOF. Что касается измерений первого коэффициента Таунсенда α_ion в неоне, то здесь имеется достаточное количество экспериментальных данных [21–24], в том числе и в области больших значений E/N.

Для сравнения с экспериментальными данными были рассчитаны величина первого коэффициента Таунсенда α_ion и дрейфовой скорости υ_SST для условий эксперимента “steady-state Townsend” (SST), а также дрейфовая скорость υ_TOF и отношения D_L/μ_e и D_T/μ_e для условий эксперимента “time-of-flight” (TOF). Расчеты проводились в той же постановке, что и для гелия [3], формулы для расчета величин α_ion, υ_SST, υ_TOF, D_L/μ_e и D_T/μ_e приведены в той же работе. Значения рассчитанных коэффициентов приведены в табл.3.

На рис. 5 рассчитанная зависимость коэффициента Таунсенда сравнивается с экспериментальными данными [21–24]. Рассчитанные значения α_ion очень хорошо согласуются с данными всех измерений в области относительно небольших значений E/N, в отличие от гелия, для которого в области малых полей экспериментальные значения α_ion заметно превосходили теоретические [3], что, скорее всего, связано с влиянием реакций Пеннинга, − если энергия возбуждения атома гелия превышает энергию ионизации атомарной или молекулярной примеси A, то могут протекать реакции ${\text{He}}{\kern 1pt} * + \;A \to {\text{He}} + {{A}^{ + }} + {{e}^{ - }}$ с ионизацией примеси. В области средних и больших E/N рассчитанные значения коэффициента Таунсенда хорошо согласуются с данными [22] во всем исследованном диапазоне значений E/N, и с данными [21] в области E/N < 600 Тд, но при больших E/N экспериментальные значения из [21] заметно превосходят рассчитанные значения. Подобный результат был получен и для гелия [3]. Это расхождение можно объяснить тем, что в эксперименте [21] в области полей E/P > > 100 В/(см · торр) (288 Тд) не достигалось равновесие электронов с полем.

На рис. 6 сравниваются рассчитанные и измеренные зависимости дрейфовой скорости и коэффициентов D_L/μ_e и D_T/μ_e. Теоретические значения υ_TOF очень хорошо согласуются с измеренными в работе [19] в диапазоне 15 ≤ E/N < 200 Тд, но в области больших полей экспериментальные значения заметно превышают теоретические, и относительная разница составляет примерно 10%, что превышает заявленную в [19] точность измерений 5%. Подобное поведение наблюдается и для отношения D_L/µ_e. В области E/N ≥ 200 Тд экспериментальные значения превышают теоретические, и относительная разница также больше заявленной точности измерений 15%. Отметим, что для гелия наблюдалось хорошее согласие расчетных и экспериментальных данных во всем исследованном в эксперименте [19] диапазоне значений E/N. Поскольку для неона (кроме единственной работы [19]) отсутствуют экспериментальные данные для величин υ_TOF и D_L/μ_e, то мы не можем ответить на вопрос о причинах наблюдаемого расхождения расчетных и экспериментальных значений в области E/N ≥ 200 Тд: вопрос о том, связано ли это с погрешностью эксперимента, которая выше заявленной в [19], или обусловлено недостатками представленной здесь теоретической модели, остается открытым.

На рис. 6 также видно, что расчетные значения D_T/µ_e хорошо согласуются с данными [20] в области E/N < 850 Тд, но в более сильных полях теоретические значения заметно превосходят экспериментальные и качественное поведение кривых отличается. Это связано с тем, что в эксперименте не достигалось равновесия электронов с полем в области больших значений E/N, о чем авторы [20] упоминают в своей статье. Отметим, что экспериментальные данные [20] приведены в книге [14] в качестве рекомендованных значений D_T/µ_e в диапазоне 1.41 < E/N < 5650 без каких-либо замечаний о том, что в области сильных полей эти экспериментальные данные могут быть недостоверны.

ЗАКЛЮЧЕНИЕ

Методом Монте-Карло вычислены зависимости от приведенной напряженности поля E/N в диапазоне от 15 до 1500 Тд средней энергии электронов, дрейфовой скорости электронов, частоты ионизации, а также коэффициентов продольной и поперечной диффузии, которые рекомендуются для использования в расчетах эволюции концентрации электронов путем решения диффузионно-дрейфового уравнения. Выполнено сравнение коэффициентов, рассчитанных методом Монте-Карло (МК), и из решения кинетического уравнения Больцмана по программе BOLSIG+. Показано, что основное отличие в поведении зависимостей, рассчитанных по обеим методикам, связано с неприменимостью в области больших значений E/N используемого в BOLSIG+ двучленного приближения для угловой части ФРЭ. Отмечена также ограниченность применимости диффузионно-дрейфового приближения к расчету пространственно временной эволюции электронов в сильных электрических полях, где при очень высоких значениях E/N необходимо учитывать появление большого число быстрых электронов.

С целью верификации модели МК вычислены значения дрейфовой скорости υ_TOF, отношений D_T/µ_e и D_L/µ_e, а также коэффициента ионизации Таунсенда α_ion/N, которые сравнивались с данными измерений. Рассчитанные значения коэффициента Таунсенда согласуются с экспериментальными данными во всем исследованном диапазоне значений E/N. Рассчитанные значения дрейфовой скорости и отношения D_L/µ_e в области E/N < 200 Тд хорошо согласуются с данными измерений [19], но в области больших значений E/N относительная разница между расчетными и измеренными значениями немного превышает заявленную в [19] экспериментальную погрешность. Расчетные значения величины D_T/µ_e в области E/N < 850 Тд хорошо согласуются с экспериментальными данными [20], но заметно превышают их в более сильных полях, скорее всего, в связи с нарушением равновесия ансамбля электронов с полем.