Журнал физической химии, 2021, T. 95, № 8, стр. 1225-1231

МЕТОДОЛОГИЧЕСКАЯ ЧАСТЬ

Катион ${{{\text{H}}}_{{\text{3}}}}{\text{SO}}_{4}^{ + }$, образующийся при автопротолизе, имеет структуру с равнозначными связями S–O–H вида [1, 4]:

Согласно [3], при изомеризации возникает катион вида:

При отщеплении от него молекулы воды образуется катион ${\text{HSO}}_{3}^{ + }$.

Таким образом, в качестве объектов моделирования выбраны:

• исходный катион (${{{\text{H}}}_{{\text{3}}}}{\text{SO}}_{4}^{ + }$);

• его изомер (H₂O–SO₂OH⁺);

• переходное состояние, соединяющее их (TS);

• отдельные частицы ${\text{HSO}}_{3}^{ + }$ и H₂O.

В качестве функционалов применены наиболее часто используемые из каждой группы: PBE [5] (GGA функционалы), M-06L [6] (meta-GGA- функционалы), B3LYP5 [7], PBE0 [8] (глобальные гибридные GGA-функционалы), M06, M06-2X [9] (глобальные гибридные meta-GGA-функционалы), LC-BLYP [10], ωB97x [11] (пространственно-разделенные гибридные функционалы) и B2T-PLYP [12] (дважды гибридные функционалы).

Для учета сил дисперсионного взаимодействия выбрана коррекция Гримме по методам D3BJ [131] (PBE, B3LYP5, PBE0, ωB97x), D3Zero [14] (M-06L, M06, M06-2X) и D2 [15] (LC-BLYP).

Во всех случаях использованы базисные наборы семейства Даннинга с добавлением диффузионных функций aug-cc-pVXZ [16, 17]. Для методов DFT применены наборы до трехэкспоненциальных включительно (X = D, T), для методов теории возмущений и связанных кластеров – до пятиэкспоненциальных (X = D, T, Q, 5).

При использовании метода связанных кластеров выбрана его реализация DLPNO-CCSD(T) [18]. Также использовано приближение к полному базисному набору по методу F-12 [19] (DLPNO-CCSD(T)-F12D) в сочетании с базисом cc-pVDZ-F12 [20].

Учет влияния растворителя осуществлен по континуальной модели COSMO [21]. В качестве растворителя рассмотрена безводная серная кислота в оценочных расчетах и, дополнительно, 5 и 15 (мас.)%-ные растворы SO₃ в Н₂SO₄ при исследовании механизма образования ${\text{HSO}}_{3}^{ + }$. Диэлектрическая проницаемость растворителя составляла: ε = 101 (0% SO₃), 259 (5% SO₃), 359 (15% SO₃) [22]. Эффективный радиус растворителя принят равным 1.3 Å для всех концентраций [23].

На уровне теории CCSD(T) проведены только одноточечные расчеты, геометрии для которых получены методом MP2/aug-cc-pV5Z.

Все расчеты осуществлены в программном комплексе ORCA [24], визуализация расчетов – в программном комплексе wxMacMolPlt [25].

Для оценки точности метода использованы два параметра: значение энергии изомеризации и длина связи S1–O5 в изомере H₂O–SO₂OH⁺. В качестве эталона сравнения использованы значение энергии, полученное на уровне теории DLPNO-CCSD(T)/aug-cc-pV5Z, и геометрия, полученная на уровне теории MP2/aug-cc-pV5Z. В качестве критерия адекватности метода принято выполнение условий “химической точности” (погрешность менее 1 ккал/моль ≈ 4.19 кДж/моль) и наименьшая погрешность расчета длины связи S1–O5.

Связь S1–O5 выбрана в качестве критерия оценки по причине сильной зависимости ее длины от метода расчета, как будет показано ниже.

Оценка энергии проведена для стадии изомеризации катиона ${{{\text{H}}}_{{\text{3}}}}{\text{SO}}_{4}^{ + }$:

для стадии отщепления молекулы воды: и суммарной реакции:

Дополнительно оценены профили энтальпии и свободной энергии Гиббса для реакций:

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

Геометрические параметры исходного катиона и его изомера приведены на рис. 1.

Рассчитанное расстояние S1–O5 (1.922 Å в вакууме) существенно больше длины ковалентных одинарных связей S–O (1.665 Å для пиросерных кислот [26], 1.54 Å для серной кислоты [27]). Можно предположить, что данная связь образована за счет электростатического взаимодействия частичных зарядов на атомах серы и кислорода. Значение ее длины свидетельствует о том, что эта связь слабее, чем ковалентная, следовательно, можно ожидать ее относительно легкого разрыва с образованием катиона ${\text{HSO}}_{3}^{ + }$ и молекулы воды.

При учете растворителя геометрия исходного катиона сохраняется (рис. 1), для изомера наиболее существенное изменение наблюдается в длине связи S1–O5 (уменьшение на 0.145 Å).

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

Из табл. 1 видно, что при использовании метода теории возмущений второго порядка даже применение широких базисов не позволяет получить оценку энергии процесса с приемлемой точностью. Применение методов связанных кластеров позволяет приблизиться к точной оценке энергии на базисах не ниже QZ. Использование малых базисов приводит к завышенному значению длины связи S1–O5 и заниженному значению энергии изомеризации.

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

Из табл. 2 видна общая тенденция уменьшения ошибок при переходе к более широкому базису. В целом, можно выделить лучшие результаты для пространственно-разделенных функционалов: ωB97x попадает в условие “химической точности”, LC-BLYP близок к нему. Функционал M06-2X демонстрирует изменение знака ошибки оценки энергии, что может указывать на случайное совпадение полученного значения с эталонным. Из остальных функционалов погрешность, близкую к допустимой, демонстрирует только PBE0.

