Моделирование реактора парового риформинга в Aspen Plus с использованием кинетики Xu & Froment

Настоящий отчет представляет собой исчерпывающее, пошаговое руководство по созданию цифрового двойника реактора парового риформинга метана (Steam Methane Reforming, SMR) типа GHR (Gas-Heated Reformer) в среде Aspen Plus. Обратите внимание: данный материал предоставляется исключительно в образовательных и информационных целях. Параметры кинетического моделирования являются теоретическими и не составляют профессиональную инженерную рекомендацию; перед любым промышленным масштабированием, проектированием или изменением операционных режимов они должны быть строго валидированы на основе закрытых спецификаций поставщиков катализаторов и эмпирических данных пилотных установок.

Интеграция кинетической модели Ленгмюра-Хиншелвуда-Хоугена-Уотсона (LHHW), первоначально разработанной Xu и Froment в 1989 году, требует глубокой математической адаптации для корректной работы решателей Aspen Plus [cite: 1]. Данное руководство описывает не только алгоритм действий («Click-by-Click»), но и фундаментальную физико-химическую логику, стоящую за каждым вводимым параметром. Особое внимание в отчете уделено преодолению главного ограничения интерфейса Aspen Plus при работе с данной кинетикой — проблеме деления на ноль из-за присутствия парциального давления водорода в знаменателях кинетических уравнений [cite: 2]. Представленные ниже расчеты и настройки переведены в строгий базис СИ (Паскали), что является стандартом по умолчанию для блока LHHW в Aspen Plus, исключая ошибки размерностей на несколько порядков [cite: 3]. Отчет сформирован с учетом необходимости интеграции высших углеводородов, расчета внутридиффузионного торможения и проектной оптимизации под заданный технологический критерий (мольная доля метана $\le$ 0.4%).

Executive Summary (Управляющее резюме)

Для успешного моделирования GHR-реактора с целевым содержанием $CH_4 \le 0.4\%$ в Aspen Plus, инженер должен реализовать следующие ключевые решения:

• Термодинамика: Использование уравнения состояния Peng-Robinson (PENG-ROB) для строгого расчета фугитивностей неидеальной газовой смеси при высоких температурах и давлениях.
• Управление сырьем (C2-C6): Интеграция каскада из адиабатического реактора минимизации энергии Гиббса (RGibbs) для симуляции предриформинга и блока нагревателя (Heater) для точной фиксации входной температуры трубного пространства (368.8 °C).
• Кинетическая адаптация: Строгий пересчет предэкспоненциальных множителей (фактор $A$) из базиса bar в базис Pa и математический перенос экспоненты водорода из знаменателя оригинального уравнения в числитель с отрицательной степенью во вкладках Aspen Kinetic Factor и Adsorption Expression.
• Аппаратная конфигурация (RPLUG): Использование модели противоточного многотрубного реактора с обязательным введением фактора эффективности катализатора ($\eta = 0.05 - 0.1$) для компенсации жесткого внутридиффузионного торможения в промышленных пеллетах.

---

1. ПОДГОТОВКА СИСТЕМЫ И ТЕРМОДИНАМИКА

Первый этап моделирования закладывает фундамент для физико-химических свойств смесей. Неправильный выбор термодинамического пакета или некорректная обработка сырья приведут к нефизичным результатам конверсии и тепловых балансов.

1.1. Выбор компонентов и обработка высших углеводородов (C2-C6)

Оригинальная кинетическая модель Xu & Froment (1989) описывает взаимодействие исключительно метана и пара на никелевом катализаторе [cite: 4]. При температурах риформинга высшие углеводороды (этан, пропан, бутан и др.) подвергаются мгновенному, практически необратимому крекингу и предриформингу [cite: 5, 6].

