Chemical reactions rates in the TEFIS database on thermophysical properties of dense plasma

Abstract

Since the 1960s at the Keldysh Institute of Applied Mathematics team led by N.N. Kalitkin creates the TEFIS database on the thermophysical properties of substances. In recent years, this database has included a section on the rates of chemical reactions that describe several relevant processes. For this, new mathematical methods for regression of experimental data measured with significant errors were developed. This paper presents the current state of work in this area.

Full Text

ВВЕДЕНИЕ

Моделирование задач химической промышленности проводят по программам газодинамики с химической кинетикой. Примерами таких задач являются расчеты горения водорода в топливных элементах и двигателях, проектирование и оптимизация газофазных реакторов, расчеты вредных выбросов и многие другие. Результаты моделирования существенно зависят от используемых данных по скоростям K(T) химических реакций, зависящих от температуры T. Поэтому важной задачей является получение достоверных аппроксимаций для K(T).

Одним из источников таких данных является регрессия результатов экспериментальных измерений. В таких экспериментах непосредственно измеряют концентрации реагирующих веществ в зависимости от времени, и по этим профилям определяют скорости реакций при фиксированной температуре T и давлении p. Затем экспериментальную зависимость аппроксимируют формулой Аррениуса

$K\left(T\right)=A\exp (-E/T)$ (1)

либо ее обобщением

$K\left(T\right)=A{T}^m\exp (-E/T)$. (2)

Формула (1) основана на следующем рассуждении. Пусть среда является однородной, изотропной и термодинамически равновесной. Тогда частота актов реакции пропорциональна средней частоте соударений частиц и вероятности преодолеть энергетический порог реакции, описываемой в рамках распределения Больцмана. В формуле (2) множитель T^m является эвристическим. Он вводится для учета зависимости частоты столкновений от температуры.

Как правило, по каждой реакции опубликовано много экспериментальных работ, причем диапазон условий в них частично перекрывается, частично различается. Аппроксимации вида (1), (2) собраны в компендиумы [1, 2].

Из-за погрешностей эксперимента результаты разных авторов отличаются, причем зачастую весьма существенно. Чтобы нивелировать индивидуальные отклонения отдельных экспериментов, проводят совместную обработку большого количества данных различных источников. Такие обработки публикуются в компендиумах [1, 2] и в тематических сборниках (см., например, [3–5]).

Однако при этом остается неясным, с какой точностью получены коэффициенты формулы Аррениуса. Во многих случаях для аппроксимирующей кривой (1) приводят экспертную оценку доверительного коридора. Как правило, он составляет ± (0.3 – 0.5) для lgK, т.е. истинное значение K до 2–3 раз больше либо меньше приведенного в сборнике.

С 1960-х гг. в Институте прикладной математики им. М.В. Келдыша РАН коллектив, основанный Н.Н. Калиткиным, создает базу данных ТЕФИС по теплофизическим свойствам веществ [6]. Термодинамическая часть этой базы представлена в работе [7] В последние годы в эту базу включен раздел, посвященный скоростям химических реакций.

Для создания этого раздела были разработаны новые математические методы регрессии экспериментальных данных, измеренных со значительными погрешностями. Отобраны системы реакций, которые описывают ряд актуальных процессов: водородо-воздушное и водород-кислородное горение, термическое разложение этана, образование оксидов серы и сероводорода в воздухе. Для этих реакций cобран компендиум первичных экспериментальных измерений зависимости K(T). Проведена уточненная обработка этих данных с помощью предложенных нами методов. Точность этих аппроксимаций в десятки раз превосходит известные базы данных. В данной публикации представлено текущее состояние работ по этому направлению.

МЕТОД РЕГРЕССИИ

Сделаем замену переменных x = 1 / T, u = lgK. Тогда закон Аррениуса превращается в прямую. Пусть имеется набор экспериментальных точек: аргументов x_j и значений функции u_i ± δ_i, где δ_i – абсолютные погрешности измерений. Для обработки таких данных был разработан новый математический метод регрессии экспериментальных данных, измеренных со значительными погрешностями [8]. Он основан на использовании полиномов, ортогонализованных на множестве экспериментальных точек. Обозначим их через φ_n(x), где n – степень полинома. Введем скалярное произведение