Для сравнения приведены значения энергии изомеризации ΔE₁, рассчитанные в работе [3] различными методами, кДж/моль:

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

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

По результатам, представленным в табл. 3, можно сделать вывод, что при учете растворителя погрешность практически всех методов снижается. Для методов связанных кластеров достаточным становится использование трехэкспоненциального базисного набора. Исключение составляет метод DLPNO-CCSD(T)-F12D, погрешность которого, наоборот, сильно возросла. Методы теории возмущений по-прежнему не обеспечивают требуемую точность расчета.

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

Для методов теории функционала плотности также наблюдается тенденция к снижению ошибки энергии при переходе к растворителю. Для трехэкспоненциального базисного набора погрешность в пределах допустимой достигается при использовании функционалов PBE0, M06-2X, LC-BLYP и ωB97x. Однако увеличение погрешности на функционале M06-2X при переходе к более широкому базису также может рассматриваться в пользу скорее случайного попадания ошибки в требуемый диапазон.

Таким образом, имеется один функционал, удовлетворяющий условиям адекватности как для вакуума, так и для растворителя – ωB97x, который был использован для дальнейших расчетов.

Следующий этап – расчет геометрии (рис. 2) и энергии переходного состояния, а также отдельных частиц ${\text{HSO}}_{3}^{ + }$ и H₂O (табл. 5).

Погрешность метода ωB97x/aug-cc-pVTZ при расчете энергий ΔE^≠ и ΔE₂ в вакууме оказалась несколько выше требуемой (5.04 и 6.62 кДж/моль), однако погрешность расчета суммарной энергии ΔE₃ (3.69 кДж/моль) не превышает допустимой. Метод MP2/aug-cc-pV5Z не обеспечивает требуемой точности расчетов, кроме оценки ΔE₂. Суммарная энергия образования катиона ${\text{HSO}}_{3}^{ + }$, рассчитанная с использованием функционала ωB97x (229.75 кДж/моль), достаточно хорошо коррелирует с экспериментальными данными (225 ± ± 18 кДж/моль). В работе [3] приведены следующие значения энергии ΔE^≠ газофазной реакции, кДж/моль:

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

Учет растворителя снижает энергетические барьеры реакции изомеризации (ΔE^≠) на ~9 кДж/моль, реакции отрыва молекулы воды (ΔE₂) на 18 (CCSD(T)) и 12 (ωB97x) кДж/моль, суммарной реакции (ΔE₃) на 32 (CCSD(T)) и 33.5 (ωB97x) кДж/моль. Однако изменение концентрации SO₃ в растворителе мало влияет на энергетические характеристики, для всех характеристик разность энергии при варьировании концентрации составляет менее 0.5 кДж/моль.

Для метода связанных кластеров расчет гессиана не проводился, и термодинамические характеристики оценивались через электронную энергию CCSD(T) и термохимические поправки по методу MP2. Оценка термодинамических характеристик исследуемых реакций приведена в табл. 6 и 7.

Согласно данным табл. 6, учет термодинамических характеристик несколько снижает энергию активации реакции (1) (ΔH^≠) и суммарный тепловой эффект (ΔH₃), однако расчетные данные по-прежнему коррелируют с экспериментальными.

Сравнение значений изменения энергии Гиббса и энтальпии реакции (табл. 6, 7) показывает, что энтропийный фактор оказывает существенное влияние на реакцию (2) (и соответственно суммарную реакцию), для реакции (1) он незначителен.

Энтропия суммарной реакции положительна, следовательно, увеличение температуры должно смещать равновесие в сторону образования катиона ${\text{HSO}}_{3}^{ + }$. Для оценки данного влияния была рассчитана величина ΔG₃ в диапазоне температур 25–100°С (298–373 K) (рис. 3).

Наибольшие значения константы равновесия суммарной реакции для выбранного диапазона температур равны 1.6 × 10^–23 (вакуум) и 8.3 × 10^–19 (растворитель), что крайне мало для обеспечения сколь-нибудь заметной концентрации катиона ${\text{HSO}}_{3}^{ + }$ при данных температурах.

Учет влияния растворителя в явном виде [28] показывает, что в олеуме частицы ${\text{HSO}}_{3}^{ + }$ существуют преимущественно в виде ассоциатов ${\text{HSO}}_{3}^{ + }$ ⋅ H₂SO₄ (точнее HO₃S⋅⋅⋅O–SO₃H₂, что можно рассматривать и как протонированную форму пиросерной кислоты). Отсюда можно сделать вывод, что катион ${\text{HSO}}_{3}^{ + }$ в свободном виде при температурах ниже 100°С не будет оказывать какого-то заметного влияния на процессы сульфирования, несмотря на то, что он – более активный сульфирующий агент по сравнению с SO₃ [29].

Таким образом, проведено сравнение адекватности методов теории функционала плотности и теории возмущений для исследования реакции изомеризации катиона ${{{\text{H}}}_{{\text{3}}}}{\text{SO}}_{4}^{ + }$, выбран наиболее точный метод (ωB97x/aug-cc-pVTZ). Методами CCSD(T)/aug-cc-pV5Z, MP2/aug-cc-pV5Z и ωB97x/aug-cc-pVTZ проведено исследование реакции образования катиона ${\text{HSO}}_{3}^{ + }$ в вакууме и при неявном учете растворителя по модели COSMO.

Произведен расчет термодинамических характеристик реакции образования катиона ${\text{HSO}}_{3}^{ + }$ и показано, что при температурах ниже 100°С его концентрация как в вакууме, так и в растворителе слишком мала, для того, чтобы рассматривать его как истинный сульфирующий агент.