Учет сольватации в квантовохимических расчетах

Учет сольватации в квантовохимических расчетах

КЛАССИФИКАЦИЯ МОДЕЛЕЙ

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

В основе современных представлений о взаимодействии растворенных соединений с растворителем лежат ставшие классическими работы Борна, Дебая, Леннард-Джонса, Лондона, Кирквуда и Онзагера. Эти работы, их дальнейшее развитие и роль в химии подробно проанализированы в монографии [84]. Однако использовать непосредственно приведенные там результаты для учета сольватации в квантовохимических расчетах не удается, поэтому были разработаны более простые модели, которые позволяют рассчитать энергию сольватации с достаточно высокой степенью точности. Прежде чем перейти к описанию этих моделей, рассмотрим возможные типы сольватации.

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

Несколько особняком лежит случай сольватации малых одно-, двух- и некоторых трех- и четырехатомных ионов. Для этих соединений характерны все перечисленные особенности, но связь оказывается необычайно прочной. Так, энергия гидратации ОН - составляет 445 кДж/моль, Н3О+ 382 кДж/моль, NH3 344 кДж/моль. Такие структуры, по-видимому, было бы правильнее рассматривать как координированные.

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

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

В некоторых работах рассматривались дискретные системы, содержащие 200-300 молекул растворителя (метод Монте-Карло с атом-атомными потенциалами). Такие большие кластеры пригодны для описания как специфического, так и неспецифического взаимодействия растворенного соединения с растворителем. Однако расчеты для них связаны с очень большими затратами машинного времени. Поэтому для учета неспецифической сольватации широкое распространение получили континуальные модели (сольватонная модель, модель самосогласованного поля реакции). В них не учитывается микроскопическая структура растворителя, поэтому отпадает необходимость определять строение сольватационной оболочки. Благодаря этому расчет становится существенно проще.

Полную энергию молекулы в растворителе можно представить как сумму двух отдельных вкладов:

Еполн= ЕМ + ЕС ,

где ЕМ -- энергия изолированной молекулы; ЕС -- энергия сольватации. В предположении отсутствия специфического взаимодействия между растворенным соединением и растворителем, т.е. существенного перераспределения заряда между ними (этот случай будет рассмотрен отдельно), ЕС является суммой трех вкладов:

ЕС = Еэл + Едисп + Екав.

где Еэл - энергия электростатического взаимодействия между собственными и наведенными зарядами растворенного соединения и молекулами растворителя;

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

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

СОЛЬВАТОННАЯ МОДЕЛЬ

Эта модель была предложена Джано [85] и Клопманом [86]. В ее основе лежит предположение, что с каждым атомом растворенного соединения связан индуцированный в растворителе электрический заряд - сольватон. Заряд сольватона по абсолютной величине равен заряду на атоме, который его индуцировал, но имеет противоположный знак. Сольватоны не взаимодействуют между собой. Энергия взаимодействия между сольватоном С и атомом А вычисляется по следующей формуле:

EСА=qCqAе2(е - 1)/(2rе)

где r -- радиус атома А, если сольватон С индуцирован зарядом данного атома, или расстояние между атомами А и В, если сольватон индуцирован зарядом атома В;

qC и qА -- заряды сольватона С и атома А в единицах заряда электрона;

е - заряд электрона; е - диэлектрическая проницаемость среды.

Чаще для оценки ЕСА используют аналогичное выражение с кулоновскими интегралами, взятыми из квантовохимического расчета, т.е. е2/r полагается равным гАА, если сольватон С индуцирован зарядом на атоме А, или гАB, если сольватон С индуцирован зарядом атома В (гАA и гАB -- одно- и двухцентровые кулоновские интегралы соответственно; см.разд. 1.2). Для энергии сольватации в этом приближении можно записать следующее выражение, где сумма берется по всем атомам растворенной молекулы:

В качестве примера использования сольватонной модели для изучения механизмов органических реакций рассмотрим результаты работы [87], в которой изучена Z,E-изомеризация б,в-дизамещенных б,в-непредельных альдегидов.

Расчеты были проведены методами ППДП/2 и МЧПДП/3. Эта реакция связана с внутренним вращением вокруг двойной связи С=С и в газовой фазе идет с очень высоким активационным барьером 200 - 250 кДж/моль. Преодолеть такой активационный барьер при нормальных условиях очень трудно, поэтому скорость газофазной реакции будет ничтожно мала.

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

