Теплофизика высоких температур, 2023, T. 61, № 1, стр. 68-74

ВВЕДЕНИЕ

Формирование сажи при горении и пиролизе углеводородных топлив является комплексным многостадийным процессом, многие этапы которого до сих пор остаются слабо изученными. Одной из основных трудностей является ограниченный перечень экспериментальных методов, позволяющих проводить прямой анализ процессов, протекающих в пламенах при температурах в диапазоне 1000–2000 К и характеризуемых достаточно быстрой кинетикой. Традиционно процесс образования сажи подразделяют на три этапа: на первом в ходе горения и пиролиза топлива происходит образование полиароматических углеводородов (ПАУ), на втором молекулы ПАУ объединяются в зародыши сажевых частиц (как правило, имеющие размер 1–3 нм), на третьем эти зародыши увеличиваются и объединяются в более сложные фрактальные комплексы микрометровых размеров. Подробный обзор актуальных представлений о механизмах сажеобразования дан в [1]. На текущий момент в научном сообществе сформировалась относительно цельная картина начальных стадий пиролиза простейших углеводородов и строения уже сформированных зрелых сажевых частиц. Механизмы же, ответственные за промежуточные процессы (нуклеацию и рост зародышей сажи), являются предметом наиболее бурных дискуссий.

По мере роста сажевая частица претерпевает ряд существенных морфологических изменений, которые принято называть “старением”. Сам по себе термин “старение” не имеет однозначного определения [2] и служит для обобщенного описания целого ряда физико-химических процессов, протекающих в сажевых частицах при высоких температурах. К таковым относятся графитизация, карбонизация, рост (из газовой фазы) и коалесценция (с другими подобными частицами). В совокупности эти процессы приводят к изменениям оптических и термодинамических свойств сажевых частиц, увеличению соотношения [C]/[H] и плотности, упорядочению структуры на наноразмерных масштабах. Степень зрелости сажевых частиц сказывается также и на температуре их сублимации, которая может изменяться в широких пределах: от 2500 до 4500 К [3, 4].

Исторически же термин “старение” был введен в первую очередь для описания эффекта снижения реакционной способности поверхности сажевой частицы в процессе ее эволюции в условиях высоких внешних температур. Предполагается, что снижение реакционной способности в ходе старения связано преимущественно с увеличением отношения [С]/[H] и графитизацией поверхности – оба этих процесса снижают концентрацию потенциальных точек связывания для молекул газового окружения. В численных моделях роста сажевых частиц для учета данного эффекта вводится стерический параметр α, выступающий в составе предэкспоненциального множителя в уравнениях химических реакций на поверхности частицы. Сравнительный анализ, проведенный в работе [5], показывает, что в разных моделях значения параметра α различаются более чем на порядок. Подобный разброс подчеркивает наличие проблемы в детальном понимании процессов, ответственных за старение сажевых частиц и за сопутствующие ему изменения в морфологии поверхности.

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

ОПИСАНИЕ МОДЕЛИ И МЕТОДИКА МОДЕЛИРОВАНИЯ

Потенциал межчастичного взаимодействия. Все молекулярно-динамические расчеты проводились с использованием программного пакета “LAMMPS” [6, 7] и потенциала ReaxFF [8]. Шаг численного интегрирования уравнений движения составлял 0.1 фс. Для вывода систем на заданные температуры использовался термостат Нозе–Хувера [9] с константой связывания 50 фс.

Выбор потенциала межатомного взаимодействия – ключевой аспект любого МД-расчета, от которого зависит достоверность полученных данных. Описание поведения углерода при высоких температурах требует учета сложной химии этого элемента, атомы которого могут изменять степень гибридизации, образовывать ковалентные связи различной кратности, а также проявлять и более сложные многоатомные эффекты, такие как π-сопряжение. Описание всего этого перечня явлений – крайне сложная задача для аналитических МД-потенциалов, в особенности для аморфных [10, 11] и жидких фаз углерода [12]. В последние несколько лет эта проблема начала устраняться благодаря развитию машинно-обучаемых потенциалов, таких как GAP [13] или SNAP [14], однако на текущий момент они позволяют моделировать только углеродные системы без добавления каких-либо других атомов [15].