$〈{\phi }_j,{\phi }_k〉=\sum _{i=1}^I\frac{{\phi }_j\left(x_i\right){\phi }_k\left(x_i\right)}{{\delta }_i^2}$. (3)

Здесь I – число экспериментальных точек. Положим φ₀(x) = 1. Выберем такую линейную функцию φ₁(x) = x – c, которая ортогональна j₀(x) в смысле скалярного произведения (2): $〈{\phi }_0, {\phi }_1〉=0$. Отсюда нетрудно найти $c=\overline{x}$, где

$\overline{x}=\frac{1}{\Delta }\sum _{i=1}^I{\delta }_i^{-2}{x}_i, \Delta =\sum _{i=1}^I{\delta }_i^2$. (4)

Аналогично, но более громоздко строятся полиномы второй и более высоких степеней [8].

Составим аппроксиманту φ(x) = a₀φ₀(x) + a₁φ₁(x). Коэффициенты a₀, a₁ найдем методом наименьших квадратов

$\sum _{i=1}^I\frac{{(u_i-\phi (x_i\left)\right)}^2}{{\delta }_i^2}\to \underset{a}{\min }$. (5)

Подставляя в (4) явный вид φ(x) и проводя минимизацию, получим

$a_j=\frac{⟨u,{\phi }_j⟩}{⟨{\phi }_j,{\phi }_j⟩}, j=\mathrm{0,1}$. (6)

Для коэффициентов a_j и аппроксимирующей кривой были найдены и строго обоснованы статистически достоверные оценки доверительных интервалов. Доказано

Утверждение. Доверительные коридоры коэффициентов δa_j и аппроксиманты δφ(x) равны

$\delta a_j=\pm {⟨{\phi }_j,{\phi }_j⟩}^{-1/2}, \delta \phi \left(x\right)=\pm {\left(\sum _{j=1}^n\delta a_j^2{\phi }_j^2\left(x\right)\right)}^{1/2}$. (7)

Коэффициент a_j является достоверным, если |δa_j| < a_j. Слагаемые φ_j с недостоверными коэффициентами следует отбросить. В [9] проведено сравнение описанного метода и регрессии по неортогональной системе степеней j_j(x) = x^j. Показано, что оценки (7) существенно (в 10 и более раз) точнее таковых для неортогональных полиномов.

ПРИМЕР ОБРАБОТКИ

В качестве примера рассмотрим реакцию N₂O + M → N₂ + O + M. По ней в базе данных NIST [1] приведены ссылки на 18 оригинальных работ, содержащих в общей сложности ~640 точек. Эти данные в переменных 1 / T – lgK приведены на рис. 1. Массив точек похож на размытую прямую, хотя местами это размытие значительно.

Рис. 1. Реакция N₂O + M → N₂ + O + M. Точки – эксперименты. Сплошная линия – аппроксимация ТЕФИС, пунктир – границы ее доверительного интервала

Аппроксимация, найденная по описанному выше методу, имеет вид $y=a_0+a_1(x-\overline{x})$, lga₀ = 15.040 ± 0.007, a₀ = 1.1081 ± 0.006, $\overline{x}$ = 5.4580. Она показана на рис. 1. Полученная погрешность составила 1%. На рис. 1 показаны границы доверительного коридора (6). Видно, что этот коридор значительно уже массива экспериментальных точек. С визуальной точностью все три линии сливаются, поэтому в целях иллюстрации мы привели небольшой участок графика более крупно.

На рис. 2 приведен десятичный логарифм отношения нашей аппроксимации к нескольким широкодиапазонным аппроксимациям [10–15] из базы NIST. Видно, что уклонение последних достигает от 0.2 до 0.5 для lgK, т.е. формулы для K(T) из NIST от 1.6 до 3 раз больше либо меньше нашей аппроксимации. Такое расхождение согласуется с оценками точности, приведенными в базе NIST. Таким образом, наша кривая в 60200 раз точнее ранее известных аппроксимаций.

Рис. 2. Реакция N₂O + M → N₂ + O + M. Отношение аппроксимации ТЕФИС к аппроксимациям из базы NIST. Цифры около линий – номер ссылки по списку литературы

РЕДУКЦИЯ СИСТЕМЫ РЕАКЦИЙ