1. Если растворитель обладает кислотными свойствами, то становится возможной протонизация исходного соединения, после чего по данным квантохимических расчетов активационный барьер реакции должен уменьшиться до 105 - 110 кДж/моль.

Учет влияния растворителя с использованием сольватонной модели (величина е полагалась равной 5 и 10) приводит к дополнительному понижению высоты активационного барьера, и он становится равным 75 - 85 кДж/моль. Реакции с таким активационным барьером могут проходить при нормальных условиях. Таким образом, если растворитель обладает кислотными свойствами, то он может существенно ускорить реакцию Z,E-изомеризации б,в-дизамещенных б,в-непредельных альдегидов.

2. Если растворитель способен катализировать перенос атома Н от группы СН3 к СО, то становится возможным второй механизм увеличения скорости изомеризации:

Высота активационного барьера изомеризации, проходящей по второму механизму, определяется высотой активационного барьера переноса протона от группы СН3 к карбонильному кислороду. Эта реакция, согласно результатам квантовохимических расчетов в газовой фазе, требует затраты около 65 кДж/моль. Учет сольватации повышает эту величину до 80 кДж/моль. Если растворитель обладает хорошей способностью к переносу протона (такой же, как вода), то высота активационного барьера будет близка к тепловому эффекту реакции, т.е. активационный барьер будет немногим выше 80 кДж/моль. С таким активационным барьером реакция может проходить при нормальных условиях.

Сольватонная модель проста и экономична. Однако она обладает рядом серьезных недостатков. Самый неприятный из них заключается в том, что сольватоны локализуются в непосредственной близости от атомов растворенной молекулы, которые их индуцировали, т.е. фактически проникают внутрь молекулы. В случае ионов сольватонная модель учитывает изменение энергии сольватации только за счет изменения степени делокализации заряда, хотя в действительности определяющим фактором в данном случае является их объем. Если изменение степени делокализации заряда коррелирует с изменением объема иона, то сольватонная модель будет давать правильные результаты. Из этого становятся ясными ограничения области ее применения. В частности, сольватонная модель не может учесть эффекты экранирования заряда неполярными группами или молекулами реагентов. По этой причине ее нельзя применять в комбинации с приближением "супермолекулы" при необходимости учесть эффект специфической сольватации.

В заключение настоящего раздела отметим еще одно важное обстоятельство. Выше было указано, что заряды на атомах, вычисленные разными методами, сильно различаются. Из-за этого сольватонная модель, будучи примененной в рамках разных методов, дает разные результаты. Она неплохо работает с методом ППДП/2. Мы пробовали использовать ее также в расчетах методами МЧПДП/2, МЧПДП/3 и МПДП, но результаты получились гораздо хуже.

МОДЕЛЬ САМОСОГЛАСОВАННОГО ПОЛЯ РЕАКЦИИ

В рамках макроскопической теории (шарик с фиксированным распределением зарядов погружен в среду с диэлектрической проницаемостью е, заряды локализованы внутри шарика) энергия взаимодействия с растворителем задается электростатическим классическим потенциалом

,

где U0=(q2/2б)(1- 1/е); U1=(м2/2б3)[2(е - 1)/(2е + 1)]

и т.д. Здесь q, м и б - заряд, дипольный момент и радиус шарика соответственно. Первый член этого ряда называют борновским, второй - онзагеровским. Следующие члены ряда зависят от квадрупольного и более высоких мультипольных моментов шарика. Для электрически нейтральных шариков Uo=0 и первым ненулевым членом ряда будет U1. В работах Тапиа с сотр. [88 - 92] предложено добавить аналог члена U1 к гамильтониану изолированной молекулы и учесть таким способом взаимодействие с растворителем. В результате ими было получено следующее выражение для вычисления гамильтониана молекулы, находящейся в растворителе:

H=H0 -

где g - тензор, который определяется восприимчивостью электростатического поля молекулы средой. В работах [88 - 92] он задавался в параметрической форме.

