CC BY
19
УДК 544.431.2+519.688
И. И. Митричев, А. В. Женса, Э. М. Кольцова*
Российский химико-технологический университет им. Д.И. Менделеева, Москва, Россия 125480, Москва, ул. Героев Панфиловцев, д. 20 , корп. 1 * e-mail: kolts@muctr.ru
ПОСТРОЕНИЕ ТЕРМОДИНАМИЧЕСКИ НЕПРОТИВОРЕЧИВЫХ
ПОВЕРХНОСТНЫХ МИКРОКИНЕТИЧЕСКИХ МЕХАНИЗМОВ
В работе сообщается о создании программного обеспечения, которое обеспечивает построение термодинамически непротиворечивых микрокинетических механизмов. Рассматривается пример усовершенствования термодинамически несогласованного механизма реакции NO+CO на платине до высокой степени согласованности по изменению энергии Гиббса с одновременным улучшением совпадения степени превращения веществ с данными эксперимента. Предлагаемый подход также может быть использован при построении новых механизмов по литературным данным при полностью или частично неизвестных параметрах отдельных стадий.
Ключевые слова: термодинамическая непротиворечивость; микрокинетическое моделирование; компьютерное моделирование; каталитическая реакция, взаимодействие NO и CO.
Во второй половине XX века был сделан большой шаг по направлению к пониманию явлений, происходящих на поверхности при протекании каталитических реакций. В 1985 году К0гекоу и ЗшЬге показали возможность использования кинетических параметров, найденных в условиях сверхвысокого вакуума, в моделировании промышленного реактора синтеза аммиака [1]. Они использовали полную многостадийную кинетическую схему,
основанную на работах в. Ег1 Бите8ю [2] отказался от дополнительных предположений о лимитирующей стадии и положил начало микрокинетическому моделированию. В классическом микрокинетическом
моделировании, развитом в работах X Бите8Ю, Ь.8сЬт1&, О. БеШзсЬтапп [3,4] каждая стадия имеет свои параметры и соответствует элементарному акту на поверхности катализатора. Преимущество микрокинетического подхода -возможность внедрения кинетического механизма в макроскопический расчет химического реактора.
Кинетический механизм (здесь под механизмом мы будем понимать совокупность отдельных стадий - схему каталитической реакции, вместе с численными значениями параметров, входящих в константу скорости каждой стадии) должен обладать двумя основными свойствами [5]:
1) стехиометрической согласованностью;
2) термодинамической непротиворечивостью.
Первое свойство означает сохранение структурных видов (атомов или групп) в каждой стадии, второе означает то, что кинетические параметры не могут быть выбраны независимо -необходимо выполнение ограничений на термодинамическом уровне, для того чтобы
каждой газофазной реакции соответствовала единственная константа равновесия.
Как отмечают У1асЬо8 и др., существуют вычислительные программы для автоматической генерации механизмов реакций, но оценка кинетических параметров все еще остается сложной задачей [6]. Квантово-механические и ИЫ^ЕР предсказания энергий активацию имеют точность 5 ккал/моль [4], кроме того, квантовые вычисления возможны для конкретной структуры поверхности, а для моделирования промышленного катализатора требуется получить усредненные параметры. Поэтому, для построения детального механизма каталитических реакций, применимого на макроуровне, требуется модификация начальных значений параметров отдельных стадий.
Для решения задачи нахождения кинетических параметров поверхностных реакций было создано программное обеспечение на основе генетического алгоритма бинарного кодирования, который решает обратную кинетическую задачу. Отличительной особенностью созданного программного обеспечения является возможность получения термодинамически непротиворечивых механизмов - как общей (по энергии Гиббса, методика ОеШзсЬтапп [4]) согласованности, так и отдельно энтальпийной и энтропийной (методика У1асЬо8 [7]). Методика может быть применена как к вновь создаваемым механизмам, так и к ранее созданным. В данной работе мы покажем возможность создания термодинамически непротиворечивого механизма с общей непротиворечивостью на основе механизма реакции N0+00 на Р^Ш) МапШ и др. [8] с сохранением кинетической схемы процесса.
В качестве математической модели реактора использовалась модель изотермического реактора с неподвижным слоем катализатора:
dz
= a ■R
(1)
где Ci - концентрация i-го газа (кмоль/м3);
vSup - средняя скорость потока (м/с);
Ri,mol - мольная скорость изменения концентрации
i-го вещества в результате поверхностной
химической реакции (кмоль/м2/с);
z - осевая координата реактора (м);
a - отношение активной поверхности катализатора
к объему реактора (м2/м3).
Модель отличается от использованной авторами механизма [9] включением приведенной скорости потока внутрь конвективной производной, что соответствует случаю изменения плотности среды при химической реакции. Модель может быть получена упрощением общей трехмерной системы уравнений движения газа механики сплошных гетерогенных сред с учетом неподвижности катализатора и отсутствия диффузионного перемешивания.
В рамках нашего подхода используется приближение среднего поля: поверхность рассматривается как однородная, и в модели различаются не отдельные адсорбаты, а их концентрации Sj или степени заполнения 0j. Используется приближение стационарного состояния (Steady State Approach, SSA). Полностью описание модели можно найти в [9]
Скорость стадий рассчитывалась по модифицированному уравнению Аррениуса:
R = 4 (T )*'- ^7 RT П j Пс
(2)
где Яг - скорость г-ой поверхностной стадии (кмоль/м2/с);
А - предэкспоненциальный множитель (кмоль,м,с);
в - температурная экспонента;
Е - энергия активации (Дж/моль);
Я - универсальная газовая постоянная
(Дж/моль/К);
Т - температура (К);
Sj - концентрация j-го адсорбата (кмоль/м2); То - номинальная температура, равная 300 К; шг - частный порядок г-ой стадии по _]-ому веществу.
Предэкспонента адсорбции рассчитывалась с использованием коэффициента прилипания, см. подробнее [4].
Для всех адсорбционных стадий в данной работе частный порядок стадии по свободным реакционным центрам полагался равным стехиометрическим коэффициентам. Например, диссоциативная адсорбция азота на платине (стадия 4) происходит с участием двух реакционных центров, следовательно, Ш4,№ = 2.
Начальный диапазон поиска значений для предэкспоненциальных множителей и энергий активации был взят таким образом, чтобы включить данные моделирования и расчетов ряда авторов [8, 10]. Значения этих параметров из публикации, с которой производилось сравнение [8], приведены в табл. 1.
Оценкой для термодинамической
непротиворечивости была сумма квадратов отклонений гда изменения энергии Гиббса стадии, рассчитанной через константу равновесия, от той же величины, рассчитанной через значения энергии Гиббса отдельных участников стадии, по всем стадиям по 7 температурным точкам в диапазоне температур от 300 до 800 K. Температурные точки были распределены равномерно. Энергии Гиббса отдельных веществ-участников стадии задавались в формате термодинамических полиномов NASA [11].
Производилась минимизация суммы квадратов отклонений экспериментальной конверсии от расчетной по формуле: Nэксп. I
X X[(yni эксп.-Уш расч.)2] + (3)
n =1 i =1 10
где n - номер экспериментальной точки; Notch. - число экспериментальных точек; i - номер газофазного вещества; I - число газофазных веществ, по которым ведется подбор конверсии;
уэксп - экспериментальная конверсия (%);
со,
урасч. - расчетная конверсия (%); X - весовой коэффициент, определяющий значимость второго слагаемого (меньшие значения соответствуют большей значимости).
Второе слагаемое в формуле (5), включающее гда, позволяет учитывать в оценке качества подбора термодинамическую согласованность механизма. Экспериментальные данные по конверсии взяты из [12], так же как в [8]. Параметры реактора, используемые в модели, уточнены в [8].
По веществам, которые отсутствуют в исходной смеси (№0), вместо конверсии рассчитывается сходная величина, где взамен начальной концентрации используется опорная концентрация. В качестве опорной выбрана концентрация другого реагента - СО.
Случайно взятый набор параметров, где каждый параметр выбран из начального диапазона поиска, имеет значение гда порядка 1010-1011.
Сначала подбор параметров ведется при Х=20, что полностью снимает требование термодинамической непротиворечивости. Когда достигается совпадение полученных в результате моделирования и экспериментальных кривых конверсии с абсолютной ошибкой не более 5-7% во всем диапазоне экспериментальных данных, то X снижалось до 7. Снижение X ниже 5 может привести к ухудшению предсказания конверсии.
т
Таблица 1
№ Элементарная стадия кинетической схемы Параметры из оригинальной работы [7] Параметры, подобранные в данной работе
Ас, (с-1) или У0 Е, кДж/моль А0, (с-1) или У0 в Е, кДж/моль
1 СО + * ^ СО* 0,89 0 0,855 0 0
2 СО2 + * ^ СО2* 0,005 0 0,0054 0 0
3 N0 + * ^ N0* 0,6 0 0,812
4 N2 + * ^ 2№ 0,001 0 0,001067 -0,455 55,77
5 N20 + * ^ N20* 0,001 0 0,0198
6 СО* ^ С0 + * 11013 133,88 3,5141014 -0,082 136,41
7 С02* ^ С02 + * 11013 15,185 2,368^ 1012 0,00671 1,68
8 N0* ^ N0 + * 11013 109,085 3,575^ 1018 0,47 136,59
9 2№ ^ N2 + * 11011 111,276 6,506^ 1012 0 110,03
10 N20* ^ N20 + * 11013 23,43 8,97110" -0,313 31,43
11 СО2* + * ^ С0* + 0* 1-1011 77,147 2,658^ 1015 -0,077 173,96
12 С0* + 0* ^ С02* + * 1-1011 97,269 4,7924012 0 99,66
13 N0* + * ^ N + 0* 1-1011 61,0 1,6311013 0 78,02
14 N + 0* ^ N0* + * 1-1011 152,68 2,919^ 1014 0,0195 135,02
15 N0* + N ^ N20* + * 1-1011 89,066 7,330Т012 0 75,93
16 N20* + * ^ N0* + N 1-1011 15,094 2,126^ 1015 -0,0445 35,1
Было выполнено 400 итераций генетического алгоритма, из них - 200 - с начальным значением параметра X и 200 итераций со значением Х=7. При проведении последних 200 итераций верхние и нижние границы поиска по каждому из параметров были центрированы относительно наилучших значений, полученных в результате первых 200 итераций. Полученный набор кинетических параметров (табл. 1) является термодинамически непротиворечивым, гда = 5,61 •Ю5. Величина гда для исходного набора параметров [8] гда = 2,61-1010, что соответствует полностью несогласованному набору. Таким образом, в результате применения процедуры поиска параметров, обеспечивающих
термодинамическую непротиворечивость, было достигнуто улучшение общей термодинамической (по энергии Гиббса) непротиворечивости на 5 порядков.
Улучшение предсказательной способности модели по конверсии (рис. 1) в данной работе связано с тем, что авторы [8,9] использовали энергии активации, полученные расчетным путем (ИВ1^ЕР) и предэкспоненты, оцененные по теории переходного состояния, без их изменения.
Исследование выполнено при финансовой поддержке РФФИ в рамках проекта 14-07-00960.
Рис. 1. Экспериментальные [12] и смоделированные значения конверсии. Пустые точки -эксперимент; сплошные линии - моделирование с набором параметров, найденных авторам данной работы; штриховые линии - моделирование с набором параметров из [8]. Концентрации СО отмечены квадратом, N0 - кружком, N20 - треугольником.
Митричев Иван Игоревич, аспирант, ведущий программист кафедры информационных компьютерных технологий РХТУ им. Д. И. Менделеева, Россия, Москва
Женса Андрей Вячеславович, к.т.н., доцент кафедры информационных компьютерных технологий РХТУ им. Д. И. Менделеева, Россия, Москва
Кольцова Элеонора Моисеевна, д.т.н., профессор кафедры информационных компьютерных технологий РХТУ им. Д. И. Менделеева, Россия, Москва
Литература
8. Stoltze P., N0rskov J.K. Bridging the ''pressure gap'' between ultrahigh-vacuum surface physics and high-pressure catalysis // Phys. Rev. Lett. - 1985. - Т. 55: - С. 2502.
9. Dumesic J. A., Trevino A. A. Kinetic simulation of ammonia synthesis catalysis // Journal of Catalysis.
- 1989. - Т. 116. - №. 1. - С. 119-129.
10. The Microkinetics of Heterogeneous Catalysis / Washington, DC: ACS, 1993. - 315 pp.
11. Modeling and Simulation of Heterogeneous Catalytic Reactions: From the molecular process to the technical system. / John Wiley & Sons [Deutschmann O. (ed.)]. - 2011. - 370 pp.
12. Stoltze P. Microkinetic simulation of catalytic reactions // Progress in surface science - 2000. - Т. 65.
- №3: - С. 65-150.
13. Raimondeau S., Vlachos D.G. Recent developments on multiscale, hierarchical modeling of chemical reactors // Chemical Engineering Journal. - 2002. - Т. 90. - С. 3-23.
14. Mhadeshwar A. B., Wang H., Vlachos D. G. Thermodynamic consistency in microkinetic development of surface reaction mechanisms // The Journal of Physical Chemistry B. - 2003. - Т. 107. - №. 46. - С. 12721-12733.
15. Mantri D., Aghalayam P. Micro-kinetic study of reduction of NO on Pt group catalysts // International Journal of Chemical Reactor Engineering. - 2007. -№5:A1.
16. Mantri D., Aghalayam P. Detailed surface reaction mechanism for reduction of NO by CO // Catalysis Today. - 2007. - Т. 119. - С. 88-93.
17. Koop J., Deutschmann O. Detailed surface reaction mechanism for Pt-catalyzed abatement of automotive exhaust gases //Applied Catalysis B: Environmental. - 2009. - Т. 91. - №. 1. - С. 47-58.
18. McBride B. J., Gordon S., Reno M. A. Coefficients for calculating thermodynamic and transport properties of individual species // NASA Report TM-4513 - 1993.
19. Chambers D. C., Angove D. E., Cant N. W. The Formation and Hydrolysis of Isocyanic Acid during the Reaction of NO, CO, and H2 Mixtures on Supported Platinum, Palladium, and Rhodium // Journal of Catalysis. - 2001. - Т. 204. - №. 1. - С. 11-22.
Mitrichev Ivan Igorevich, Zhensa Andrey Vyacheslavovich, Kol'tsova Eleonora Moiseevna* D.I. Mendeleev University of Chemical Technology of Russia, Moscow, Russia. * e-mail: kolts@muctr.ru
CONSTRUCTION OF THERMODYNAMICALLY CONSISTENT SURFACE MICROKINETIC MECHANISMS
Abstract
We report the software development for thermodynamically consistent microkinetic mechanism construction. We consider an example of improvement of thermodynamically inconsistent mechanism of NO+CO reaction on platinum to the high degree of consistency in respect of free energy change with simultaneous improvement of fit for experimental conversion data. The proposed approach could be used in new kinetic mechanisms construction with the use of literature data when reaction step parameters are totally or partially undefined.
Keywords: thermodynamic consistency; microkinetic modeling; computer simulation; catalytic reaction; interaction of NO and CO