Отдельной проблемой является выбор системы реакций. Нужнo нaпиcaть дocтaтoчнo пoлную cиcтeму химичecких peaкций, учитывaющую вce реакции, существенные в рассматриваемых условиях. Ecли нeкoтopыe вaжныe peaкции нe учтeны, тo этo мoжeт пpивoдить к кaчecтвeннo нeпpaвильным peзультaтaм.

В литературе традиционно используют так называемый анализ чувствительности (sensivity analysis) системы реакций по отношению к вариации скоростей [16]. Считается, что если изменение скорости отдельно взятой реакции не приводит к значимому изменению концентраций, то эту реакцию можно откинуть.

Однако для анализа больших систем химических реакций такая процедура недостаточно надежна. Причина этого заключается в следующем. Известно [17], что обратная задача кинетики (то есть нахождение K(T) по временному профилю концентраций) некорректна. Это означает, что зависимость временных профилей концентраций от скорости отдельно взятой реакции достаточно слабая, и для определения чувствительности к скорости данной реакции требуется избыточно сильная вариация K(T), превосходящая погрешность ее нахождения. Иллюстрацией этого является возможность введения эффективной брутто-реакции вместо детальной кинетической схемы.

Нами предложена следующая процедура отбора. Выберем максимально полную систему реакций. Проведем расчеты кинетики полной системы реакций в заданном диапазоне температур. Полученные временные профили концентраций примем за эталонные. Далее исключим одну из реакций и повторим расчет. Если отличие полученных профилей концентраций от эталонных заметно меньше, чем погрешность экспериментального измерения концентраций конечных продуктов, то данной реакцией можно пренебречь. В противном случае реакцию необходимо оставить. Реакции целесообразно исключать в порядке увеличения K. Система реакций считается окончательной, когда ни одну из оставшихся реакций нельзя отбросить. Насколько нам известно, ранее такая процедура не использовалась.

Такой подход обладает рядом преимуществ по сравнению с известными методами анализа чувствительности. Во-первых, редуцированная система реакций дает те же практические результаты, что и оригинальная система. Во-вторых, предложенная процедура более экономична по сравнению с традиционными подходами. Для каждой исследуемой реакции требуется однократное решение задачи кинетики с «выключением» этой реакции, в то время как в методах анализа чувствительности требуется выборка таких решений с многократным варьированием K(T).

Проиллюстрируем описанную процедуру на примере системы реакций, описывающих водород-кислородное горение. Наиболее полная система реакций, которую удалось найти в литературе, содержит 50 реакций (25 прямых и столько же обратных) с участием 9 компонент: O, H, H₂O₂, OH, H₂, H₂O, HO₂, O₃, O₂ [18]. Эта система весьма громоздка, поэтому мы ее не приводим. Были рассмотрены температура T ∼ 3000 К и давление p ∼ 150 атм. Эти условия соответствуют горению топлива в сопле реактивного двигателя.

На рис. 3 показана зависимость важнейших компонент от времени при указанных условиях. Начальные концентрации соответствуют «гремучей смеси»: 2 части H₂ на 1 часть O₂. Расчеты проводились с использованием пакета GACK [19] с автоматическим формированием правых частей и контролем фактической точности в ходе расчета. Видно, что концентрации исходных компонент монотонно убывают. Конечным продуктом является вода, концентрация которой возрастает. Математическая точность расчета была не хуже 1%.

Рис. 3. Концентрации важнейших компонент при горении H₂ + O₂

Применение предложенного метода позволило описать этот процесс системой 10 реакций с участием 7 компонент. Эти реакции перечислены в табл. 1. Отклонение профилей концентраций, внесенное редукцией системы реакций, показано на рис. 4. Оно не превышает 6%, что соответствует характерной точности экспериментального измерения концентраций.

Таблица 1. Реакции с участием H, O

Реакция	lgA	E	dK, %
H+O2+M → HO2+M	16.54 ± 0.04	0	11
H+O2+M ← HO2+M	16.5 ± 0.2	2.03 ± 0.04	70
H+OH+M → H2O+M	17.04 ± 0.11	0	30
H+OH+M ← H2O+M	14.12 ± 0.03	4.01 ± 0.04	9
H+HO2 → H2+O2	0.1372 ± 0.0001	0	0.02
H+HO2 ← H2+O2	12.96 ± 0.04	1.97 ± 0.15	12
O+OH → O2+H	13.08 ± 0.01	0	2
O+OH ← O2+H	13.989 ± 0.004	0.643 ± 0.003	1
OH+H ← H2+O	13.85 ± 0.04	0.479 ± 0.007	13
OH+H2 → H+H2O	13.31 ± 0.03	0.229 ± 0.003	11