Этот метод учета взаимодействия с растворителем получил название модели самосогласованного поля реакции. Основные его недостатки очевидны. Они заключаются в необходимости учета зависимости тензора g от значения диэлектрической проницаемости растворителя и размеров молекулы. Это нельзя сделать без очень грубых дополнительных предложений. Обычно полагают, что диагональные элементы g равны 2(е - 1)/[(2е + 1)б3], а недиагональные - нулю. Однако перенесение макроскопической формулы на микроскопический уровень, вообще говоря, неверно, диэлектрическая проницаемость при этом теряет свой физический смысл, а выбор значения параметра б (радиуса молекулы) достаточно произволен. Поэтому модель самосогласованного поля реакции имеет ограниченную область применения. Практически ею можно пользоваться при рассмотрении реакций изомеризации и переноса протона. В остальных случаях ее использование нежелательно.

ПРИБЛИЖЕНИЕ СУПЕРМОЛЕКУЛЫ

С точки зрения построения модели наиболее простым способом учета сольватации является включение большого числа молекул среды в систему, для которой проводится квантовохимический расчет. Все ее электроны (сольватируемой молекулы и молекул среды) включаются в электронный гамильтониан. Этот способ является непосредственным обобщением квантовохимических методов, развитых для отдельных молекул на случай больших систем, состоящих из нескольких или даже большого числа отдельных молекул. Если нам удастся учесть таким способом взаимодействие растворенной молекулы с большим числом молекул cреды, рассчитать энергетически наиболее выгодную конформацию растворенной молекулы и конфигурацию молекул растворителя и получить для этой конфигурации электронную волновую функцию, то мы сможем объяснить или даже предсказать практически все интересующие нас свойства молекулы в растворе. При этом следует с особой осторожностью подходить к выбору квантовохимического метода, к которому в этом случае предъявляются повышенные требования, а именно он должен быть пригоден для изучения межмолекулярных взаимодействий.

Из предыдущей главы было видно, что многие квантовохимические методы, которые успешно используются для изучения реакционной способности органических соединений дают неправильные результаты при расчете параметров, характеризующих межмолекулярные взаимодействия. Например, метод МПДП дает неудовлетворительные результаты при расчете систем с водородными связями. Поэтому при использовании приближения "супермолекулы" для учета сольватации приходится выбирать достаточно совершенный квантовохимический метод. При этом возникают трудности, связанные с очень большим размером супермолекулы, которая состоит из растворенного соединения и молекул среды. Кроме того, в супермолекуле обычно неизвестно большое число геометрических параметров, которые определяют строение сольватационной оболочки. Поэтому выполнить реально такие расчеты удается только для систем с небольшим числом молекул среды.

Такой учет сольватации, при котором система из сольватируемой молекулы и некоторого ограниченного числа молекул растворителя рассчитывается квантовохимическим методом как одна молекула, получил название "приближения супермолекулы" [93, 94]. Расчеты в этом приближении широко распространены. К ним относятся работы, в которых вычисляются параметры для комплексов, состоящих из растворенной молекулы и одной молекулы растворителя, а также расчеты для комплексов с двумя, тремя и т.д. молекулами среды. В этих работах впервые удалось получить данные о строении и энергии взаимодействия для гидратных оболочек ионов ОН- и Н3О+ [65]. Несколько позже аналогичными методами исследованы гидратационные оболочки простейших ионов: Li+, Be2+, Na+, Mg2+, Al3+, К+, Ca2+, F-, Сl-, NH4+ [95 - 102], CH5+, CH5- [103], алкиламмониевых ионов [104].

Методика проведения перечисленных выше работ заложила основу подхода к изучению сольватации в приближении супермолекулы. Суть этих подходов заключается в следующем. Вначале рассчитывается энергия взаимодействия иона с одной молекулой воды и определяется наиболее энергетически выгодная конформация, потом в систему добавляется еще одна молекула воды и вновь рассчитываются энергия взаимодействия и структура комплекса и т.д. В результате такого расчета получается набор величин, которые являются энергиями гидратации иона А каждой последующей молекулой воды (En,n-1):

А*(n - 1)Н2O + Н2O > А*nН2O + En,n-1.

Для простейших ионов эти энергии были сопоставлены с экспериментальными данными, полученными методами масс-спектрометрии высоких давлений и циклотронного резонанса. Хорошее согласие с экспериментом подтвердило широкие возможности квантовохимических расчетов для изучения сольватации. Кроме того, были получены данные о числе молекул воды в первой гидратационной сфере и их конфигурация.

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