Если направить поток, содержащий $C_2-C_6$, напрямую в блок кинетического реактора (RPLUG) с настроенной кинетикой только для метана, Aspen Plus проигнорирует эти компоненты (они выступят как инертные газы), что приведет к грубейшей ошибке в расчете материального и теплового балансов. Попытка использовать стехиометрический реактор (RStoic) для этих целей является ошибочной практикой: RStoic требует точного задания стехиометрии и фракционной конверсии для каждой реакции крекинга, которые нелинейно зависят от состава и температуры.

Практическое решение (Best Practice): Интеграция каскада предриформинга Для корректного учета высших углеводородов необходимо смоделировать стадию адиабатического предриформинга с последующей температурной коррекцией до основного аппарата GHR. В адиабатическом предриформинге паровой риформинг тяжелых углеводородов протекает как сильно эндотермический процесс, в то время как сопутствующие реакции сдвига фаз (WGS) и метанирования экзотермичны. Для классического природного газа баланс этих реакций формирует специфический профиль: изначальное падение температуры смеси сменяется ее небольшим ростом, что в итоге приводит к общему снижению температуры на выходе из слоя по сравнению с температурой на входе [cite: 5, 7, 8].

1. Перейдите во вкладку Properties -> Components -> Specifications.
2. Введите список компонентов (Component ID): CH4, H2O, CO, CO2, H2, N2 (если есть в сырье), а также весь спектр высших углеводородов (C2H6, C3H8, NC4H10, IC4H10 и т.д.).
3. В среде Simulation, перед блоком основного реактора RPLUG, установите реактор типа RGibbs (назовите его PRE-REF). Выбор RGibbs элегантно решает проблему сложной кинетики C2-C6 путем минимизации свободной энергии Гиббса, полностью симулируя адиабатический 100% распад высших углеводородов [cite: 6].
4. Настройте RGibbs:
   • Вкладка Setup: укажите Adiabatic (Duty = 0). Давление задайте равным давлению сырья (45.9 кгс/см² изб).
   • Вкладка Products: выберите опцию Calculate phase equilibrium and chemical equilibrium.
   • Вкладка Restricted Eq: выберите компоненты $C_2-C_6$ и спецификацию конверсии равную 100%.
5. Выходящий поток из PRE-REF теперь будет содержать только равновесную смесь CH4, H2O, CO, CO2 и H2. Из-за эндотермической природы распада природного газа, температура потока упадет относительно исходной [cite: 8].
6. Температурная корректировка: Поскольку исходные данные строго регламентируют вход в трубное пространство GHR при температуре 368.8 °C, установите блок теплообменника типа Heater (назовите его ADJ-TEMP) сразу после PRE-REF. В настройках Heater задайте Temperature = 368.8 °C и Pressure = 45.9 кгс/см² изб. Именно этот стабилизированный поток направляется в трубное пространство RPLUG.

1.2. Выбор термодинамического пакета (Property Method)

Для процессов SMR, протекающих при высоких температурах (360–1000 °C) и давлениях (30–50 бар), газовая смесь ведет себя неидеально, однако отклонения связаны преимущественно с летучестью и фугитивностью компонентов при высоких давлениях [cite: 4, 9].

Обоснование метода: Уравнение состояния Пенга-Робинсона (Peng-Robinson, PENG-ROB) является индустриальным стандартом для углеводородных систем и синтез-газа при таких параметрах [cite: 10, 11]. Метод PRMHV2 (Peng-Robinson with Modified Huron-Vidal mixing rules) часто рекомендуется для смесей с высоким содержанием полярной воды, однако при температурах выше критической точки воды (374 °C) и низких плотностях газа классический PENG-ROB обеспечивает превосходную точность без излишней вычислительной жесткости.

Пошаговая настройка:

1. Перейдите во вкладку Properties -> Methods -> Specifications.
2. В поле Base method выберите PENG-ROB.
3. Перейдите во вкладку Parameters -> Binary Interaction -> PRKBV-1, чтобы убедиться, что бинарные коэффициенты взаимодействия (особенно для пары H2O-CH4 и H2O-CO2) загружены из баз данных Aspen.

2. ДЕТАЛЬНЫЙ АЛГОРИТМ ВВОДА КИНЕТИКИ LHHW