Ранее авторы работы сравнили нескольких наиболее распространенных моделей аналитических реакционных потенциалов в задаче моделирования углерода при высоких температурах [16] и продемонстрировали хорошее согласие с экспериментальными и квантово-химическими данными для потенциала ReaxFF с набором параметров CHO-2017 [17]. Результаты для этого набора параметров показали хорошее согласие с экспериментальными данными при моделировании оксида графена [18]. Также данный набор параметров был применен для анализа изменения морфологии сажевых частиц при высоких температурах [19]. Многочастичный потенциал ReaxFF способен описывать акты образования и разрыва химических связей и требует для расчета энергии только координаты атомов, без задания в явном виде ковалентных связей. Помимо задач, связанных с высокотемпературной химией углеводородов, данный потенциал активно используется для моделирования высокоэнергетических соединений [20].

Классическая молекулярная динамика не является единственным численным подходом, используемым для атомистического моделирования сложных углеродных систем. Ряд квантово-химических методов, например метод теории функционала плотности, также активно применяются для изучения свойств углеродных материалов и наночастиц [21–24]. Преимуществом этих методов является то, что они не требуют внедрения в модель дополнительных эмпирических параметров и представляются, таким образом, более универсальными, чем классическая МД. Это, однако, обходится ценой больших вычислительных затрат – моделирование в рамках теории функционала плотности систем, содержащих даже сотню атомов, на наносекундных временах является практически невыполнимой задачей.

Классическая МД тоже является достаточно ресурсоемким подходом и не позволяет выходить далеко за пределы наноразмерных систем. Для описания кинетики процессов образования и роста сажевых частиц на макроскопических временных и пространственных масштабах применяются кинетические модели, в основе которых лежат системы дифференциальных уравнений, описывающих протекающие в исследуемой системе реакции [25–27]. При таком подходе все потенциально возможные механизмы течения реакций между реагентами должны быть учтены в модели заранее. Атомистическое же моделирование позволяет выявлять в том числе и неочевидные механизмы, влияющие на кинетику реакций, которыми потом могут быть дополнены уже существующие кинетические модели. В данной работе рассматривается самый простой сценарий взаимодействия сажевой частицы с молекулой газового окружения, ограниченный моделированием налетания одиночной молекулы ацетилена. Такой упрощенный подход позволяет продемонстрировать влияние наноструктуры поверхности сажевой частицы на кинетику ее роста. В дальнейшем планируется распространить исследование на более широкий спектр молекул, включающий не только ацетилен, но и более тяжелые продукты пиролиза, наблюдающиеся в газовых смесях. Используемый потенциал [17] позволяет описывать системы, состоящие из атомов углерода, водорода и кислорода, т.е. дает возможность моделировать в том числе и процессы окисления.

Генерация моделей сажевых частиц. Одной из основных подзадач данной работы является попытка сформировать морфологически реалистичную атомистическую структуру для зародыша сажевой частицы. Такой зародыш, как правило, имеет размер порядка нескольких нанометров и обладает менее упорядоченной и менее графитированной структурой, чем более взрослые сажевые частицы.

Вопрос о том, в рамках каких именно механизмов происходит формирование и рост зародышей сажевых частиц, является по-прежнему открытым. Долгое время предполагалось, что в формировании участвуют исключительно переконденсированные ПАУ (например, пирен [28–31], коронен [30–33]), однако даже при сравнительно высоких температурах (Т = 1500–2000 К) они оказываются недостаточно реакционно способными [1], чтобы объяснить наблюдаемую скорость нуклеации сажи. В качестве альтернативных объяснений выдвигались, например, гипотезы о резонансно-стабилизированных димерах ПАУ [34, 35].

Модель молодой сажевой частицы. Для создания модели “молодой” сажевой частицы в расчетную ячейку с размерами 7 × 7 × 7 нм³ помещались молекулы ацетилена и молекулы ПАУ, как в экспериментальной работе [36] (рис. 1). Основным аргументом в пользу моделирования именно гетерогенной смеси из крупных ПАУ и легких молекул ацетилена стали результаты двух недавних экспериментальных исследований [36, 37]. В работе [36] с помощью методов электронной микроскопии была установлена атомистическая структура ПАУ, обнаруженных в областях образования сажи. По мнению авторов, эти молекулы, обладающие массами 180–600 Д, и выступают начальными реагентами в процессе нуклеации сажи. Однако результаты [37] указывают на то, что концентрация ПАУ в типичных пламенах недостаточна для объяснения скорости нуклеации сажевых частиц – эту роль авторы, наоборот, отводят бензолу и более легким молекулам. Подобные противоречия могут объясняться участием в формировании зародышей сажевых частиц как более крупных полициклических продуктов, так и более легких неароматических молекул.