1) детально изучается взаимодействие сольватируемой молекулы с одной молекулой растворителя и получается предварительное представление о строении сольватационной оболочки;

2) изучается влияние второй, третьей и т.д. молекул растворителя на строение сольватационной оболочки и уточняются данные, полученные на первом этапе;

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

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

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

На втором этапе увеличивают число молекул растворителя в супермолекуле. При этом используют данные о положениях локальных минимумов, полученные на первом этапе, и лишь уточняют их положение. На третьем этапе в супермолекулу включают максимально возможное количество молекул растворителя (исходя из возможностей ЭВМ) и проводят расчет того или иного параметра.

В качестве примера рассмотрим результаты работы [105], в которой изучено строение гидратационной оболочки мочевины. Это соединение обладает рядом уникальных свойств. Так, в его присутствии резко снижается температура, при которой в водных растворах происходит денатурация белков и других биополимеров, повышается растворимость неэлектролитов и т.д. Принято говорить, что в присутствии мочевины разрушаются водородные связи между молекулами воды, другими словами, молекулы воды в водном растворе становятся более подвижными. Это приводит к перечисленным выше химическим и биохимическим эффектам. Однако механизм этого явления на микроскопическом уровне установлен не был. Для решения этого вопроса было изучено [105] строение гидратационной оболочки мочевины с помощью квантовохимических расчетов методами ППДП/2 и ОСТ-ЗГФ.

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

В этом комплексе молекула воды связана двумя водородными связями с двумя атомами Н связей N--Н, которые находятся в транс-положении по отношению к карбонильной группе. Расчет также показал, что в этом положении молекула воды обладает необычайно высокой подвижностью. Она может смещаться на расстояния до 0,15--0,20 нм от положения равновесия, и при этом энергия взаимодействия в комплексе почти не меняется. Таким образом, было обнаружено, что в исследованном комплексе молекулы воды и мочевины частично сохраняют свободу независимых перемещений в пространстве.

На втором этапе было увеличено количество молекул воды в комплексе, однако оказалось, что вторая молекула Н2О не может связаться с молекулой мочевины аналогично первой и образовать комплекс типа

Таким образом, первая молекула воды, присоединяясь к мочевине, попадает в широкую потенциальную яму (именно этим объясняется ее высокая подвижность), но вторая молекула воды в эту потенциальную яму попасть не может, поэтому при увеличении количества молекул воды в комплексе высокая подвижность одной из них сохраняется. По-видимому, наличие этой достаточно широкой и глубокой потенциальной ямы обусловливает способность мочевины "разрушать структуру воды".

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

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

У метилмочевины наиболее устойчивой является следующая конформация:

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

У диметилмочевины возможны два изомера

У первого из них две связи N--Н в трансположении по отношению к карбонильной группе, а у второго - лишь одна. Следовательно, первый изомер будет повышать термодинамическую температуру воды и приводить к тем же эффектам, что и мочевина, а второй - нет. Это теоретическое предсказание получило экспериментальное подтверждение.

Приближение супермолекулы позволяет решать многие задачи, связанные с влиянием сольватации на реакционную способность органических соединений. Однако при его использовании необходимо самое серьезное внимание уделить выбору метода расчета. При этом следует руководствоваться двумя требованиями: 1) метод расчета должен быть достаточно точным и хорошо передавать как основные свойства растворенной молекулы, так и строение сольватационных оболочек; 2) расчет в приближении супермолекулы связан с вычислением электронной волновой функции для очень большой системы, и он должен быть практически реализуем исходя из возможностей имеющихся ЭВМ, их быстродействия и памяти. До недавнего времени расчеты с учетом сольватации в приближении супермолекулы проводились методами ППДП/2 и ОСТ-ЗГФ. Именно эти методы были использованы в приведенном выше примере изучения особенностей строения гидратационной оболочки мочевины. В настоящее время предпочтительнее использовать метод МПДП/Н и неэмпирические расчеты в базисе 3-21ГФ.

МЕТОД АТОМ-АТОМНЫХ ПОТЕНЦИАЛОВ

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