Кинетика Xu & Froment (1989) описывается тремя обратимыми реакциями [cite: 12]:

1. SMR 1: $CH_4 + H_2O \rightleftharpoons CO + 3H_2$
2. WGS (Сдвиг фаз): $CO + H_2O \rightleftharpoons CO_2 + H_2$
3. SMR 2 (Глобальная): $CH_4 + 2H_2O \rightleftharpoons CO_2 + 4H_2$

Математическая структура модели LHHW в Aspen Plus требует приведения уравнения скорости к следующему виду [cite: 1, 2]: $r = \frac{(\text{Kinetic Factor}) \times (\text{Driving Force})}{(\text{Adsorption Term})^n}$

2.1. Решение проблемы деления на ноль ($p_{H_2}$)

В оригинальной статье скорости реакций имеют вид $r_1 = \frac{k_1}{p_{H_2}^{2.5}} (...)$, что вызывает фатальную ошибку деления на ноль в начале трубы, если в сырье нет водорода [cite: 2]. Меры по устранению:

1. В интерфейсе мы зададим отрицательные экспоненты для $H_2$ напрямую во вкладке Kinetic Factor.
2. Аналогичная проблема возникает в адсорбционном члене: слагаемое $K_{H_2O} \cdot \frac{p_{H_2O}}{p_{H_2}}$ требует отрицательной степени для водорода во вкладке Adsorption Expression [cite: 9].
3. Обязательное условие: В сырьевой поток, подаваемый в систему, необходимо искусственно добавить небольшое количество $H_2$ (если оно отсутствует после предриформера), чтобы обеспечить инициализацию расчета и избежать сингулярности в первой ячейке реактора [cite: 2]. Обычно достаточно мольной доли ~0.001 (0.1%).

2.2. Настройка блока Reactions

Перейдите в Simulation -> Reactions -> Reactions. Создайте новую реакцию (назовите XF-KINET), выберите тип LHHW. Нажмите New и добавьте три реакции (R1, R2, R3) согласно стехиометрии. Обязательно укажите тип Reversible для всех трех.

А) Вкладка Kinetic Factor (Адаптация размерностей)

Aspen Plus по умолчанию использует базис парциального давления в Паскалях (Pa). Оригинальные константы из статьи рассчитаны для базиса давления в bar и скорости реакции в $kmol / (kg_{cat} \cdot h)$. Если ввести оригинальные константы $A_{bar}$ в интерфейс без пересчета, скорости реакций будут занижены на порядки, и симуляция покажет нулевую конверсию [cite: 3]. Расчет предэкспоненциального фактора $A$ для базиса СИ (Pa) выполняется по формуле пересчета, учитывающей степени парциальных давлений в числителе и знаменателе для каждой реакции, а также перевод часов в секунды (деление на 3600) [cite: 13].

Ниже представлена детальная таблица сравнения констант для обоих базисов. В интерфейс Aspen Plus необходимо вводить значения исключительно из колонки "Базис Pa".

Ввод в Aspen: Во вкладке Kinetic, выберите реакцию. В секции Factor введите $A$ (из колонки Pa) и $E$. Убедитесь, что единицы энергии активации выбраны как J/mol. В секции Concentration basis выберите компонент H2 и укажите соответствующую отрицательную степень (Exponent).

Б) Вкладка Driving Force

Движущая сила в Aspen рассчитывается как $(\text{Term 1} - \text{Term 2})$.

• Term 1 соответствует прямым компонентам (константа = 1, параметры $A=0, B=0, C=0, D=0$).
• Term 2 соответствует обратным компонентам. Константа для Term 2 равна $1/K_{eq}$. Aspen вычисляет ее как $\exp(A + B/T)$. Поскольку нам нужно ввести именно $1/K_{eq}$, знаки $A$ и $B$ инвертированы по отношению к стандартным константам химического равновесия [cite: 13].

Таблица ввода для Term 1 (одинакова для базисов Pa и bar):

