CC BY
16
УДК 539.19
МОДЕЛИРОВАНИЕ РЕАКЦИЙ В ВОДНЫХ КЛАСТЕРАХ МЕТОДАМИ КВАНТОВОЙ ХИМИИ
Б. Л. Григоренко, М. А. Князева, А. В. Немухин
(кафедра физической химии)
Неэмпирическими методами квантовой химии исследовано влияние конечного числа молекул воды на энергетические профили реакций ОН" + С02 ^ НС03- и 2^0 ^ Н30+ + ОН-. Для молекул воды гидратных оболочек использовано представление эффективными фрагментами по схеме, реализованной в программе ОАМЕ88. В обоих случаях вдоль выбранной координаты реакции выполнена полная оптимизация геометрических параметров систем, включая и положения окружающих молекул воды.
Неэмпирические методы квантовой химии достаточно активно применяются при моделировании строения и динамики молекулярных систем, причем акценты все более смещаются в сторону расчетов элементарных химических реакций в конденсированных средах. По понятным причинам реакциям в водной среде уделяется особое внимание. Теоретические подходы к расчетам энергетических профилей реакций в полярных растворителях, в частности в воде, основываются либо на континуальных, либо на кластерных моделях. Данная работа принадлежит второму направлению, и молекулы растворителя рассматриваются явно в рамках квантово-механического описания системы.
В последние годы авторским коллективом известной программы вЛМБ88 [1] сделан серьезный шаг на пути моделирования явлений сольватации, а именно реализованы алгоритмы построения и учета в традиционных квантовохимических схемах так называемых эффективных фрагментов (ЭФ) [2]. Согласно этому подходу молекулы растворителя, в частности воды, включаются в анализируемую квантовую систему, но их влияние (как эффективных фрагментов) на частицы растворенного вещества передается посредством потенциалов, параметры которых определяются на предварительных стадиях по результатам неэмпирических квантовохимических расчетов. Потенциалы ЭФ включают в себя дальнодействующие мультиполь-ные и поляризационные вклады, а также вклады короткодействующего обменного отталкивания. Подобранные параметры для молекул воды позволили с высокой точностью передать структуру кластеров (И20)п вплоть до п = 20 [3].
Способность приближения ЭФ столь хорошо описывать межмолекулярные комплексы с водородными связями позволяет надеяться на успешное рассмотрение важных химических реакций в водных кластерах, и прежде всего на оценку влияния молекул воды из гидратных оболочек на энергетические профили реакций. Одним из наиболее узких мест при моделировании реакций в кластерах молекул растворителя является необходимость подстраивать их положения для каждого значения координаты реакции, т.е. многократно решать задачу об оптимизации геометрии системы с большим числом степеней свободы. При этом следует принять во внимание чрезвычайную трудоемкость квантовохимических расчетов взаимодействий молекул растворителя между собой, а также с час-
тицами растворенного вещества на достаточном уровне, требующем и базисов большой размерности, и адекватного учета эффектов электронной корреляции. Описание молекул растворителя эффективными фрагментами при полной оптимизации геометрических параметров всей системы вдоль пути реакции позволяет значительно ускорить расчеты энергетических профилей.
В данной работе приведены результаты расчетов двух реакций при явном включении молекул растворителя (воды). Рассмотрены достаточно важные процессы, протекание которых в газовой фазе принципиально отличается от условий водной среды,
OH-
+ CO2 ^ HCO3-
и 2H2O ^ H3O+ + OH-.
Основная цель сообщения - демонстрация возможностей подхода ЭФ для построения энергетических профилей реакций в кластерах. Все расчеты выполнены с использованием программы PC GAMESS [4] на вполне доступных рабочих станциях на основе Pentium III.
Реакция OH- + CO2 ^ HCO3- относится к классу наиболее значимых для биохимических приложений процессов, ответственных, в частности, за регулирование pH в крови и за перенос диоксида углерода в живых системах. Надежные квантово-химические расчеты [5, 6] показывают, что на энергетическом профиле пути реакции в отсутствии растворителя барьер отсутствует, в то время как экспериментальные исследования данной реакции в водном растворе [7] свидетельствуют о наличии активацион-ного барьера. Таким образом, именно водная среда ответственна за существенную модификацию энергетического профиля. В работе Пенга и Мерца [6] этот вывод был подтвержден модельными расчетами методами молекулярной динамики с эмпирическими потенциальными функциями. В данной работе описывается первый шаг на пути квантовохимического моделирования рассматриваемой реакции в водных кластерах.
На рис. 1 изображен рассчитанный энергетический профиль реакции OH- + CO2 ^ HCO3- в кластере из 30 молекул воды. В качестве координаты реакции выбрано расстояние R между кислородом гидроксил-аниона и углеродом молекулы CO2. Энергии и силы в подсистеме OH- + CO2 ^ HCO3- вычислены методом Хартри-Фока с двухэкспонентными базисами. Описание молекул воды как эффективных фрагментов проведено с параметрами, приведенными в работе [2]. Геометрическая конфигурация
молекул Н20 предполагается фиксированной с межъядерными расстояниями 0-Н, равными 0,9572 А, и углом Н—0—Н, равным 104,52°. Для каждого значения координаты реакции Я все остальные геометрические параметры системы, включая положения окружающих молекул воды, получены при минимизации полной энергии.
Рельеф рассчитанного потенциала показан на рис. 1. В полном согласии с данными прежних неэмпирических расчетов [5, 6] глобальный минимум отвечает продукту реакции НС03 при Я = 1,45 А. Глубина потенциальной ямы 44 ккал/моль, отсчитанная от энергии реагентов при Я = 6 А, также согласуется с литературными данными (48 ккал/моль без поправки на присутствие молекул растворителя [6]). Видно, что при движении вдоль пути реакции в водном кластере энергия сначала понижается, а затем проходит через барьер, т.е. результаты расчета качественно согласуются с экспериментальными данными [7] и данными моделирования методами молекулярной динамики [6].
Расположение частиц системы в области реагентов (Я = 6 А) и продукта (Я = 1,45 А), показанное на рис. 2, дает основание предложить более детальное качественное объяснение эффекта. В начале пути (рис. 2) молекулы растворителя занимают положение между реагирующими частицами ОН и С02 и до определенного момента сближение реагентов приводит к понижению энергии за счет образования эффективной сетки водородных связей. В дальнейшем принудительное уменьшение Я должно сопровождаться выталкиванием молекул воды из пространства между 0Н- и С02 (рис. 2, б), что и отвечает за энергетический барьер (рис. 1). Для более определенных заключений требуются дополнительные расчеты, в частности с варьированием числа молекул воды в кластере, а также с привлечением методов молекулярной динамики, однако уже первые результаты весьма интересны.
Второй пример, рассмотренный в данной работе, также относится к системе, вызывающей повышенное внимание, а именно к реакции, в результате которой из двух молекул воды образуется ионная пара Понятно, что без стабилизирующего влияния среды
Рис. 1. Энергетический профиль пути реакции 0Н- + С02 ^ НС03- в кластере 30 молекул воды вдоль координаты Я, соответствующей расстоянию от кислорода гидроксил-аниона до углерода. Для каждого значения Я геометрические параметры всей системы отвечают минимуму полной энергии
Рис. 2. Равновесная конфигурация системы: а - (0Н + С02) @ (Н20)30 для Я = 6 А (11,3 ат. ед.); б - (НС03-) @ (Н20)30 для Я = 1,45 А (2,7 ат. ед.). Для наглядности атомы кислорода реакционной подсистемы выделены темным цветом
подобная система не может быть устойчивой. В связи с этим достаточно важен вопрос о том, какое минимальное количество молекул растворителя может обеспечить стабилизацию ионной пары. Авторы работы [8] показали, что на поверхности потенциальной энергии пентамера воды (Н20)5 кроме основного изомера с наименьшей энергией, отвечающего циклической структуре (I), можно найти точки локальных минимумов, которые интересны для сформулированной проблемы (II). Эти две конфигурации изображены на рис. 3. Структура I составлена из 5 молекул воды, а структура II содержит частицы Н30 (3) и 0Н- (1), разделенные слоем трех молекул воды. В работе [8] был проведен анализ найденных стационарных точек, а вычислением гармонических частот колебаний было подтверждено, что они отвечают точкам истинных минимумов.
В данной работе предлагается вариант оценки энергетического профиля на пути реакции при переходе от структуры I к структуре с ионной парой II. Для этого была выбрана подсистема из трех молекул воды (1, 2, 3), между которыми предполагается перенос протона от частицы 1 к частице 3. Расчет энергий и сил в этой подсистеме проводили методом Хартри-Фока с двухэкспонентными базисами, в то время как ассистирующие процессу две молекулы воды (4, 5) описываются как эффективные фрагменты. Точки вдоль координаты реакции задаются следующим образом. На первом этапе проведено растяжение связи 0-Н в молекуле 1 , направленной вдоль водородной связи с молекулой 2. При этом всем геометрическим параметрам системы, включая положения молекул 4 и 5, описываемых
а
б
11 ВМУ, Химия, № 2
жесткими ЭФ, разрешено подстраиваться под минимум полной энергии. При достижении расстояния 1,03 А между мигрирующей водородной частицей и кислородом молекулы 2, характерного для катиона Н30+, происходит переключение на новую степень свободы. А именно, последовательно увеличивается расстояние О-Н вдоль водородной связи частиц 2 и 3, в то время как все остальные геометрические параметры варьируются.
На рис. 4 изображен полученный таким образом профиль пути реакции. Видно, что основные затраты энергии происходят на первом этапе при переносе протона между двумя соседними молекулами воды. На втором этапе практически сразу же достигается точка на вершине барьера (промежуточная геометрическая конфиграция показана на рис. 4), и профиль приходит к точке локального минимума, отвечающего структуре рис. 3, б. Анализ электронных заселенностей показывает, что мигрирующая водородная частица несет заряд меньше, чем +1, т.е. не столь однозначно процесс можно классифицировать как перенос чистого протона. Нельзя исключить, что другие варианты пути реакции для переноса водородной частицы между структурами I и II будут более предпочтительными, нежели предлагаемый.
Результаты данной работы демонстрируют возможности современных квантовохимических подходов для эффективных расчетов энергетических профилей путей реакций при явном учете молекул растворителя в рамках кластерного
Рис. 4. Энергетический профиль пути реакции 2Н2О ^ ^ Н3О++ ОН вдоль координаты, соответствующей движению водорода на первом этапе от молекулы 1 к молекуле 2 и на втором этапе к молекуле 3. Для каждого значения координаты реакции геометрические параметры всей системы (Н2О)5 отвечают минимуму полной энергии
приближения. Построение таких профилей представляет первую часть общей задачи моделирования реакций в растворах. Дальнейшим развитием могут быть молекулярно-динамические расчеты различных процессов, протекающих либо на одной, либо на нескольких адиабатических потенциальных поверхностях всей системы. Скорость вычисления энергии и градиентов энергии при описании молекул растворителя в рамках метода эффективных фрагментов оказывается вполне достаточной для того, чтобы молекулярно-динамическое моделирование стало осуществимым.
При написании данной статьи использованы работы, поддержанные РФФИ (проект № 99-03-33178), а также программой «Университеты России» (проект № 990939). Авторы благодарят А. А. Грановского и А.В. Боченкову за помощь в проведении квантовомеханических расчетов.
Рис. 3. Равновесная конфигурация кластера: а - (Н2О)5, соответствующая одному из локальных минимумов пентамера воды; б - (ОН- + Н3О+)- @ (Н2О)3
СПИСОК ЛИТЕРАТУРЫ
1. Schmidt M.W., Baldridge K.K., Boatz J.A. et al. // J. Comput.
Chem. 1993. 14. P. 1347.
2. Day P.N., Jensen J.H., Gordon M.S. et al. // J.Chem.Phys. 1996.
105. P. 1968.
3. Day P.N., Pachter R., Gordon M.S., Merrill G.N. // J.Chem.Phys.
2000. 112. P. 2063.
4. Granovsky A.A. URL http://classic.chem.msu.su/gran/gamess/
index.html.
5. Johnsson B., Karlstrom G., Wennerstrom H. // J.Am er.Chem.Soc.
1978. 100. P. 1658.
6. Peng Z, Merz K.M. // J. Amer. Chem. Soc. 1978. 114. P. 2733.
7. Pinsent B.R.W., Pearson L., Roughton F.J.W. // Trans. Far. Soc.
1956. 52. P. 1512.
8. Lee C, Sosa C, Novoa J.J. // J. Chem. Phys. 1995. 103. P. 4360.
а
б
Поступила в редакцию 26.12.00