Различные эмпирические потенциалы, которые предлагались различными авторами для описания межмолекулярных взаимодействий, не удовлетворяют точности, необходимой при учете сольватации. Отсутствие достаточно надежных потенциалов делало невозможным использование такого подхода для изучения сольватации. Существенный прогресс в этой области был достигнут благодаря работам Клементи с сотрудниками. Ими была выдвинута идея использовать неэмпирические квантовохимические расчеты для определения потенциалов межмолекулярных взаимодействий [106--109]. Первоначально потенциалы генерировались в численной форме путем неэмпирического расчета энергии взаимодействия сольватированной молекулы и молекулы растворителя. Но, так как дальнейшее использование численного потенциала для построения сольватационной оболочки затруднительно, была подобрана аналитическая функция для его аппроксимации, которая представляла собой сумму атом-атомных потенциалов. При этом все атомы в молекуле были разбиты на классы в зависимости от того, к каким функциональным группам и в каких положениях в группах находится данный атом. В результате число классов во много раз превышало число различных атомов. Для атомов каждого класса подбирались свои атом-атомные потенциалы. Аналитическая форма, в которой производился поиск атом-атомных потенциалов, выбиралась разной и зависела от базиса, использованного в расчете. При вычислении потенциала взаимодействия между молекулами в небольших базисах обычно использовалась относительно простая аналитическая функция

Uij = -Aij/r6ij + Bij/r12ij + Cijqiqj/rij

где Uij - энергия взаимодействия между атомами i и j; rij - расстояние между этими атомами; qi и qj - заряды на атомах; Аij, Вij и Сij - эмпирические параметры, зависящие от того, к каким классам принадлежат атомы i и j.

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

Число классов атомов у молекул типа оснований ДНК и аминокислот составляет несколько десятков, число неизвестных параметров в атом-атомных потенциалах достигает нескольких сотен. При определении значений этих параметров приходится варьировать относительное положение и взаимную ориентацию молекул в достаточно широких пределах, практически для вычисления каждого параметра приходится делать 15 - 20 расчетов. Таким образом, для расчета потенциала межмолекулярного взаимодействия молекул среднего размера типа оснований ДНК и молекул воды необходимо сделать несколько десятков тысяч расчетов полной энергии системы неэмпирическим методом. Поэтому процедура подбора параметров связана с очень большими затратами машинного времени. Но следует отметить одно благоприятное обстоятельство: по мере накопления наборов параметров и создания их банка для каждого нового соединения объем вычислений сокращается, так как оказывается возможным отнести большинство атомов к уже известным классам, для которых все параметры аналитических потенциалов известны из расчетов других молекул [110].

После нахождения потенциалов расчет строения сольватационной оболочки и энергии взаимодействия между растворителем и растворенным соединением становится относительно простой задачей, аналогичной задачам конформационного анализа. Аналитическая форма, в которой в настоящее время найдены потенциалы для описания взаимодействия молекул среднего размера и молекул воды, также совпадает с наиболее широко распространенными потенциалами, которые используются в конформационных расчетах. Однако параметры в потенциалах Клементи для межмолекулярных взаимодействий имеют совершенно иную природу. В конформационном анализе потенциалы типа 6 - 12 описывают ван-дер-ваальсовы взаимодействия между атомами, а у Клементи - электронодонорные и электроноакцепторные взаимодействия. Третий член в формуле для атом-атомных потенциалов соответствует кулоновскому взаимодействию. Для электрически нейтральных молекул значение коэффициента Сij в потенциалах Клементи близко к единице. Однако для ионов оно не превышает 0,5; это связано, по-видимому, с эффектами экранирования и перераспределения заряда.

Использование потенциалов Клементи позволяет рассматривать гидратацию весьма сложных молекул большим числом молекул воды. Пока число молекул воды не превышает 10 - 15, каких-либо существенных трудностей при расчете строения гидратационной оболочки не возникает. Однако при дальнейшем увеличении числа молекул воды появляется ряд новых проблем. Для достаточно точного описания гидратационной оболочки даже небольшого соединения количество молекул воды желательно увеличить до 200 - 300. При расчете строения такой огромной гидратационной оболочки основная трудность заключается в существовании большого числа структур с близкими энергиями. Задача сводится к нахождению всех таких структур, определению вероятности реализации каждой из них и усреднению по всем найденным структурам. В таком расчете приходится учитывать температурную зависимость.