Таблица ввода для Term 2 (Базис Паскали - Pa):

Синтез логики: Параметр $A$ в базисе Паскалей включает в себя поправку $\ln(10^{10})$ для реакций с изменением числа молей [cite: 13]. Это критически важный шаг; ошибка здесь приведет к сдвигу равновесия (driving force) в неправильную сторону.

В) Вкладка Adsorption Expression

Знаменатель в кинетике Xu & Froment имеет вид: $DEN = 1 + K_{CO}p_{CO} + K_{H_2}p_{H_2} + K_{CH_4}p_{CH_4} + K_{H_2O}\frac{p_{H_2O}}{p_{H_2}}$

В интерфейсе Aspen Plus это вводится в виде суммы слагаемых (Term 1 ... Term 5). Константы адсорбции $K_i$ описываются уравнением $\ln(K_i) = A + B/T$ [cite: 13]. Во вкладке Adsorption Expression задайте Adsorption Exponent = 2 для всех трех реакций.

Обратите внимание: Для Term 5 компонент H2 имеет степень -1, что отражает дробь $p_{H_2O}/p_{H_2}$ [cite: 9]. Термодинамически константа адсорбции воды в этой формуле безразмерна, поэтому параметр $A$ остается неизменным при переходе от бар к Паскалям (12.084) [cite: 13].

3. ТОЧНАЯ КОНФИГУРАЦИЯ БЛОКА RPLUG

Переходим к конфигурации самого аппарата теплообменного риформинга (GHR). Здесь кинетика встречается с массо- и теплопереносом.

3.1. Вкладка Setup

• Reactor Type: Выберите Reactor with thermal fluid. Это позволяет моделировать теплообмен через стенку трубы.
• Fluid direction: Выберите Countercurrent (противоток). Конвертированный газ (теплоноситель) традиционно подается снизу-вверх, нагревая опускающийся по трубам сырьевой газ.

3.2. Вкладка Configuration

• Multitube Reactor: Поставьте галочку Multitube reactor.
• Number of tubes: Для начала задайте ориентировочное количество труб промышленного масштаба, например, 500 (этот параметр мы будем оптимизировать позже).
• Length / Diameter: Введите характерные для SMR значения (например, Length = 10 m, Diameter = 0.1016 m [cite: 14]).

3.3. Вкладка Catalyst

Риформинг метана протекает исключительно на поверхности катализатора, поэтому корректное задание плотности загрузки критично для расчета истинного времени контакта.

• Phase: Установите Vapor-Only (в аппарате нет жидкости).
• Particle density (Плотность частиц): Стандартное значение для никель-алюминиевых катализаторов составляет около 2350 kg/m3.
• Bed voidage (Порозность слоя): Введите значение 0.45 - 0.5 [cite: 11]. Aspen Plus автоматически вычислит Bulk density (Насыпную плотность) как $\text{Particle density} \times (1 - \text{Bed voidage})$.

3.4. Вкладка Pressure

Для точного моделирования распределения давлений (от которого напрямую зависит кинетика), активируйте расчет перепада давления.

• В секции Pressure drop выберите Ergun Equation.
• Particle diameter: Введите диаметр частиц катализатора, типичный для колец Рашига или цилиндров, например, 0.015 m (15 мм) [cite: 14].

3.5. Вкладка Thermal Fluid

Теплоноситель в межтрубном пространстве GHR обеспечивает эндотермические реакции теплом.

• Specification: Выберите T-P-Flow-Frac или укажите отдельный поток теплоносителя, если он связан с остальной схемой.
• Температура: 1005 °C (исходные данные).
• Давление: 38.5 kgf/cm2 (убедитесь, что единицы выбраны верно — Aspen имеет опции kg/sqcm для kgf/cm²).
• Heat Transfer (U): Коэффициент теплопередачи от межтрубного к трубному пространству. Задайте Constant значение в диапазоне 60 - 100 W/m2-K.

3.6. Вкладка Kinetics / Effectiveness Factor