Рис. 1.

Схема формирования молодой сажевой частицы в МД-расчете: (а) ‒ исходная конфигурация молекул в системе, t = 0 нс; (б) ‒ наиболее крупный кластер, t = 5 нс; (в) ‒ наиболее крупный кластер, t = 7 нс.

Используемый здесь подход отличается от большинства молекулярно-динамических расчетов, в которых формирование сажевых частиц моделируется в условиях первоначально однокомпонентной газовой фазы. В качестве такой фазы принимается либо набор ПАУ (например, пирен [28–31] или коронен [30–33] – из предположений, что все молекулы первоначальной углеводородной смеси уже трансформировались в ароматические соединения), либо, наоборот, исходные молекулы топлива, например ацетилен [38, 39], метан [40], н-декан [41] и т.п. В первом варианте из-за высокой стабильности и малой реакционной способности ПАУ процессы их димеризации и роста оказываются существенно затруднены (в условиях отсутствия радикалов или более активных молекул), что вынуждает проводить моделирование при температурах Т = 2500–3000 К [42], заметно превышающих реальные значения в экспериментах. Во втором варианте МД расчет требует значительно больших временных масштабов моделирования: сначала требуется дождаться формирования ПАУ из первоначальной смеси (что само по себе может требовать субмикросекундных времен) и только потом рассматривать процесс формирования зародышей сажевых частиц.

В настоящей расчетной ячейке пропорции смеси подобраны таким образом, чтобы общее количество атомов углерода в молекулах ацетилена было равно общему количеству атомов углерода в ПАУ (по 800 в обеих подгруппах). Суммарное число атомов в расчетной ячейке составляло 1600 для углерода и 1166 для водорода, а плотность равна 0.01 г/см³.

Во время моделирования в системе поддерживалась постоянная температура Т = 2000 К. Выбор канонического ансамбля (NVT) обусловлен тем, что из-за происходящих экзотермических реакций моделируемая система может существенно отклоняться от начальных условий. Полное время моделирования составляло 7 нс, по истечении которого в расчетной ячейке была сформирована частица, насчитывающая 620 атомов углерода и 227 атомов водорода. На рис. 2 представлена динамика роста частицы, а также количество пяти-, шести- и семичленных колец в ее составе. По мере роста частицы соотношение [С]/[H] увеличивается, достигая в конце значения 2.7, которое лишь немного выходит за пределы диапазона, характерного для молодой сажи ([C]/[H] = 1.4–2.5 [43]).

Рис. 2.

Процесс формирования молодой сажевой частицы: (а) ‒ зависимость количества атомов углерода (1) и водорода (2) в составе частицы от времени; (б) ‒ количество углеродных колец в составе частицы с течением времени: 3 – гексагоны, 4 – пентагоны, 5 – гептагоны.

Модель “взрослой” сажевой частицы. В то время как молодые сажевые частицы характеризуются достаточно разупорядоченной аморфной структурой, взрослая сажа является более графитированной и преимущественно состоит из упорядоченных графеноподобных слоев. Среднее межслойное расстояние уменьшается по мере роста и старения частицы, достигая 3.3–3.5 Å [44], близких к аналогичным значениям для пиролитического графита. Для имитации таких графеноподобных слоев рассмотрены плоскости, замощенные молекулами циркумкоронена (рис. 3). Параметры такой системы соответствуют параметрам более взрослой сажи: среднее межслойное расстояние составляет d ≈ 3.5 Å, а плотность ρ ≈ 1.9 г/см³ [45]. В недавней работе было показано, что импульсный нагрев сажевых частиц до температур Т = 2500‒3500 К также может приводить к повышению их степени графитированности [46].