В работах Клементи показано, что для нахождения строения гидратационных оболочек можно успешно использовать метод Монте-Карло, с помощью которого были проведены расчеты строения гидратационных оболочек ряда простых ионов с учетом их взаимодействия с 200 - 250 молекулами воды. При этом возникла еще одна проблема. Распределение молекул воды в расчетах методом Монте-Карло носит вероятностный характер, поэтому перед исследователями встала задача перехода к таким простым и наглядным характеристикам гидратации, как число молекул воды в первой гидратационной сфере и ее радиус. Для получения этой информации было предложено вычислить зависимость плотности атомов водорода или кислорода от расстояния до центра иона. На таких графиках получается ряд четко выраженных максимумов. Их положение для атомов кислорода обычно связывают с радиусами гидратационных оболочек, а площадь под кривыми - с количеством молекул воды в оболочке. Ниже приведены радиусы первых гидратационных оболочек (R) и число молекул воды в них (N), вычисленные таким способом [108, 109].

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

1. Используется приближенный квантовохимический метод (для молекул среднего размера весьма грубый) для вычисления параметров атом-атомных потенциалов. В случае небольшого числа молекул растворителя ошибки могут быть невелики, но по мере увеличения их количества они будут накапливаться.

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

3. Атом-атомные потенциалы, которые обычно используют для изучения сольватации, не являются аддитивными функциями, а взаимодействия трех тел учесть довольно сложно и этого почти никогда не делают, хотя эти коллективные взаимодействия существенно влияют на результаты расчета (это, по-видимому, самый большой недостаток метода Клементи).

В ряде работ подход Клементи к учету сольватации был использован для изучения влияния растворителя на поверхности потенциальной энергии органических реакций [111, 112]. Эти работы будут более подробно рассмотрены ниже. Здесь же мы только отметим, что проведение таких расчетов требует очень больших затрат машинного времени. Их порядок таков: 1) вычисляют полную энергию для какой-либо точки на поверхности потенциальной энергии газофазной реакции; 2) в этой точке рассчитывают параметры атом-атомных потенциалов, описывающих взаимодействие реагентов с молекулой растворителя; 3) с помощью атом-атомных потенциалов, полученных на предыдущем этапе расчета (см. пункт 2), методом Монте-Карло вычисляют энергию сольватации.

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

МОДЕЛЬ ТОЧЕЧНЫХ ДИПОЛЕЙ

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

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

Последовательная электростатическая модель, в которой учитывалось отталкивание в потенциале, описывающем межмолекулярные взаимодействия, была предложена в работе Уоршела [113]. Молекулы растворителя в ней моделировались шариками с фиксированным дипольным моментом и ван-дер-ваальсовым радиусом. Если молекул растворителя немного, то использование метода Уоршела не вызывает затруднений. Однако при увеличении количества молекул растворителя возникают трудности, связанные с поиском оптимальной структуры сольватной оболочки. Использование для этой цели метода Монте-Карло связано с очень большим объемом вычислений и не годится для решения прикладных задач, а методы оптимизации геометрии типа наискорейшего спуска позволяют найти лишь один из многочисленных локальных минимумов (не обязательно самый глубокий). Поэтому метод Уоршела в настоящее время применяется редко. В качестве примера его использования приведем работу [114], в которой с его помощью было учтено влияние гидратации на потенциальный профиль реакции OH- + CO2 > Н3СО-. Было показано, что гидратация существенно уменьшает тепловой эффект реакций. Кроме того, было установлено, что в газовой фазе эта реакция идет без активационного барьера, а в растворе - с активационным барьером.

Другой вариант модели точечных диполей был предложен в работе [115]. Чтобы избежать вычисления оптимальной структуры сольватной оболочки, дипольные моменты молекул растворителя были раздроблены и почти непрерывно и равномерно "размазаны" по объему растворителя, т.е. каждая молекула растворителя была заменена большим количеством точечных диполей с малыми дипольными моментами. Такое ''размазывание" дипольного момента молекул растворителя ранее было использовано Уоршелом при построении модели ланжевеновских диполей [113].