В кинетике Xu & Froment получены внутренние (intrinsic) константы скорости на измельченном порошке [cite: 15]. Однако промышленные пеллеты катализатора SMR (диаметром 15 мм) подвержены мощнейшему внутридиффузионному торможению, при котором реакция протекает только в тонком макропористом поверхностном слое.

Если не учесть этот фактор в симуляции, реактор завершит конверсию на первых 50 сантиметрах трубы. В промышленной практике эффективность гранулы (Effectiveness Factor) для парового риформинга эмпирически оценивается всего в 5–10% [cite: 15, 16].

• Во вкладке Kinetics добавьте созданный ранее набор реакций XF-KINET.
• В столбце Effectiveness factor выберите Constant для всех трех реакций.
• Задайте значение $\eta$ в диапазоне 0.05 – 0.1 (рекомендуемое стартовое значение: 0.07). Это означает, что эффективно работает лишь 7% физической массы катализатора [cite: 15].

4. ПРОЕКТНАЯ ОПТИМИЗАЦИЯ И СХОДИМОСТЬ

Цель моделирования — достичь жесткого целевого показателя: остаточная мольная доля $CH_4 \le 0.4\%$ на выходе трубного пространства. Нарушение этого критерия несет тяжелые технологические последствия: непрореагировавший метан начинает выступать как балластный инерт в последующих контурах синтеза (например, аммиака или метанола), что катастрофически снижает парциальное давление активных реагентов, вынуждает увеличивать сброс продувочного газа (purge rate) и обрушивает общую энергоэффективность завода [cite: 16].

4.1. Настройка Design Specs (Проектные спецификации)

Мы поручим решателю Aspen Plus автоматически подобрать количество труб (или их длину), чтобы достичь требуемой конверсии.

1. В навигаторе перейдите в Flowsheeting Options -> Design Specs и создайте новую спецификацию (например, CH4-TARG).
2. Вкладка Define: Создайте переменную XCH4.
   • Type: Mole-Frac
   • Stream: Выберите поток, выходящий из RPLUG.
   • Component: CH4.
3. Вкладка Spec:
   • Spec: XCH4
   • Target: 0.004 (это 0.4% в мольных долях).
   • Tolerance: 0.0001.
4. Вкладка Vary:
   • Type: Block-Var
   • Block: RPLUG
   • Variable: Выберите NTUBES (количество труб).
   • Limits: Задайте разумные границы (например, от 100 до 1500).

4.2. Анализ профилей реактора (Отклонение от Sensitivity)

Модуль Sensitivity в Aspen Plus отлично подходит для варьирования граничных условий (например, расхода сырья или температуры), однако он не способен нативно строить профили вдоль длины (координаты Z) внутри самого реактора без написания сложных Fortran-вставок (Calculator blocks), извлекающих массив данных B_RPLUG_Z.

Практическое решение (Использование Profiles): Для визуализации профилей температуры, давления, конверсии $CH_4$ и выхода $H_2$ по длине реактора, необходимо использовать встроенный инструмент блока RPLUG:

1. Убедитесь, что расчет сошелся без ошибок (Status: Results Available).
2. В дереве навигации разверните блок вашего реактора (RPLUG) и выберите папку Profiles.
3. Перейдите во вкладку Plot.
4. В разделе X-axis выберите Length.
5. В разделе Y-axis отметьте галочками Temperature (для оценки теплопереноса) и Mole-Frac (выберите CH4 и H2 для оценки кинетики). Нажмите кнопку Plot. Это создаст точный график изменения параметров от 0 до 10 метров трубы.

4.3. Рекомендации по устранению ошибок сходимости (Convergence)

Кинетика LHHW для SMR является жесткой (stiff) математической системой дифференциальных уравнений. Возможные ошибки и пути решения:

• Сингулярность в начале интегрирования: Убедитесь, что сырье содержит следы $H_2$. Если ошибка повторяется, увеличьте долю $H_2$ до 1% за счет уменьшения $H_2O$.
• Сбой интегратора (Integration Error): Перейдите в блок RPLUG, вкладка Setup -> Integration. Измените метод интегрирования на Gear (специально для stiff-систем) и уменьшите Step size (например, Initial step = $10^{-5}$ м).
• Температурный кроссовер: Если температура трубного пространства превышает температуру теплоносителя, возникает нефизичный обратный теплоперенос. Проверьте коэффициент $U$ или увеличьте расход теплоносителя.

5. ЧЕК-ЛИСТ ВЕРИФИКАЦИИ МОДЕЛИ

После того как расчет завершен, инженер должен валидировать модель по 6 ключевым точкам проверки:

1. Проверка инициализации H2 (Отсутствие сингулярности): Проверьте лог решателя (Control Panel). Там не должно быть предупреждений о делении на ноль или обнулении знаменателя Adsorption term [cite: 2].
2. Температурный подход к равновесию (Approach to Equilibrium, AT): Рассчитайте равновесную температуру для состава на выходе. AT должен составлять от $5$ до $20 °C$. Если AT стремится к нулю в середине трубы, значит катализатор избыточен (или Effectiveness Factor завышен).
3. Энергетический баланс (Тепловая нагрузка): В результатах RPLUG сравните выделенное тепло теплоносителем (Shell Side Heat Duty) с поглощенным теплом реакционным потоком. Они должны совпадать с учетом знака.
4. Профиль температуры (Эндотермический шок): Откройте график из вкладки Profiles. На первых 10–20% длины трубы (входная зона) должно наблюдаться резкое падение температуры реакционной смеси на 30–60 °C. Это происходит потому, что скорость поглощения тепла химическими реакциями значительно опережает теплопередачу через стенку. В промышленности этот феномен жестко мониторится поверхностными термопарами (skin thermocouples) для предотвращения локальных переохлаждений [cite: 17]. После прохождения "дна" температура плавно возрастает.
5. Проверка перепада давления (Ergun validation): Откройте профиль давления. Падение давления в трубах $\Delta P$ должно составлять реалистичные $1–3$ бара (1-3 кгс/см²). Если падение превышает 10 бар, необходимо пересмотреть диаметр частиц или количество труб.
6. Отсутствие проскока высших углеводородов: Проверьте выходящий поток. Мольная доля $C_2-C_6$ должна быть строго равна нулю (благодаря корректной работе блока PRE-REF, предшествующего GHR).

6. ВЕКТОРЫ РАЗВИТИЯ МОДЕЛИРОВАНИЯ (FUTURE OUTLOOK)

Классический подход с адаптацией макрокинетики Xu & Froment (1989), описанный в данном руководстве, остается золотым стандартом для базового инженерного проектирования. Однако современная парадигма моделирования риформинга смещается в сторону двух инновационных направлений:

1. Микрокинетическое моделирование (Wei & Iglesia): Переход на фундаментальные механизмы, такие как кинетика Вей и Иглесиа, которые описывают реакции на уровне элементарных стадий активации C-H связей на кластерах рутения/никеля. Данные механизмы изначально лишены математической проблемы деления на ноль из-за отсутствия парциального давления водорода в знаменателе [cite: 16].
2. Интеграция CFD и Aspen Plus: Для точного учета радиальных градиентов температуры (которые в трубках SMR могут достигать 50 °C от стенки к центру) все чаще применяется со-симуляция. Одномерная кинетика Aspen Plus связывается с 3D-моделями тепло- и массопереноса в программных комплексах вычислительной гидродинамики (например, ANSYS Fluent), что позволяет отказаться от константного фактора эффективности катализатора в пользу прямого расчета диффузии в грануле [cite: 16].

Sources:

1. mdpi.com
2. aidic.it
3. youtube.com
4. acs.org
5. ecn.nl
6. tno.nl
7. csic.es
8. nlr.gov
9. nih.gov
10. combustion-institute.it
11. rjwave.org
12. lu.se
13. Link
14. aidic.it
15. ijret.org
16. ulisboa.pt
17. ucla.edu