Рис. 4. Отклонения концентраций от эталонных при редукции системы реакций

БАЗА ДАННЫХ

В базе ТЕФИС рассмотрены следующие системы реакций.

1. 20 реакций с участием H, N, O, у которых характерное время протекания при T < 1000 K, p ∼ 1 атм. не превышает 0.01 с. Эти реакции имеют значение для задач плазмохимии.
2. 15 реакций термического разложения этана [20] при давлениях более 1 атм.
3. 10 реакций водород-кислородного горения, существенных при T < 3000 K, p ∼ 150 атм.
4. 16 реакций с участием S, H, O, N, у которых характерное время протекания при T < 300 K, p ∼ 1 атм. не превышает 0.1 с. Эти реакции имеют значение при расчете распространения загрязнений от химических производств.

Для указанных реакций были собраны результаты первичных экспериментальных измерений скоростей в зависимости от температуры. Насколько нам известно, ранее такие компендиумы не создавались, а экспериментальные данные были рассеяны по отдельным публикациям. Эти данные были обработаны с помощью разложения по ортогонализованным полиномам. Во всех случаях оказалось достаточно постоянного и линейного членов разложения. Уже квадратичный член разложения оказывался недостоверным в пределах имеющихся экспериментальных погрешностей.

Результаты для первой и второй систем реакций опубликованы в [8] и [21], соответственно. Для систем 3 и 4 результаты приведены в табл. 1 и 2. Ранее они не публиковались. Здесь стрелка вправо → соответствует экзотермической реакции, а стрелка влево ← показывает эндотермической. Величины lgA даны в единицах см³/(моль с), энергии E – в эВ. Для ряда реакций из табл. 2 имелась всего по одной экспериментальной точке. Для них регрессия не проводилась. В остальных случаях полученная точность аппроксимаций составляет 1–10% для хорошо изученных реакций и 30–60% для малоизученных. Это в десятки раз точнее известных баз данных.

Таблица 2. Реакции с участием S, H, N, O

Реакция	lgA	E	dK, %
HOSO2+O2 → HO2+SO3	11.46 ± 0.02	0	4
HSO+NO2 → HSO2+NO	12.76	0	–
HSO+O3 → продукты	10.69 ± 0.05	0	12
HSO2+O2 → HO2+SO2	13.26	0	–
O+H2S → OH+SH	13.18 ± 0.15	0	35
O+SH → SO+H	14.3 ± 0.2	0	60
O3+SO2 → SO3+O2	3.78	0	–
OH+H2S → SH+H2O	12.459 ± 0.005	0	1
OH+S → H+SO	13.5957	0	–
OH+SO → H+SO2	13.43 ± 0.08	0	20
S+O2 → SO+O	12.55 ± 0.04	0	11
SH+NO2 → HSO+NO	13.56 ± 0.04	0	7
SH+O2 → HSO+O2	12.48 ± 0.02	0	5
SO+NO2 → SO2+NO	12.905 ± 0.005	0	1
SO+O2 → SO2+O	10.76 ± 0.07	0	17
SO+O3 → SO2+O2	6.08 ± 0.08	0	20

ЗАКЛЮЧЕНИЕ

Выбор значений скоростей реакций оказывает существенное влияние на результаты расчетов кинетики химических процессов. Приведем пример. Ранее мы проводили расчеты изотремического горения водород-кислородной смеси, пользуясь более ранней версией базы ТЕФИС [22]. Исходные реагенты представляли собой «гремучую смесь» при атмосферном давлении. При температуре 2000 К сгорание смеси (т.е. выход концентраций на стационарные значения) происходило за ~ 5 мкс. Аналогичные расчеты проводила Е.Г. Синько [22], пользуясь аппроксимациями из базы NIST [1]. В том расчете сгорание происходило за ~0.1 мкс, т.е. в 50 раз быстрее! Это сравнение выразительно показывает важность использования достоверных данных по скоростям химических реакций.

Работа поддержана Программой стратегического академического лидерства РУДН.