Точечные диполи с малыми дипольными моментами могут быть расположены в узлах любой упорядоченной решетки. Из них следует выделить ту часть точечных диполей, которые вносят существенный вклад в энергию сольватации. Положение в пространстве этой части точечных диполей должно удовлетворять следующим двум условиям. Во-первых, точечный диполь не может находиться ближе некоторого критического расстояния ни к одному из атомов растворенной молекулы. Во-вторых, энергия взаимодействия между точечными диполями и растворенной молекулой должна превышать некоторую пороговую величину. Физический смысл первого условия очевиден (точечные диполи не должны попасть внутрь растворенного соединения). Суть второго условия заключается в учете взаимодействия растворенного соединения лишь с теми точечными диполями, которые находятся в области достаточно сильного электрического поля растворенного соединения. Их дипольные моменты будут ориентированы в направлении этого поля, и поэтому будут вносить основной вклад в энергию сольватации. Дипольные моменты остальных точечных диполей из-за взаимодействия молекул растворителя между собой будут ориентированы в произвольных направлениях, поэтому энергия их взаимодействия с электрическим полем растворенного соединения будет мала. Этим вкладом в работе [115] предложено пренебречь.

Таблица 1 Энергии гидратации простых ионов и супермолекул, вычисленные с помощью модели точечных диполей (кДж/моль) [115]

Результаты расчета энергий гидратации простых ионов и некоторых супермолекул (комплексов, образованных ионом и одной или двумя молекулами воды) методом МПДП с использованием модели точечных диполей [115] приведены в табл. 2.1, из которой видно, что получается хорошее согласие с экспериментом. Однако не следует переоценивать возможности модели точечных диполей в описанном выше варианте. Она годится лишь для достаточно грубого качественного учета сольватационных эффектов. Основное ее достоинство - малые затраты машинного времени, поэтому ее легко применять в прикладных расчетах для получения качественно правильных результатов для реакций в полярных растворителях.

ИЗУЧЕНИЕ МЕХАНИЗМОВ ОРГАНИЧЕСКИХ РЕАКЦИЙ В РАСТВОРАХ МЕТОДАМИ КВАНТОВОЙ ХИМИИ

Диссоциация молекул в полярных растворителях.

Рассмотрим вопрос о механизме диссоциации молекул в полярных растворителях. Химики-экспериментаторы очень часто имеют дело с такими процессами при изучении многих органических реакций, поэтому необходимо уметь их моделировать в квантовохимических расчетах. На первый взгляд эта задача кажется очень простой, однако при более близком знакомстве с ней выясняется, что механизм этих процессов весьма сложен и общепринятые представления о нем являются неточными. В квантовохимических работах образование ионов в полярных растворителях наиболее часто рассматривалось на примере следующей модельной реакции [116]:

FH + NH3 > F- + NH4+. (I)

FH - кислота, a NH3 - основание, поэтому в водных растворах эта система должна существовать в ионной форме F- + NH4+. Однако в квантохимических расчетах удавалось получить лишь локальный минимум для ионной пары и то только при фиксированном и достаточно большом расстоянии между атомами фтора и азота. Включение этого расстояния в число варьируемых параметров неизменно приводило к исчезновению локального минимума для ионизированной системы. В качестве примера на рис. 2.1 приведены результаты такого расчета методом МПДП/Н, сольватация моделировалась в приближении супермолекулы шестью молекулами воды.

Рис. 2.2. Сечение ППЭ для реакций переноса протона 1 - реакция (I); 2 - Н2О...Н2О > ОН- + Н3О+. R - расстояние между атомами F--Н или О--Н. Расчет выполнен методом МПДП/Н; для учета сольватации использована модель точечных диполей

Е, кДж/моль

Рис. 2.3. Сечение ППЭ для реакции переноса протона (I) Расчет выполнен методом МПДП/Н, для учета сольватации использована модель точечных диполей м(1030*D/м3): 1 - 0; 2 - 0,06; 3-0,12; 4 - 0,15.

При дальнейшем исследовании этой системы [116] было установлено, что минимум для электрически нейтральной системы является локальным, а глобальный минимум соответствует разделенной ионной паре, в которой расстояние между ионами составляет около 1 нм (рис. 2.2). Для контактных ионных пар минимума обнаружено не было, что свидетельствует о невозможности их существования в водных растворах.

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

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

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