Методика моделирования соударений с молекулами газовой фазы. Для моделирования процесса взаимодействия молекулы ацетилена с поверхностью сажевой частицы выбрана ячейка постоянного объема с периодическими граничными условиями и размерами 10 ×10 ×10 нм³ для “молодой” частицы и 5 × 1.5 × 3 нм³ для “взрослой”. В ходе расчетов по направлению центра масс частицы запускалась молекула С₂Н₂ (рис. 4). Чтобы исключить влияние термостата на процесс взаимодействия между С₂Н₂ и сажевой частицей, расчеты производились в микроканоническом ансамбле (полная энергия системы оставалась постоянной). При этом как сами сажевые частицы, так и налетающие молекулы С₂Н₂ были предварительно выведены на равновесную температуру Т = 2000 К в ходе релаксации длительностью 0.1 нс. Полная кинетическая энергия молекулы С₂Н₂ в системе отсчета неподвижной частицы составляла KE = = 107.6 кДж/моль. При этом поступательная скорость движения центра масс С₂Н₂, выбранная в соответствии с максимумом распределения Максвелла при Т = 2000 К, равна $v = \sqrt {\frac{{2RT}}{M}} = 1130$ м/с, где M – масса молекулы. Остальная часть KE распределялась по внутренним степеням свободы молекулы. Для набора статистики для каждого типа поверхности (“молодой” и “взрослой”) моделировалось 400 актов столкновения, различающихся начальным положением молекулы ацетилена (рис. 4).

Рис. 3.

Модель поверхности взрослой сажевой частицы, состоящая из слоев плотно расположенного циркумкоронена.

Рис. 4.

Схема выполнения расчетов с молодой (а) и взрослой (б) частицами; стрелки – векторы направления поступательной скорости молекулы С₂Н₂.

РЕЗУЛЬТАТЫ И ОБСУЖДЕНИЕ

В описанных выше молекулярно-динамических расчетах вычислялись моментальные расстояния между ближайшими атомами молекулы ацетилена и сажевой частицы d_min, а также длительность контакта между молекулой и частицей. Предполагалось, что частица и налетающая молекула находятся в контакте при выполнении условия d_min < 7 Å.

Для молодой сажевой частицы при этом наблюдалось три основных сценария взаимодействия ацетилена с поверхностью: 1) моментальный отскок, 2) образование короткоживущего невалентного комплекса, 3) ковалентное присоединение молекулы к частице (рис. 5). Для взрослой же частицы ни одного акта образования ковалентных связей в рамках проведенных расчетов зафиксировано не было, и все сценарии взаимодействия сводились либо к отскоку, либо к образованию невалентного комплекса.

Рис. 5.

Изменение расстояния d_min между С₂Н₂ и поверхностью сажевой частицы (сплошные кривые) в ходе столкновения для трех характерных сценариев: (а) ‒ моментальный отскок, (б) ‒ образование короткоживущего невалентного комплекса, (в) ‒ образование ковалентной связи между С₂Н₂ и сажевой частицей; стрелки – момент образования связи; штриховые прямые – период контакта молекулы и частицы; на вставке ‒ момент образования ковалентной связи увеличено.

На рис. 6 приведено распределение плотности вероятностей длительности контакта между налетающей молекулой ацетилена и обоими типами рассматриваемых сажевых частиц. Распределения имеют схожий вид для времен до 2 пс, т.е. для сценариев, когда происходит практически моментальный отскок (рис. 5а). Однако распределение для молодой частицы имеет гораздо больший “хвост” на длительных временах (времена контакта более 5 пс). Иными словами, для молодой частицы, в отличие от взрослой, наблюдается ряд ситуаций, в которых налетающая молекула проводит вблизи поверхности частицы более 5 пс. Это объясняется попаданием молекулы ацетилена в карманы на поверхности частицы сложной формы (см. рис. 1в). Экспериментальные данные HRTEM свидетельствуют о том, что для молодых частиц характерные размеры таких межслойных пустот могут превышать 0.5 нм [4], что также наблюдается и в настоящих расчетах. Молекулы газовой фазы, попав в такие межслойные пустоты, оказываются в условиях так называемого наноконфайнмента – в замкнутом или в почти замкнутом наноразмерном пространстве. Известно, что кинетика химических реакций в таких условиях может существенно изменяться [47].

Рис. 6.

Распределение длительностей контакта между молекулой ацетилена и молодой (а) и взрослой (б) частицами.

Данный сценарий образования ковалентной связи в результате захвата в наноразмерную полость представлен на рис. 5в. В интервале времени t = 0–5 пс молекула летит в направлении частицы, а в интервале t  = 6–9 пс находится в контакте с ее поверхностью, но взаимодействует с ней исключительно невалентно. И только в момент t ≈ 9 пс между молекулой и частицей образуется ковалентная связь: расстояние d_min падает до 1.5–1.6 Å, что соответствует типичной длине С‒С-связи.

ЗАКЛЮЧЕНИЕ

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

Работа выполнена при поддержке Российского научного фонда (соглашение № 20-79-00245).