В модель точечных диполей дискретные свойства сред такого типа можно ввести, например, путем ограничения количества точечных диполей, которые могут взаимодействовать с растворенной молекулой. На рис. 2.3 показаны результаты расчета сечения поверхности потенциальной энергии для реакции переноса протона в предположении, что количество точечных диполей объемом Vт.д= 10-3 нм3 не может превышать 200. Видно, что при таком ограничении количества точечных диполей глобальный минимум соответствует контактным ионным парам. Получить его без введения в использованную сольватационную модель элемента дискретности невозможно.

Реакции нуклеофильного замещения, протекающие по механизму SN2

За последние, годы существенно расширились наши знания о детальных механизмах реакций нуклеофильного замещения, широко распространенных в химии [117, 118]. Методами масс-спектрометрии высоких давлений [119], ионного циклотронного резонанса [120] и квантовой химии [121 - 123] показано, что в газовой фазе первоначально образуется устойчивый предреакционный комплекс (K1). Ему соответствует первый минимум на профиле потенциальной энергии. При дальнейшем движении вдоль координаты реакционная система преодолевает активационный барьер, после чего образуется второй комплекс (К2), который распадается на конечные продукты:

Таблица 2.2 Энтальпия (ДН) и активационный барьер (Е#) реакции (III) в водном растворе (кДж/моль) [125]

Экспериментально эти реакции обычно изучаются в растворах. Влияние сольватации на профиль потенциальной энергии было установлено на основе результатов квантовохимических расчетов.

В работе [124] рассчитано сечение поверхности потенциальной энергии для реакции Сl- + Н3ССl > СlСН3 + Сl- неэмпирическим методом 3-21ГФ. Для моделирования сольватации использовано приближение супермолекулы (к реагентам добавлено две молекулы воды). Было показано, что сольватация существенно повышает активационный барьер реакции. Однако вычисленное значение энергии активации оказалось сильно заниженным по сравнению с экспериментом. Это связано, во-первых, с недостаточным количеством молекул воды, включенным в гидратационную оболочку реагентов, и, во-вторых, с использованием базиса 3-21 ГФ, расчеты в котором существенно переоценивают стабильность комплексов K1 и К2 и переходного состояния (ПС) реакции (II) в газовой фазе.

В работе [125] модифицированным методом ППДП/2 вычислен профиль потенциальной энергии для реакций [126],

X- + Н3СY > XСН3 + Y- (X,Y=H,F,Cl) (III)

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

В работе [127] профиль потенциальной энергии для реакции

Сl- + Н3ССl > СlСН3 + Сl- (IV)

ОН- + Н3СВr > НОСН3 + Вr- (V)

рассчитан методом МПДП. Для учета гидратации использована модель точечных диполей. Полученные результаты приведены на рисунке 2.4 и в таблице 2.3. Видно, что гидратация существенно изменяет профиль потенциальной энергии исследованных реакций, однако траектория движения реагентов, которая проходит по дну долины на поверхности потенциальной энергии, при этом остается практически неизменной. Аналогичные результаты были получены в работах [128, 129], выполненных методом МПДП/Н с использованием приближения супермолекулы для расчета гидратации. Приблизительно такая же форма профиля потенциальной энергии для этих реакций получена и в работе [111], в которой расчеты проводились неэмпирическим методом в большом гауссовом базисе и методом Монте- Карло с атом-атомными потенциалами. Влияние сольватации на траекторию движения реагентов в последней работе не рассматривалось.

Рис. 2.4. Сечение ППЭ для реакций (IV) (a) и (V) (б)

1 - в газовой фазе; 2 - в водном растворе. Д = R1 - R 2; R1 - расстояние между атомом углерода и уходящим атомом хлора (брома); R2 - расстояние между атомом углерода и атакующим атомом хлора (брома).

Таблица 2.3 Характеристики траектории движения реагентов в газовой фазе и в водном растворе

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

Нуклеофильные реакции карбонильных соединений.

Реакции нуклеофильного присоединения и замещения карбонильных соединений играют важную роль в химии и биологии. В частности, к ним относятся реакции образования пептидных связей в условиях биосинтеза белка. Обычно эти реакции идут в водной среде. При этом растворитель играет роль резервуара, с которым происходит обмен протонами. Поэтому при подходе нуклеофильного реагента NuH к карбонильному соединению XCOY возможно образование интермедиатов с разным характером протонирования: