Молекулярное моделирование с программой PC GAMESS: от двухатомных молекул до ферментов Текст научной статьи по специальности «Нанотехнологии»

УДК 539.19

МОЛЕКУЛЯРНОЕ МОДЕЛИРОВАНИЕ С ПРОГРАММОЙ PC GAMESS: ОТ ДВУХАТОМНЫХ МОЛЕКУЛ ДО ФЕРМЕНТОВ

А.В. Немухин, Б.Л. Григоренко, А.А. Грановский

(кафедра физической химии)

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

Введение

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

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

альных пакетов программ. В иерархии квантовохи-мических программ комплекс GAMESS (The General Atomic and Molecular Electronic Structure System) занимает достаточно почетное место. Не последнюю роль в ее популярности играет то обстоятельство, что так называемый американский вариант пакета GAMESS (US) [1] относится к числу свободно распространяемых (некоммерческих) программ, предоставляя пользователю почти те же самые возможности, что и многие другие программные комплексы квантовой химии. В настоящее время функционирование пакета GAMESS (US) поддерживается усилиями группы профессора М. Гор -дона (M. Gordon, Ames Laboratory / Iowa State University, USA). В сети Интернет существует официальная страница GAMESS(US) [2], на которой представлена наиболее важная информация для пользователей, а также сформулированы правила получения различных вариантов программы. Английский вариант программы, GAMESS (UK), поддерживается группой проф. М.Ф. Геста (M.F. Guest, Daresbury Laboratory, UK) [3].

Не наша прерогатива излагать здесь историю создания и развития GAMESS - для этого нужно посмотреть официальную документацию на соответствующих страницах в Интернет [2, 3]. В России программа GAMESS (US) впервые появилась на химическом факультете Московского государственного университета им. М.В. Ломоносова в 1991 г. при любезном содействии профессора Ф. Вейнхольда (F. Weinhold, University Wisconsin-Madison, USA) и с разрешения администратора GAMESS (US) доктора М. Шмидта (M. Schmidt). После решения ряда технических проблем на химическом факультете МГУ появилась возмож-

ность выполнять квантово-химические расчеты с помощью этой программы. Практически сразу же в работу сначала по адаптации, а затем по усовершенствованию программы включился А.А. Грановский. В результате его многолетних усилий был создан вариант программного комплекса, получивший широкую известность как PC GAMESS [4].

Первая стадия работы над созданием программы PC GAMESS заключалась в адаптации оригинальных текстов GAMESS (US) к платформе Intel 386/ 32-bit extended DOS, затем к Intel 386/OS/2 и Intel 386/Windows NT. Поддержка операционной системы Linux была добавлена существенно позже. На следующем этапе практически все участки программы, определяющие скорость ее работы, были переписаны и оптимизированы для процессоров Intel 386 (сначала для Intel Pentium, а на более поздних стадиях и для Intel Pentium Pro/II/III). Были разработаны эффективные алгоритмы работы с оперативной памятью, алгоритмы упаковки и распаковки данных, хранящихся на жестком диске, в режиме реального времени; была добавлена прямая поддержка больших (>2 ГБ) файлов.

После этого началась стадия создания новых и усовершенствования существующих алгоритмов для таких расчетных методов квантовой химии, как HF, MCSCF, TDHF, MP2, CI и MCQDPT. Была добавлена возможность расчетов на уровне MP3 и MP4-SDQ, а затем и на уровне MP4-SDTQ. Одновременно проводилась работа по эффективной поддержке многопроцессорных систем (NT, OS/2, Linux) и вычислительных кластеров (Win32, Linux), были разработаны многопотоковые и параллельные алгоритмы для методов HF, TDHF, MP2, MP3, MP4. Была улучшена масштабируемость многих существующих вычислительных методов, а также добавлена поддержка современных процессоров Intel (Pentium 4) и AMD (Athlon).

В последующих разделах будут приведены некоторые результаты наиболее интересных исследований, выполненных с помощью PC GAMESS, по мере усложнения рассматриваемых молекулярных систем. В обзоре приведены в основом конкретные данные, полученные при участии авторов настоящей работы, но дан также краткий анализ результатов использования PC GAMESS другими исследователями в России.

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

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

Исследования с использованием программы PC GAMESS

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

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

Начнем с примеров построения картины потенциальных кривых двухатомных молекул в основном и возбужденных электронных состояниях с помощью программы PC GAMESS, выполненных для последующих расчетов более сложных систем методом двухатомных фрагментов в молекулах (DIM). Этот метод [5] можно рассматривать как полуэмпирический вариант приближения валентных схем, позволяющий комбинировать парные двухатомные потенциалы в строго определенных электронных состояниях в потенциальные поверхности многоатомных молекул. Особую ценность потенциальным поверхностям метода DIM обеспечивает сравнительно невысокая стоимость вычислений энергии и градиентов энергии для многоатомных систем.

В работе [6] были вычислены потенциальные кривые молекулы OH, использованные далее для моделирования методом DIM поверхности энергии молекулы H2O [7]. На рис. 1 изображены рассчитанные зависимости энергии OH от расстояния R между ядрами для дублетных состояний симметрии П и Х+. Рассмотренные потенциальные кривые коррелируют с указанными на рис. 1 справа пределами

Т а б л и ц а 1

Расшифровка общепринятых сокращений в квантовой химии, используемых в данном обзоре

Англоязычный термин Русскоязычный эквивалент Расшифровка

АО АО атомная орбиталь (часто эквивалент базисной функции)

МО МО молекулярная орбиталь

SCF ССП самосогласованное поле

HF ХФ Хартри-Фока (метод)

RHF ОХФ ограниченный (метод) Хартри-Фока

UHF НХФ неограниченный (метод) Хартри-Фока

GVB ОВС обобщенный (метод) валентных схем

MCSCF МК ССП многоконфигурационный (метод) самосогласованного поля

CASSCF ССП ППАО (метод) самосогласованного поля в полном пространстве активных орбиталей

MP2 МП2 (метод теории возмущений) Меллера-Плессе 2-го порядка

MP4-SDQ МП4-ОДК (метод теории возмущений) Меллера-Плессе 4-го порядка с учетом однократных, двукратных и четырехкратных возбуждений

MP4-SDTQ МП4-ОДТК (метод теории возмущений) Меллера-Плессе 4-го порядка с учетом однократных, двукратных, трехкратных и четырехкратных возбуждений

CI КВ конфигурационное взаимодействие

SOCI КВВП конфигурационное взаимодействие 2-го порядка

CISD+Q КВОД+К конфигурационное взаимодействие с включением в волновую функцию всех одно- и двукратно возбужденных конфигураций и с поправкой Давидсона на квартичные возбуждения

MCQDPT МК КВТВ многоконфигурационная квазивырожденная теория возмущений

SBK СБК (псевдопотенциалы) Стивенса-Баша-Краусса

NBO НСО натуральные связевые орбитали

CAS/NBO ПАП/НСО полное активное пространство с натуральными связевыми орбиталями

DIM ДФМ двухатомных фрагментов в молекулах (метод)

Рис. 1. Потенциальные кривые ОН для состояний 2П и 2£+, рассчитанные в работе [6]. Начало отсчета энергии соответствует 75,0 ат. ед.

диссоциации электронных термов на атомы при больших значениях Я, среди которых есть и ион-парный предел О"Н+. Для получения этой картины были использованы необходимые приемы вычислительной квантовой химии, реализованные в комплексе РС ОЛМБ88.

Для расчетов потенциальных кривых в столь широкой области значений Я требуется многоконфигурационное представление волновых функций. Точки на кривых получены методом конфигурационного взаимодействия с молекулярными орбиталя-ми, оптимизированными в рамках многоконфигурационного приближения 8СБ. В программе (ОЛМБ88(Ш)) и соответственно РС ОЛМБ88 реализован один из наиболее эффективных вариантов МС8СБ, который называется приближением самосогласованного поля в полном пространстве активных орбиталей (СЛ88СБ) [8, 9]. В данном приложении в волновых функциях СЛ88СБ 7 валентных электронов молекулы ОН распределялись всеми возможными способами, совместимыми с заданной

пространственной и спиновой симметрией электронного состояния, по активным орбиталям 2а, 3а, 4а, 1п. Состояния 2£+ попадают в блок симметрии Л1, а состояния П в блок В1 (или В2) наивысшей возможной группы симметрии С2у, в рамках которой можно рассматривать в РС ОЛМБ88 гетеро-ядерные двухатомные молекулы. Поскольку целью расчета было построение потенциальных кривых двухатомной системы для нескольких электронных состояний предпочтительно с одинаковой точностью, в процессе оптимизации орбиталей и коэффициентов перед конфигурационными функциями СЛ88СБ проводили процедуру усреднения по состояниям [10]. С равными весами усредняли матрицы плотности для трех первых состояний (Х2П, 2 П, 3 П) в блоке В1 и для первых четырех состояний (Л2Х+, В2£+, С2£+, 2 А) в блоке 2Л1. Полученные наборы молекулярных орбиталей были затем использованы в приближении конфигурационного взаимодействия второго порядка с выбором в качестве исходных всех конфигураций, задействованых в предварительных расчетах СЛ88СБ (опция 8ОС1 в программе РС ОЛМБ88). Всего было учтено 117408 конфигураций в блоке В1и 118003 в блоке 2Л1. Все вычисления выполнены с использованием базиса атомных орбиталей, который считается в настоящее время вполне приемлемым для неэмпирических расчетов с учетом эффектов электронной корреляции (а^-сс-рУТ2).

В работе [6] тщательно проанализированы погрешности вычисленных молекулярных постоянных для рассмотренной совокупности электронных состояний молекулы ОН. Сделан вывод, что равновесные межъядерные расстояния для связанных состояний рассчитывают с погрешностями, не превышающими 0,01 А, в то время как значения энергии диссоциации могут иметь погрешности до 0,15 эВ. Значительная часть погрешностей при описании потенциальных кривых относится к диссоциационным пределам. После введения поправок с использованием экспериментальных значений для атомных термов вычисленные адиабатические энергии возбуждения из основного состояния X П достаточно хорошо согласуются с оценками по спектральным

данным: 4,06 эВ (эксперимент 4,05); 8,56 (экспери-

2 +

мент

8,65); 10,82 эВ (эксперимент 11,09) для Л £ , В2£+, С2£+ соответственно.

На рис. 2 воспроизведена также рассчитанная в работе [6] картина потенциальных кривых аниона ОН , коррелирующих с низшим диссоциационным

Рис. 2. Потенциальные кривые состояний 1П и 1Х+ аниона ОН , рассчитанные в работе [6]. Показано положение потенциала основного состояния Х2П нейтральной молекулы ОН (начало отсчета энергии соответствует 75,0 ат. ед.)

пределом О(2Р)+Н(ь). Показано также положение потенциала основного состояния Х2П нейтральной молекулы ОН. Этот пример иллюстрирует определенную проблему квантовой теории, связанную с расчетами автоионизационных состояний. В данном случае потенциальная кривая состояния П аниона пересекает кривую состояния той же пространственной симметрии Х2П нейтральной молекулы, и при малых межъядерных расстояниях система ОН-нестабильна относительно распада ОН + е. Для того чтобы построить соответствующую потенциальную кривую, были выполнены такие же расчеты, как и для нейтральной молекулы, но в базис была добавлена еще одна диффузная функция на водородном центре с показателем экспроненты 0,0001.

Как уже отмечалось выше, рассчитанные с помощью программы РС ОЛМБ88 потенциальные кривые ОН и ОН- были успешно использованы (наряду с хорошо известными энергиями Н2 и Н2+) для построения потенциальной поверхности молекулы Н2О. В работе [6] было показано, что профили по-

тенциальной поверхности молекулы воды, вычисленные методом DIM, согласуются с реперными данными высокой спектроскопической точности [11] с различиями в пределах 50 см . Трудно разделить ответственность приближения DIM и надежность двухатомных потенциальных кривых при анализе окончательных результатов для трехатомной системы, но очевидно, что точность неэмпирических расчетов двухатомных потенциалов достаточно высока.

В работе [12] была поставлена более амбициозная задача построить методом DIM потенциальные поверхности малых кластеров воды (H20)n (n < 6). Для реализации этого проекта потребовалось множество потенциальных кривых двухатомных фрагментов, и использование неэмпирических квантово-химических расчетов представляло практически единственную возможность решить проблему. Особые трудности были связаны с необходимыми для схемы DIM потенциалами ионных состояний. На рис. 3, 4 изображены рассчитанные по программе PC GAMESS потенциальные кривые всех дублетных и квартетных состояний аниона 02 , коррелирующие с низшим по энергии диссоциационным пределом 0(2P) + 0(3P). Был использован ориентированный на расчеты с учетом корреляционных эффектов базис aug-cc-pVTZ, молекулярные орбитали найдены при оптимизации методом CASSCF в варианте полного пространства активных орбиталей с усреднением по всем рассматриваемым состояниям в каждом блоке симметрии (здесь классифицированных в рамках группы D2h). Расчеты методом конфигурационного взаимодействия в данном случае являются слишком дорогостоящими, поэтому корреляционные вклады более высоких порядков по отношению к уровню CASSCF были учтены с помощью процедуры многоконфигурационной квазивырожденной теории возмущений (MCQDPT) [13].

Полученная картина термов двухатомных систем хорошо согласуется с немногочисленными известными экспериментальными оценками молекулярных постоянных. Применение рассчитанных неэмпирическими методами потенциальных кривых систем OH, OH , 0H+, 02 в рамках оригинальной схемы DIM, учитывающей баланс ионных и нейтральных составляющих [14], позволило с хорошей точностью характеризовать кластеры (H20)n (n < 6). Полная картина потенциальных кривых аниона 02 представляет и самостоятельную ценность.

Еще один пример расчетов картины потенциальных кривых двухатомных молекул, выполненных

Е, эВ

1,4 1,8 2,2 2,6 3,0 3,4

Ло-О, А

Рис. 3. Потенциальные кривые дублетных состояний О2 [12]

Т а б л и ц а 2

Вертикальные энергии возбуждения (эВ) от точки минимума основного состояния Х2П молекулы 8И

Ссылка Л21+ 2Л B2I+ 22П 321+ 32П

[15], базис а^-сс-рУТ2 3,98 6,94 7,85 9,69 10,53 9,99

[15], базис АНО 3,96 - 7,73 9,74 - -

Эксперимент 3,85 - - - - -

Расчет, [19] 4,07/3,96 7,04 7,95 8,05 - -

Расчет, [20] 4,19 7,04 7,99 7,61 - -

Расчет, [21] 4,08 7,16 8,30 - - -

была успешно применена для молекулярно-динами-ческого моделирования спектров люминисценции радикалов 8Н в криптоновых матрицах [16, 17]. На рис. 5 изображены рассчитанные энергетические кривые дублетных состояний молекулы 8Н, коррелирующие с несколькими низшими диссоциацион-ными пределами, включая и ион-парный 8"Н+. Использовали стандартную вычислительную схему, а именно, оптимизировали молекулярные орбитали методом СЛ88СБ с усреднением по состояниям, а затем строили волновые функции конфигурационного взаимодействия второго порядка на основе исходных конфигураций СЛ88СБ. Функции СЛ88СБ строили при полном распределении 7 валентных электронов по 9 активным орбиталям; в методе 8ОС1 учитывали до 500 000 конфигурационных функций. Были рассмотрены два варианта базисных наборов АО, а именно, набор а^-сс-рУТ2, а также так называемые атомные натуральные орбитали (АНО) [787р4^ для серы и [683р] для водорода [18].

В данном случае интересно сравнить результаты расчетов с данными, известными из литературы [19-21]. Равновесные межъядерные расстояния для связанных состояний совпадают с экспериментальными с точностью до 0,01 А. Рассчитанные значения энергии диссоциации для низших связанных состояний имеют погрешности менее 0,1 эВ. В табл. 2 значения вертикальной энергии возбуждения, вычисленные в работе [15], сопоставлены с литературными данными.

Рис. 5. Потенциальные кривые молекулы SH [15]

для использования в рамках схемы DIM, относится к проблеме построения потенциальной поверхности трехатомной системы SH...Kr [15]. В последующих исследованиях эта поверхность приближения DIM

Рис. 6. Потенциальные кривые молекулы Кгё [15]

В отличие от 8Н для другого двухатомного фрагмента, необходимого при построении потенциальной поверхности 8Н ...Кг, а именно, для Кгё наши расчеты [15] предоставляли данные для построения потенциальных кривых ряда электронных состояний. Качественно полученная картина (рис. 6) согласуется с данными неэмпирических квантово-хи-мических расчетов родственных соединений ХеО и Хе8 [22]. В данном приложении мы использовали эффективный потенциал электронного остова для тяжелого атома Кг, построенный с учетом релятивистских поправок, а также рекомендованный базисный набор для валентных орбиталей [23]. Для серы был использован базис сс-рУБ2 с дополнительными наборами диффузных функций. В методе СЛ88СБ было учтено распределение 12 валентных электронов по 11 активным орбиталям. Функции приближения 8ОС1 насчитывали до 270 000 конфигураций. Сравнение вычисленных молекулярных постоянных для связанных состояний Кгё с литературными данными для Хе8 [22] подтвердило адекватность выбранной методики.

Более трудную задачу представляют расчеты двухатомных систем с участием атомов редкоземельных элементов. В работах [24, 25] были проведены оценки оптических спектров молекул Но2 и AgHo, и результаты напрямую сопоставлены с данными экспериментальных исследований в аргоновых матрицах. Для гольмия были использованы псевдопотенциалы штуттгартской группы [26] (с включенными в остов ^электронами) и соответствующими базисными наборами [584р3^. По данным расчетов методом С^Б для гомоядерной молекулы Но2 было найдено, что нижней по энергии потенциальной кривой отвечает электронная конфигурация 3^ 3си 2пи с формальным триплетным состоянием по спину. Это состояние является связанным (равновесное межъядерное расстояние 3,32 А) с энергией диссоциации 23,6 ккал/моль, в то время как экспериментальные оценки дают 19±4 ккал/мол. По расчетным данным электронным переходам в Но2 с наибольшими значениями дипольного момента перехода должны соответствовать возбуждения пи ^ ^ (500 нм) и пи ^ с* (524 нм). В экспериментальных спектрах наблюдаются интенсивные дублетные полосы при 498/504 нм и 558/563 нм. При расчетах гетероядерной системы AgHo [25] для гольмия были использованы те же псевдопотенциалы и базисы, а для серебра - достаточно популярный псевдопотенциал и базис Стивенса-Баша-Краусса (8ВК) [27]. Было показано, что смешанный кластер из ё- и ^элементов AgHo является стабильной частицей с равновесным межъядерным расстоянием 3,22 А, хотя и с очень пологой потенциальной кривой (гармоническая частота колебания 83 см х). По результатам расчетов наиболее интенсивная полоса в оптическом спектре находится при 436 нм. По экспериментальным данным при облучении смеси двух металлов в аргоновых матрицах наблюдается рост новой полосы при 430 нм, что свидетельствует об образовании соединения. Таким образом, согласие экспериментальных и рассчитанных спектров достаточно сложных двухатомных систем более чем удовлетворительное.

Примеры расчетов двухатомных молекул, выполненных с помощью программы РС ОЛЫБ88, показывают, что задача построения картины потенциальных кривых для взаимодействия пары атомов в любых электронных состояниях неэмпирическими методами квантовой химии практически решена. Адекватными средствами являются многоконфигурационные методы (СЛ88СБ и С1). Использование

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

2. Свойства и спектры малоатомных молекул

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

Расчеты равновесных геометрических конфигураций, колебательных спектров малоатомных молекул, а также анализ их электронной структуры составляют, по-видимому, наиболее обширную область применения программы РС ОЛМБ88 в настоящее время. Именно при решении таких задач использование программы требует наименьших усилий по заданию входной информации. Задача исследователя сводится к выбору методики (чаще всего НБ или МР2), выбору базиса и заданию некоторых опций программы, если параметры, заложенные по умолчанию, по каким-либо причинам не устраивают пользователя. После завершения вычислений достаточно однозначно интерпретируются результаты. Значительно реже рассматриваются свойства малоатомных молекул в возбужденных электронных состояниях и электронные спектры.

По-видимому, исторически первой публикацией результатов квантово-химических расчетов с применением программы, которая несколько позже стала называться РС ОЛМБ88, была работа [28], посвященная анализу строения молекул перлитиевых сульфидов Б^ и Ы48. Методом КНБ были определены равновесные геометрические конфигурации, оценены энергии диссоциации по различным каналам в приближении МР2/6-3Ш*//3-2Ш с учетом

вкладов нулевых колебаний. Показано, что в структурах сульфидов и оксидов гиперлитиевых систем, рассмотренных ранее [29, 30], есть как общие черты, так и отличия. Четко выделяется структурный фрагмент Б13Х с атомом неметалла в центре треугольника. Добавление четвертого атома лития приводит к комплексам симметрии С2у с электронным распределением [Б13Х]5+...Б15 , но для оксида существует более стабильный изомер симметрии Тй. В работе [28] была использована методика анализа электронной плотности, основанная на представлении натуральных связевых орбиталей (КВО) [31], которая является очень удобной для интерпретации электронной структуры в рамках концепции Льюиса [32]. Проведенные по схеме КВО рутинные расчеты (эффективные заряды на атомах, состав связывающих и несвязывающих орбиталей, а также орбита-лей неподеленных пар) позволяют для большинства молекул сформулировать результаты квантово-хими-ческих вычислений в привычных для химиков терминах. В частности, сделан вывод о преобладающем ионном характере связи для гиперлитиевых сульфидов, причем заряд на сере близок к формальной степени окисления -2.

Работа [33], выполненная с помощью первой версии программы РС ОЛМБ88, была посвящена структуре молекулы Л12О4 и фактически завершала серию расчетов строения и колебательных спектров малоатомных соединений алюминия с кислородом [34, 35]. Частицы состава Л1пОт (п + т < 5) были зафиксированы по ИК-спектрам как продукты реакции алюминия с кислородом в низкотемпературных матрицах [36]. По результатам неэмпирических расчетов методами КНБ [33], молекула Л12О4 является стабильной, может образовываться в экспериментальных условиях, и наблюдаемая интенсивная ИК-полоса при 590 см относится именно к этому соединению. По найденной равновесной геометрической конфигурации и анализу электронной структуры соединение Л12О4 можно характеризовать как устойчивый димер Л1О2. Согласно расчетам методом МР2, энергия димеризации Л1О2 составляет 128 ккал/моль.

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

ем, чтобы оно было выполнено с PC GAMESS), предпринятых с целью интерпретации экспериментов по матричной изоляции. Проиллюстрируем подобную деятельность на примере серии работ по расчетам магний-органических соединений. Образование молекул реактивов Гриньяра CH3MgX при низкотемпературных матричных реакциях магния с алкилгалогенидами CH3X (X=Cl, Br) является хорошо известным фактом. По результатам спектральных исследований [37, 38] продуктов подобных реакций была высказана достаточно интригующая гипотеза об образовании более сложных соединений состава CH3MgMgX. Квантово-химические расчеты, выполненные, в частности, и по программе PC GAMESS методами RHF и MP2 [39-41], подтвердили, что бимагниевые соединения CH3MgMgX (X=F, Cl, Br) должны быть стабильными частицами -энергия CH3MgMgX ниже на 7-8 ккал/моль (результаты расчетов в приближении MP2/6-31G*// RHF/6-31G*), чем сумма энергий CH3MgX + Mg. ИК-спектры моно- и бимагниевых соединений должны заметно различаться. Анализ электронной структуры с натуральными связевыми орбиталями может быть интерпретирован таким образом, что в бимагниевых молекулах образуется связь Mg-Mg. В случае фтора стабильными должны быть также комплексы присоединения одного или двух атомов магния к молекуле CH3MgF со стороны фтора.

Взаимодействие магния с тетрахлоридом углерода Mg + CCl4 исследовали экспериментально при низких температурах в работе [42]. О протекании реакций свидетельствовало появление новых полос в ИК-спектрах, однако уверенно отнести спектральные полосы к возможным металло-органическим соединениям не удалось. В сравнительно недавней работе [43] квантово-химические расчеты возможных продуктов реакции Mg + CCl4 были выполнены с помощью программы PC GAMESS на вполне современном уровне. На первом этапе были найдены равновесные геометрические конфигурации молекул CCl4, MgCl2, CCl2, CH3Cl, CH3MgCl, CCl3MgCl, CCl2MgCl2, CCl3MgMgCl при полной оптимизации параметров методом MP2 с корреляци-онно-корректированными базисами aug-cc-pVTZ. Расчеты частот гармонических колебаний подтвердили, что полученные структуры соответствуют истинным минимумам на потенциальных кривых. На следующем этапе для найденных конфигураций вычислены энергии в приближении теории возмущений четвертого порядка MP4-SDTQ с теми же бази-

Рис. 7. Равновесные геометрические конфигурации молекул CCl2MgCl2 (а) и CClзMgMgCl (б). Межъядерные расстояния указаны в ангстремах, углы - в градусах

сами а^-сс-рУТ2. Для уточнения значений относительной энергии были введены поправки на энергию нулевых колебаний. Точность расчета можно характеризовать при сравнении с известными экспериментальными справочными данными [44], в частности, для ^^ получена конфигурация симметрии Тй с равновесным расстоянием ^О 1,77 А (справочные данные: 1,766 ± 0,005А), для молекулы MgCl2 получена линейная структура симметрии с расстоянием Mg-Cl 2,19 А (2,186 ± 0,011А [44]). Аналогом реактива Гриньяра является соединение CCl3MgCl симметрии Наиболее интересные результаты связаны с предсказанием соединений необычного состава, а именно, комплексов CCl2MgCl2 и CCl3MgMgCl (рис. 7).

В табл. 3 сопоставлены значения энергии связи соединений CCl3MgCl, CCl2MgCl2, CCl3MgMgCl по отношению к различным каналам распада, вычисленные в приближении МР4/а^-сс-рУТ2//МР2/а^-сс-рУТ2 с учетом поправок на энергии нулевых колебаний. Видно, что образование всех трех продуктов реакции Mg + CCl4 энергетически выгодно. Для сравнения приведем энергии реакций CH3MgCl ^ Ш3а + Mg (52,4 ккал/моль в работе [40]) и CH3MgMgCl ^ Ш^а + Mg (6,2 ккал/моль

Т а б л и ц а 3

Энергии связи CClзMgCl, СС12М§С12, CClзMgMgCl по отношению к различным каналам распада, вычисленные в приближении МР4/а^-сс-рУТ2//МР2/а^-сс-рУТ2 с учетом поправок на энергии нулевых колебаний

Реакционный канал Энергия, ккал/моль

CCbMgCl ^ CCI4 + Mg 67,0

CCl2MgCl2 ^ CCI4 + Mg 70,7

CCl2MgCl2 ^ CCbMgCl 3,7

CC^MgCh ^ CCl2 + MgCl2 21,0

CCbMgMgCl ^ CCl4 + Mg + Mg 80,6

CCbMgMgCl ^ CCbMgCl + Mg 13,6

[40]). Таким образом, замещение атомов водорода на хлор в реагенте СН3С1 способствует более энергетически выгодному образованию магний-органических соединений.

Можно ожидать, что использование методов теории возмущений (вплоть до 4 порядка) для решения задач с числом базисных функций 350-400 позволяет для малоатомных молекул получать достаточно точные оценки молекулярных постоянных - с погрешностями не более 0,01 А по равновесным межъядерным расстояниям и порядка 1 ккал/моль по энергиям связи.

Примером достаточно успешного соотнесения экспериментальных данных по спектрам матрично-изолированных магний-органических соединений с результатами квантово-химических расчетов является работа [45]. Продукты матричных реакций атомов и кластеров магния с диоксидом углерода, этиленом и их смесями были исследованы спектральными методами. Квантово-химические расчеты равновесных геометрических конфигураций и гармонических частот колебаний выполнены методом МР2/ 6-31Ш(2ё,2р) для соединений MgCO2, MgC2H4, Mg(C2H4)2, Mg(C2H4)(CO2). Для всех перечисленных систем были найдены минимумы на соответствующих поверхностях потенциальной энергии. Во всех случаях для этих конфигураций характерно наличие циклических фрагментов с участием атома магния. Согласие экспериментальных полос в ИК-спектре матрично-изолированной молекулы и вычисленных

гармонических частот колебаний можно проиллюстрировать для MgCO2: 1580, 1385 и 866 см- (эксперимент), 1700, 1385, 840 см- (расчет). Качественное (или почти-количественное) совпадение положений полос и интенсивностей с результатами кванто-во-химического моделирования позволило во всех случаях установить образование циклических комп -лексов в матрицах.

В работах ряда исследовательских групп можно найти другие примеры применения программы PC GAMES S для предсказания структуры и спектров молекул, условно относимых к малоатомным. Отметим, в частности, расчеты строения замещенных бензилхлоридов [46], расчеты строения молекулы бис-(триметилстаннил)ацетилена, а также анализ колебательных спектров его изотопомеров методами RHF и MP2 [47]. Расчеты колебательных спектров малоатомных молекул высокого уровня точности с использованием PC GAMESS выполнены в группе В.Г. Соломоника (Ивановский государственный химико-технологический университет) [48-51]. При анализе молекул состава MX3 (M=Sc, Y, La, Lu; X=F, Cl, Br, I) использованы широкие наборы базисных функций (вплоть до f-функций) атомов металла и рассмотрены варианты решения электронного уравнения с достаточно полным учетом электронной корреляции (например, варианты CISD+Q). Показана необходимость таких усилий (а также и ангармонических вкладов в колебания) для четкого отнесения экспериментальных спектров [49, 50]. Тщательно рассмотрены строение и спектры молекул солей типа M2Xa,bO4 и MReO4 (M=Li, Na, K; Xa = S, Se, Te; Xb = Cr, Mo,W) [51, 52]. В группе И.А. Годунова (химический факультет МГУ им. М.В. Ломоносова) методами HF, CASSCF, CI рассчитывали поверхности потенциальной энергии в основных и возбужденных электронных состояниях малоатомных карбонильных и нитрозосоединений [53-55]. На основании этих данных определяли равновесные геометрические конфигурации, частоты колебаний, потенциальные функции внутреннего вращения, а также вращательные контуры полос электронно-колебательных спектров.

К последним работам [49-55], конечно, не относится сформулированное в начале данного раздела утверждение о «простоте» расчетов структуры и колебательных спектров малоатомных молекул с помощью программы PC GAMES S. Если исследователя устраивает точность оценок для основного электронного состояния небольших молекул, состо-

ящих из атомов, расположенных в верхней части таблицы элементов (порядка 0,01 А по равновесным межъядерным расстояниям, 5° по валентным углам, 50 см по частотам колебаний), то достаточно простые расчеты вполне возможны с использованием методов HF, MP2. Если поставлены задачи, решение которых требует большей точности, а также если в молекуле есть тяжелые атомы или необходимы сведения о свойствах системы в возбужденных электронных состояниях, то решение требует значительно больших усилий, хотя принципиально возможно.

3. Кластеры и межмолекулярные комплексы

В этом разделе приведены результаты моделирования с помощью программы PC GAMESS свойств систем, которые относят к межмолекулярным комплексам или к атомным и молекулярным кластерам. Хотя достаточно трудно сформулировать четкое определение для такого химического образования, как кластер [56], мы будем пользоваться этим термином в тех случаях, когда частица существует как система связанных ядер и электронов, но расставить на соответствующей структурной формуле валентные связи однозначно не удается.

Самые первые случаи применения PC GAMESS к молекулярным ассоциатам сводились к квантово-химическому анализу электронного распределения в системах, геометрические конфигурации которых находили с помощью эмпирических потенциалов. В работе [57] рассматривали свойства молекулы NBr в сольватных оболочках аргона с целью интерпретации спектроскопических экспериментов в условиях низкотемпературной матричной изоляции. Колебания в системе анализировали методом молекулярной динамики, но оценки влияния матрицы на электронные характеристики молекулы NBr в основном (X3£) и возбужденном (b1^) состояниях проводили методами квантовой химии. С этой целью рассчитывали электронные распределения в кластере NBr@Arn методами HF, MP2 с базисом 6-31G*, используя преобразования к натуральным связевым орбиталям. Сделан вывод, что следствием взаимодействия с матрицей могут быть достаточно заметные сдвиги в значениях орбитальной энергии связывающих орбиталей и орбиталей неподеленных пар (вплоть до 300 см х). Похожий анализ был выполнен для димера фторида водорода (HF)2, включенного в аргоновые оболочки, содержащие до 62 атомов инертного газа [58]. Найдено, что окру-

жающие аргоновые атомы могут отдавать на соединение внедрения небольшую долю электронного заряда, которая распределяется на несвязывающие орбитали c*(H-F) димера фторида водорода.

Свойства малых кластеров металлов, как гомоя-дерных (M)m, так и гетероядерных (Mx)m(M2)n, представляют большой интерес с точки зрения фундаментального знания, а также в связи с многообещающей перспективой их практического применения [59]. В работах [60-62] рассчитывали равновесные геометрические конфигурации энергии связи кластеров (M1)m(M2)n, где Мх, М2 = Ag, Cd, Cu, Mg, Na, Pb, Sn, Zn (m + n < 4). Для каждого кластера был определен тип электронного состояния (по спину) наинизшей энергии. Для выделенной мульти-плетности методами RHF и UHF рассчитывали геометрическую конфигурацию, отвечающую минимуму на поверности потенциальной энергии. Во всех случаях были использованы эффективные потенциалы остова и соответствующие базисы для валентных орбиталей в приближении SBK. Расчеты частоты колебаний позволили проверить, действительно ли найденная конфигурация отвечает истинной точке минимума. Для найденных координат вычисляли значения энергии систем по теории возмущений MP2. Эти данные позволяют заключить, что наиболее перспективными парными композициями являются следующие: AgCu, CuPb, CuSn, AgPb, CuNa, PbSn, AgSn, NaSn с энергиями связи 20-30 ккал/моль. Исследование смешанных трех- и четырехатомных кластеров было проведено с учетом результатов, полученных для двухатомных частиц. Прежде всего интерес представляли Ag и Cu как партнеры в смешанных кластерах типа AgnMm и CunMm, именно такие комбинации и были рассмотрены. Следует отметить, что для систем с числом атомов больше двух существует большое число изомерных структур в дополнение к разнообразию спиновых состояний. Окончательно было установлено, что наиболее перспективными для образования смешанных металлических кластеров являются комбинации серебра или меди со свинцом или оловом, т.е. этим частицам соответствуют наибольшие значения энергии связи по отношению к разным каналам диссоциации.

Строение более привычных комплексов металлов с молекулярными лигандами рассматривали в работах [63, 64]. Для интерпретации экспериментов по конденсации паров гольмия в низкотемпературных матрицах диоксида углерода и азота были выполне-

Рис. 8. Равновесные геометрические конфигурации: а - Но(СО)6, б - Но(М2)6 (межъядерные расстояния указаны в ангстремах)

ны расчеты по программе PC GAMESS для ряда соединений состава HoLn (L = CO, N2; n = 1-6). Поиск равновесных геометрических конфигураций и расчет частот гармонических колебаний проводили методом HF с использованием псевдопотенциала и базиса [4s4p2d2f] для атома гольмия, предложенного в [65], и базиса 3-21G для углерода, кислорода и азота. Был использован общепринятый прием масштабирования вычисленных частот колебаний. Параметр масштабирования (0,91) был подобран в серии предварительных расчетов родственных молекул Cr(CO)6, W(CO)6 таким образом, чтобы воспроизвести известный экспериментальный колебательный спектр. По результатам оценок колебательных полос было найдено, что среди соединений типа Ho(CO)n, Ho(N2)n с экспериментальными матричными данными по ИК-спектроскопии согласуются лишь структуры гексамеров. Соответствующие конфигурации показаны на рис. 8, причем для Ho(CO)6 характерна симметрия D2d, а для Ho(N2)6 - симметрия D4h*

Много работ [66-70] посвящено моделированию свойств олигомеров фторида водорода (HF)n (n > 6). Кластеры фторида водорода представляют собой, по-видимому, наиболее простые системы с сетями водородных связей, и поэтому их интенсивно исследуют как экспериментально, так и теоретически. Для них характерно большое разнообразие изомерных структур в пределах кластера одного состава, энергии которых различаются несколькими едини-

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

Анализ этого материала в рамках данного обзора подробно не проводится, поскольку наиболее интересные результаты по равновесным геометрическим конфигурациям глобальных и локальных минимумов были получены не с помощью программы PC GAMESS, а более дешевым методом двухатомных фрагментов в молекулах с применением алгоритмов поиска стационарных точек [71], не включенных в PC GAMESS. Однако проведение этих работ без использования PC GAMESS было бы невозможно. Во-первых, множество двухатомных потенциальных кривых, необходимых для применения приближения DIM, включая и кривые ион-парных состояний (раздел 2 настоящего обзора), были рассчитаны с помощью PC GAMESS. Во-вторых, предсказания метода DIM проверяли выборочными расчетами структур надежным квантово-химическим методом MP2/6-31+G(2d,2p). В-третьих, распределения электронной плотности для всей совокупности структур анализировали по методу натуральных связевых орбиталей, что позволило сформулировать ряд заключений о системе водородных связей в данных кластерах.

Циклические структуры олигомеров хлорида водорода (HCl)n (n = 3, 4, 5), отвечающие глобальным минимумам на поверхностях потенциальной энер-

Рис. 9. Равновесные геометрические конфигурации циклического глобального (а) и некоторых локальных

(б-е) минимумов (НБ)6

гии каждого кластера данного состава, были рассчитаны методом MP2 с несколькими вариантами выбора базиса, включая aug-cc-pVTZ [72]. Найденные равновесные геометрические конфигурации, значения энергии связи и наборы гармонических частот колебаний были сопоставлены с соответсву-ющими результатами для кластеров фторида водорода. Показано, что в комплексах энергия связи, приходящаяся на одну водородную связь, монотонно увеличивается при увеличении размера кластера до значений 7,5 и 3,3 ккал/моль для (ИР)5и (HCl)5 соответственно.

Если изучение структуры гомогенных кластеров типа (HF)h рассматривается как шаг на пути моделирования свойств жидкости, то исследование молекулярных гетерокластеров вносит вклад в моделиро-

вание растворов на микроуровне. Сведения о строении сольватных оболочек, ближайших к частицам растворенного вещества, необходимы на всех стадиях изучения свойств растворов. В работе [73] были рассчитаны стационарные точки на поверхностях потенциальной энергии кластеров О (НБ)п (п = 1 -5) и проанализировано расположение молекул «растворителя» НБ относительно «растворенного вещества» О . Полная оптимизация геометрических параметров для всех кластеров проведена с помощью программы PC GAMESS методом MP2/6-311+G**. Расчетом гармонических частот колебаний подтверждено, что найденные точки отвечают истинным минимумам на потенциальных поверхностях. Энергии для выделенных конфигураций пересчитаны в приближении MP4(SDTQ)/6-311+G**. Особый инте-

Е, эВ

6-

5-

4-

3-

2-

1 -

0-

-1-

п

Рис. 10. Энергия поверхностной структуры кластера С1 (НР)п по отношению к внутренней в зависимости от размера кластера п

рес представляло сравнение двух типов структур, характерных для кластеров СГ(НР)п (п = 2, 3, 4, 5): «поверхностных», с расположением аниона хлора на поверхности кластера (НБ)п, и «внутренних», когда анион окружен молекулами растворителя. На рис. 10 показана зависимость разности энергий поверхностной и внутренней структур в зависимости от размера кластера п. Видно, что для п < 5 более выгодными являются конфигурации, в которых анион хлора расположен на поверхности, но по мере роста кластера внутренние структуры становятся более конкурентоспособными.

Понятно, что наибольший интерес вызывает исследование строения межмолекулярных кластеров с участием молекул воды. В работе [74] были выполнены оценки сольватных сдвигов в колебательных спектрах тиокомплексов золота по данным квантово-химических расчетов по программе РС ОЛМБ88. Методами Хартри-Фока с двух- и трех-экспонентны-ми базисами и использованием псевдопотенциалов для атома золота были найдены равновесные геометрические конфигурации и вычислены гармонические частоты колебаний следующих систем: Ли8Н, Ли(8Н)2-, Ли(8Н)(Н28), Ли8Н(Н20)1-8, Ли(8Н)2-(Н20)4. На основании расчетов можно проследить характерные эффекты, показывающие влияние растворителя (воды) на колебательный спектр растворенного вещества: сдвигу в положениях полос в сторону больших частот (до +33 см для валентных колебаний 8-Н) и в сторону меньших частот (до -19 см для деформационных колебаний Н-Ли-8).

В работе [75] были оценены сдвиги частот коле -баний молекулы Н28 при ее внедрении в кластеры воды. С этой целью были рассчитаны структуры систем Н28(Н20)57 с использованием метода МР2 с двухэкспонентным базисом, дополненным поляризационными функциями. По результатам расчетов кластера наибольшего размера Н^(Н20)57 влияние растворителя приводит к сдвигам в частотах +23 и -265 см 1 для симметричного и асимметричного колебаний молекулы Н28.

Структура малых кластеров воды (включая анионы и катионы как в чистом виде, так и допирован-ные малыми неорганическими молекулами) активно исследуется с помощью неэмпирических квантово-химических расчетов в лаборатории строения и квантовой механики молекул химического факультета МГУ им. М.В. Ломоносова [75].

Так же, как и при моделировании олигомеров фторида водорода [67, 68], в нашей работе [12] равновесные конфигурации изомеров малых кластеров воды (Н20)п (п = 2-6) рассчитывали методом двухатомных фрагментов в молекулах. Как описано в разделе 1, необходимые потенциальные кривые пар атомов были построены с помощью РС ОЛМБ88. Для найденных структур (рис. 11) неэмпирическим методом МР2 были достаточно точно вычислены разности энергий изомеров для каждого кластера, а также энергии связи кластеров (Н20)п по отношению к п(Н20).

В работах [77-80] рассчитывали свойства комплексов серебра с фрагментами цианобифенилов. Реакции переходных металлов с органическими ме-зогенными молекулами представляют особый интерес ввиду возможного их использования в синтезе новых материалов на основе жидких кристаллов. Подобные жидкокристаллические материалы обладают важными электрофизическими характеристиками и находят применение в оптоэлектронике. Кван-тово-химические расчеты структуры и спектров модельных комплексов серебра с центральным циано-фенильным фрагментом, входящим в состав мезо-генных цианобифенилов и ряда их производных, позволяют напрямую сопоставить результаты с экспериментальными данными по низкотемпературной соконденсации 4-пентил-4'-цианобифенила (5СВ) с атомарным серебром в инертных матрицах. Молекула цианобифенила

С5Н

11

Рис. 11. Равновесные конфигурации кластеров воды (И20)п: а - димер; б и в - структуры (ииф и (шш) тримера; г и д -структуры (и^ф и (и^ф тетрамера; е - циклический пентамер; структуры гексамера: ж - призма, з - клетка,

и - кольцо и к - книга

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

С

Ад

I

(1)

N

(2)

Р

Рис. 12. Геометрия модельной системы (PhCN)2Ag

структуры комплекса 5СВ с серебром, соответствующие п-комплексам. Согласно данным ренгено-структурного анализа молекулы производных циано-бифенилов в кристаллах упакованы в пары по принципу "голова-к-хвосту", поэтому рассматривалась модель комплекса серебра с двумя молекулами цианобензола, представленная на рис. 12.

Система обладает плоской симметрией с атомом серебра, помещенным на равных расстояниях от лигандов, а именно, от центра CN фрагмента ли-ганда (1) и от центра бензольного кольца лиганда (2). В расчетах почти все геометрические параметры были фиксированы и равны значениям аналогичных параметров в свободной молекуле РЬС^ Варьированию подлежало расстояние Я между плоскостями лигандов (1) и (2).

Расчет методом многоконфигурационного самосогласованного поля СА88СБ с двумя электронными конфигурациями показал, что два решения могут конкурировать в основном состоянии: одно соответствует (PhCN)2Ag, а другое - конфигурации с переносом заряда (РЬС^-Ag+. Красный сдвиг частот колебаний С^группы в комплексах серебра с цианобифенилом качественно понятен. При образо-

вании пленки атомы серебра входят в полости между парами органических лигандов таким образом, что достигается асимметричное окружение по отношению к CN-группам. Донирование небольшой части электронной плотности с серебра на антисвязы-вающие орбитали CN одного лиганда приводит к возрастанию соответствующего межъядерного C-N-расстояния и уменьшению кривизны потенциальной поверхности вдоль этой координаты, т.е. к уменьшению частоты колебаний. Донирование большей части заряда с Ag на п*-систему другого лиганда ведет к перераспределению электронной плотности в образовавшемся отрицательном ионе с уменьшением соответствующей частоты CN-колебаний на другую величину. Заключения качественного характера подтверждаются численными данными. В случае комп -лекса с переносом заряда вычисленные красные сдвиги в частотах колебаний CN (-150 и -175 см 1) хорошо коррелируют с экспериментально наблюдаемыми для Ag/5СВ пленок (-150 и -200 см 1).

В заключение этого раздела процитируем работу [81], в которой методами квантовой химии по программе PC GAMESS рассчитывали заряды на атомах в пептидах и комплексах пептидов с водой. Эти данные необходимы, в частности, для калибровки потенциалов молекулярно-динамических расчетов свойств пептидов.

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

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

4. Расчеты энергетических профилей химических реакций

В первой работе [82], использующей PC GAMESS для расчетов энергетических профилей химических превращений, была рассмотрена реакция

Li2 + H2S ^ Li2S + H2.

Интерес к этой системе объясняется тем, что образование Li2S является, по-видимому, одним из каналов превращения кластеров лития в реакциях с сероводородом, а также тем, что в ходе данной реакции происходит переход от соединений с ко-валентными связями (Li2, H2S) к молекуле с ионной связью (Li2S). В расчетах энергетического профиля реакции был использован эффективный как в вычислительном, так и в интерпретационном отношении прием локального учета эффектов электронной корреляции CAS/NBO [83, 84]. Согласно этой методике многоконфигурационные волновые функции CASSCF построены на натуральных связевых орбиталях, четко сопоставляемых с химическими связями. В данном случае волновые функции были построены при распределении 6 валентных электронов по орбиталям c(Li-Li), c(S-Hx), c(S-H2), c*(Li-Li), c*(S-Hx), C*(S-H2). Таким образом, учитывались заселения электронами связывающих орбиталей a (Li-Li), a(S-Hj), a(S-H2), дважды занятых у реагентов, и несвязывающих орбиталей a*(Li-Li), a*(S-Hx), a*(S-H2), что моделирует разрывы соответствующих связей. В долине реагентов такая схема достаточно адекватно передает изменения электронной структуры участников реакции. В работе [82] было показано, что вплоть до точки перевала на энергетической кривой структурные изменения, а также изменения в электронной конфигурации происходят лишь с фрагментом Li2. Использованный базис атомных орбиталей (3-21G) позволяет делать заключения преимущественно качественного характера, поэтому полученное значение акти-вационного барьера (8 ккал/моль) можно рассматривать как достаточно приближенную оценку.

Представление натуральных связевых орбиталей было также активно использовано при моделировании активации молекулы этилена структурными

фрагментами катализатора алюмофенилсилоксана [85, 86]. Электронная конфигурация молекулы этилена, описываемая конфигурацией [Остов]а (С-С)а2(С-Н1)а2(С-Н2)а2(С-Н3)а2(С-И4)л2(С-С), при сближении с катализатором претерпевает изменения. В качестве критерия активации рассматривался такой вычисляемый методами квантовой химии параметр, как заселенность орбитали п(С-С). Уменьшение заселенности от исходного значения 2 свидетельствует о готовности к разрыву п-связи в этилене. Рассмотренные фрагменты алюмофенилсилокса-на содержали центральную структурную группу А104 с различными замыкающими атомами. По результатам моделирования был сделан вывод, что наибольшее значение для активации имеет правильная тетраэдрическая геометрическая конфигурация группы А104.

Достаточно активно квантово-химические расчеты применяют для моделирования превращений органических соединений. В частности, на основе расчетов поверхности потенциальной энергии анион-радикала С6Н5С1' проанализированы возможные механизмы гетеролитического разрыва связи С-На1 [87, 88].

Примером расчета, выполненного на достаточно высоком уровне, является построение профиля пути реакции для взаимодействия атомарного магния с диоксидом углерода [89]. Выше упоминалась работа [41], в которой описаны результаты комплексного исследования MgC02 неэмпирическими методами квантовой химии и низкотемпературной матричной изоляции. Полная оптимизация геометрии MgC02 методамами ОУВ, МР2 в базисе 6-31Ш(2^2р) без явного наложения ограничений по симметрии приводила к структуре С2у типа п -0,0 с нелинейным фрагментом С02, в котором длины связей С0 незначительно увеличены. Расчеты гармонических частот колебаний, выполненные методами ОУВ и МР2, подтверждают, что найденная конфигурация действительно соответствует минимуму поверхности потенциальной энергии системы. Рассчитанные высокочастотные полосы в колебательном спектре комплекса коррелируют с результатами исследований методами спектроскопии матричной изоляции. Однако сравнение значений энергии, рассчитанных в приближении МР2/6-31 Ш(2^2р) для MgC02 в точке равновесия и для реагентов Mg + С02, показало, что данный минимум поверхности потенциальной энергии лежит выше диссоциационного предела более чем на 10 ккал/моль. Эти результаты оз-

начают, что подобная структура соответствует мета-стабильной частице, фиксируемой в матрице инертного газа и стабилизированной благодаря влиянию матрицы. Для проверки этого утверждения были проведены расчеты энергетического профиля взаимодействия Mg + CO2 на более высоком теоретическом уровне. Точки на поверхности потенциальной энергии рассчитывали с помощью многоконфигурационной квазивырожденной теории возмущений второго порядка (MCQDPT) с использованием кор-реляционно-корректированных базисов d-aug-cc-pVTZ на O, и aug-cc-pVTZ на C и Mg. Низшие по энергии 12 орбиталей были отнесены к остовным, в то время как 13 валентных орбиталей, заселяемых 10 электронами, составляли активное пространство для процедуры CASSCF, являющейся первой стадией подхода MCQDPT.

Межъядерное расстояние Mg-C было выбрано в качестве координаты реакции и для каждой точки вдоль этой координаты была проведена неградиентная оптимизация углов ZOCO и длин связей C-O по методу MCQDPT. Результаты расчета приведены на рис. 13. Внутренний минимум соответствует метастабильной структуре C2v, равновесные геометрические параметры для которой следующие: R(Mg-C) = 2,31 Á; ZOCO = 125,6°; R(C-O) = 1,26 Á. Расчеты наиболее высокого уровня, выполненные в данной работе (MCQDPT/ aug-cc-pVTZ), показывают, что энергия MgCO2 в равновесной точке выше диссоциационного предела Mg + CO2 на 14,7 ккал/моль.

Е 18 -i 161412108 6 4 2 0

-2 —i—1—i—1—i—1—i—1—i—1—i—1—i—1—i—1—г 2,0 2,5 3,0 3,5 4,0 4,5 5,0 5,5 6,0

R, Á

Рис. 13. Энергетический профиль пути реакции Mg + CO2 ^ MgCO2, рассчитанный в приближении MCQDPT/aug-cc-pVTZ

Резкое изменение характера потенциальной кривой в области 2,5 А (рис. 13) требует дополнительного обсуждения. Область, соответствующая координатам R(Mg-C) 2,0-2,5 А, коррелирует с ионным каналом диссоциации комплекса MgCO2 на Mg+ и СО 2. Часть кривой, оставшаяся после прохождения барьера, соответствует отталкивательному взаимодействию нейтральных частиц Mg и СО2. Полная диабатическая потенциальная кривая, соответствующая нейтральным реагентам Mg + СО2, может быть непрерывно продолжена в сторону меньших расстояний R(Mg-C) (<2,5 А), например, с волновыми функциями MCSCF. Однако кривая, соответствующая ионной паре Mg+ + CO 2, не имеет гладкого продолжения вправо от точки пересечения. Это объясняется автоионизацией CO 2 при определенных геометрических параметрах, где энергия ионной частицы превышает энергию CO2 и свободного электрона. Подобное поведение потенциалов обсуждалось и для двухатомной системы OH-/OH (раздел 2). Таким образом, между внутренней потенциальной ямой и

внешней областью энергетического профиля на рис. 13 не должно быть гладкого перехода.

Основные результаты работы [90], посвященной моделированию межмолекулярных комплексов HXeOH-(H2O)n (п = 1, 2, 3), можно было бы обсудить и в предыдущем разделе. Однако следует рассмотреть химические превращения этих необычных соединений, связанные с передачей протона по ориентированным цепям молекул воды. Ряд соединений благородных газов состава HRgY, где У - электроотрицательный фрагмент, были синтезированы методами матричной криохимии в университете Хельсинки (Финляндия) [91]. В частности, была зафиксирована частица HXeOH, образующаяся при лазерном облучении водно-ксеноновой смеси в аргоновых матрицах с последующим отжигом матрицы [92]. Полосы ИК-спектра, отнесенного к данному соединению, полностью согласуются с результатами квантово-химических расчетов. Система HXeOH, представляющая собой ионную пару HXe+...OH , относится к локальному минимуму на потенциаль-

Рис. 14. Рассчитанный методом MP2 путь реакции распада комплекса HXeOH-(H2O)2 ^ Xe-(H2O)3 [90]. Показаны конфигурации метастабильного состояния HXeOH-(H2O)2 (а), переходного состояния (б), продукта реакции HXeOH-(H2O)2 (д), и промежуточных структур на пути передачи протона (в) и (г)

Рис. 15. Энергетический профиль реакции ОН + С02 ^ НС03. Координатой реакции выбрано расстояние Я между кислородом гидроксила и углеродом С02; остальные внутренние координаты оптимизированы для каждого значения Я. Профиль в газовой фазе рассчитан в приближении МР2/6-311++0**. Профиль в водном окружении рассчитан для системы в окружении совокупности молекул воды, описываемых эффективными фрагментами. Положения молекул воды также оптимизированы для каждого значения координаты реакции

ной поверхности, и основному измеру соответствует комплекс Хе...Н20.

В работе [90] методом МР2 с использованием базисов а^-сс-рУТ2 на кислороде и водороде, штуттгартского псевдопотенциала и соответствующего базиса (б8брЫ)/[484рЫ] на ксеноне были найдены конфигурации локальных и глобальных минимумов на потенциальных поверхностях комплексов НХеОН с молекулами воды. Было показано, что комплексам НХе0Н-(Н20)п (п = 1, 2, 3) соответствуют связанные состояния систем. Дальнейшее увеличение числа молекул воды приводит к распаду комплекса. Расчетом частот и ИК-интенсивности колебаний было показано, что неотнесенные ранее экспериментальные полосы в спектрах матрично-изолированных частиц соответствуют комплексам с одной (НХе0Н-(Н20)) и двумя (НХе0Н-(Н20)2) молекулами воды. Наиболее интересные результаты были получены при моделировании пути распада комплексов, т. е. перехода систем от локальных минимумов к глобальным. Для этого были найдены конфигурации переходных состояний как седловых

точек с единственной мнимой частотой колебаний, и в обе стороны от этих точек построены пути минимальной энергии в долины реагентов и продуктов (опция внутренней координаты реакции (IRC) в программе PC GAMESS).

На рис. 14 показан рассчитанный энергетический профиль, соответствующий распаду комплекса HXeOH с двумя молекулами воды. Слева показана конфигурация локального минимума (метастабиль-ного состояния) HXeOH-(H2O)2 (структура а), на вершине барьера изображена конфигурация переходного состояния (структура б),справа внизу -структура продукта распада Xe + 3H2O (структура д). Промежуточные структуры (в) и (г) иллюстрируют процесс миграции протонов по системе водородных связей. При движении системы слева направо сразу же после переходного состояния протон от фрагмента XeH+ переходит на ближайшую молекулы воды (структура а). Другой протон от образовавшейся частицы H3O+ движется к следующей молекуле воды (структура г), от которой еще один протон перемещается к гидроксилу OH-. Были проведе-

Рис. 16. Геометрические конфигурации системы (ОН + С02 ^ НС03 )@(Н2О)30 для значений координаты реакции Я = 2,75 А (а), Я = 2 А (б), Я = 1,4 А (в). Стержнями обозначены молекулы воды, входящие в ближайшую оболочку сети водородных связей реагентов. Остальные молекулы воды показаны линиями

ны преобразования от молекулярных орбиталей к натуральным связевым орбиталям вдоль пути реакции. Приведенные рассуждения являются следствием анализа в терминах КБО, а не логическими построениями. Данный механизм переноса протона по ориентированным цепям молекул воды иллюстрирует молекулярную схему, которую принято называть водным (или протонным) проводником [93].

В работах [94-96] рассмотрена реакция ОН- + С02 ^ НС03 , энергетический профиль которой в газовой фазе принципиально отличается от энергетического профиля в водном окружении (рис. 15). В

отсутствие среды энергия сближающихся реагентов понижается вплоть до точки образования продукта НС03 . В водном окружении на пути реакции существует потенциальный барьер, на наличие которого указывают как экспериментальные данные, так и расчеты в рамках континуальных моделей среды. Для учета молекул среды в дискретных или кластерных приближениях в программе ОЛМБ88(и8) реализован метод потенциалов эффективных фрагментов [97]. В этом подходе каждая молекула растворителя является жесткой подсистемой (эффективным фрагментом) с фиксированными внутренни-

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

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

Если сравнить расположение частиц до барьера (рис. 16, а) и на барьере (рис. 16, б), то видно, что в ближайшей гидратной оболочке ОН происходят изменения - до барьера гидроксил-анион окружен большим количеством молекул воды, чем на барьере (5 против 4), и соответственно энергия такой системы оказывается ниже. В работе [95] был оценен барьер свободной энергии величиной 8-13

ккал/моль, что хорошо согласуется с экспериментальными данными.

Еще один пример, демонстрирующий значение учета влияния молекул растворителя на энергетические профили химической реакции, описан в работе [98]. На рис. 17 изображены сечения поверхностей потенциальной энергии кластера (HF)5 вдоль координаты, соответствующей диссоциации выделенной молекулы HF из пентамера. Остальные 4 молекулы кластера, окружающие выделенную, рассматриваются как молекулы «растворителя». Для каждого значения координаты реакции, т.е. межъядерного расстояния Rhf в диссоциирующей молекуле HF, положения всех частиц в кластере оптимизировались по минимуму энергии. Нижнее электронное состояние системы коррелирует с диссоциацией молекулы HF на основные состояния атомов H+F. Возбужденное электронное состояние коррелирует с ион-парным пределом H+ + F . Потенциальные поверхности пентамера в основном и возбужденном электронном состояниях строились методом DIM (см. раздел 3) с неэмпирическими потенциальными кривыми двухатомных фрагментов, включая и ион-парные состояния, вычисленными по программе PC GAMES S многоконфигурационными приближениями. Предсказания DIM для

Рис. 17. Потенциальные кривые диссоциации молекулы НБ в изолированном состоянии и в кластере (НБ)

Рис. 18. Структуры реагентов (а) и продуктов (б) реакции модельного фрагмента белка NCp7 с протонированной по центральному

атому азота формой 1,2-диазенендикарбоксамида

кластера (HF)5 также проверяли неэмпирическими квантово-химическими расчетами с использованием PC GAMESS.

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

Важное направление квантовой химии представляет расчет и исследование потенциальных поверхностей больших молекулярных систем биохимической значимости. Мы приведем два примера подобных расчетов, выполненных с использованием программы PC GAMESS. В работе [99] были построены энергетические профили вдоль предполагаемого пути реакции взаимодействия 1,2-диазенендикарбок-самида с цинксодержащим фрагментом белка ядер-

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

В работе [99] была предложена модельная структура для цинксодержащего фрагмента, показанная в верхней левой части рис. 18, причем для физиологических значений pH цистеиновые остатки должны быть депротонированными. Для реагента 1,2-диазе-нендикарбоксамида были рассмотрены разные варианты протонирования (по центральным или перефе-рийным атомам азота, или по кислороду). Далее методом HF с псевдопотенциалами и базисами SBK строили энергетические профили пути реакции. В качестве координаты реакции выбрали расстояние от цинка до центрального атома азота реагента, а внутренние координаты всех других атомов оптимизировали. Было найдено, что только с протонированной

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

(показаны некоторые расстояния в ангстремах)

по центральному атому азота формой реагента реакция проходит эффективно и безактивационно. Результатом реакции является отрыв группы SCH3 от цинка и присоединение ее к 1,2-диазенендикарбоксамиду. В нижней части рис. 18 показана найденная геометрическая конфигурация продуктов реакции. Видно, что структура цинксодержащего фрагмента претерпела существенные изменения - вместо тетраэдрической она стала плоской. Таким образом, результаты молекулярного моделирования подтверждают перспективность 1,2-диазенендикарбоксамида как эффективного реагента, разрушающего структуру белка КСр7. На следующих этапах результаты этих расчетов достаточно низкого уровня проверялись в более высоких

квантово-химических приближениях.

В работах [100, 101] квантово-химическими методами проведено моделирование стадии ацилиро-вания реакции гидролиза пептидной связи, катализируемой сериновыми протеазами. Катализ осуществляется боковой цепью серина, содержащей OH-группу. Сама гидроксильная группа не является сильным нуклеофилом, поэтому для протекания реакции необходим перенос протона с гидроксиль-ной группы на основание, роль которого в серино-вых протеазах выполняет гистидин (His). Также известно, что аспарагиновая кислота (Asp), находящаяся в непосредственной близости от остатка ги-стидина, принимает активное участие в процессе.

Рис. 20. Равновесные геометрические конфигурации модельного активного центра тетраэдрического интермедиата

(показаны некоторые расстояния в ангстремах)

Построение энергетической поверхности реакции для модельной системы и поиск стационарных точек (локальных минимумов и переходных состояний) проводили методом ЯНБ с применением псевдопотенциалов 8БК. Выбирали координату реакции как расстояние между кислородным атомом серина и карбонильным углеродом субстрата, при движении вдоль которой проводили оптимизацию всех остальных геометрических параметров. В локализованной области стационарной точки для выбранного значения координаты реакции проводили уточнение структуры при полной оптимизации геометрических параметров. На рис. 19 приведена рассчитанная равновесная геометрическая

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

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

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

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

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

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

Первые расчеты по программе PC GAMESS на персональном компьютере PC-AT 80486 датируются 1993 г., и за прошедшее десятилетие изменения в возможностях молекулярного моделирования поистине впечатляют. Если в первых работах речь могла идти о результатах для малоатомных молекул преимущественно качественной значимости, достижимых при расчетах методами HF, MP2, CASCSF с базисами типа 3-21G общей численностью до сотни функций, то сейчас ставятся задачи воспроизведения спектров малоатомных молекул с экспериментальной точностью. При этом используется весь арсенал квантово-химических методов с базисами достаточно больших размерностей и в задачах поиска стационарных точек на потенциальных поверхностях, а также при расчетах энергии в точке.

Перспективы дальнейшего развития молекулярного моделирования с использованием, в частности, PC GAMESS относятся к разработкам методов расчета свойств больших молекулярных систем. По-видимому, наиболее интересные пути связаны с реализацией комбинированных методов квантовой и молекулярной механики (КМ/ММ) [102]. Первые шаги в этом направлении представляются достаточно убедительными [103-106].

При написании данного обзора использованы работы, поддержанные несколькими грантами РФФИ (в частности, 01-03-32071, 02-03-39014, 02-07-90111). Мы признательны корпорации Intel за поддержку проекта развития программы PC GAMESS в 1999-2002 гг.

СПИСОК ЛИТЕРАТУРЫ

1. Schmidt M.W., Baldridge K.K., Boatz J.A., Elbert S.T., Gordon

M.S., Jensen J.H., Koseki S., Matsunaga N., Nguyen K.A., Su S.J., Windus T.L., Dupuis M., Montgomery J.A. // J. Comput. Chem. 1993. 14. P. 1347.

2. http://www.msg.ameslab.gov/GAMESS/GAMESS.html

3. http://www.dl.ac.uk/TCSC/QuantumChem/Codes/GAMESS-UK/

4. http://classic.chem.msu.su/gran/gamess/index.html

5. Tully J. // Modern Theoretical Chemistry, Semiempirical Methods

of Electronic Structure Calculations / Еd. G.A. Segal (Plenum, N.Y., 1977), V.7A, Chap. 6.

6. Nemukhin A.V., Grigorenko B.L. // Chem. Phys. Lett. 1997. 276.

P. 171.

7. Grigorenko B.L., Nemukhin A. V., Apkarian V.A. // Chem. Phys.

1998. 232. P. 321.

8. Roos B.O., Taylor P.R., Siegbahn P.E.M. // Chem. Phys. 1980. 48.

P. 157.

9. Siegbahn P.E.M.,Almlmf J.,Heiberg A., Roos B.O. // J. Chem.

Phys. 1981. 74. P. 2384.

10. Almlof J., Nemukhin A. V., Heiberg A. // Int. J. Quantum Chem.

1981. 20. P. 655.

11. Partridge H., Schwenke D.W. // J. Chem.Phys. 1997. 106. P. 4618.

12. GrigorenkoB.L., NemukhinA.V., TopolI.A., BurtS.K. //J. Chem. Phys. 2000. 113. P. 2638.

13. Nakano H. // J. Chem. Phys. 1993. 99. P. 7983.

14. Grigorenko B.L., Nemukhin A.V., Apkarian V.A. // Chem. Phys. 1998. 108. P. 4413.

15. Nemukhin A.V., Grigorenko B.L., Granovsky A.A. //Chem.Phys. Lett. 1999. 301. P. 287.

16. Grigorenko B.L., Khriachtchev L., Nemukhin A.V., Pettersson M, Isoniemi E, Rdsdnen M. // J. Chem. Phys. 1999. 110. P. 5836.

17. Firsov D.A., Grigorenko B.L., Nemukhin A.V., Khriachtchev L., Räsänen M. // Chem. Phys. Lett. 2001. 338. P. 317.

18. Pierloot K., Dumez B., Widmark P.-O, Roos B.O. //Theoret. Chim. Acta. 1995. 90. P. 87.

19. Bruna P.J., Hirsch G. // Molec. Phys. 1987. 61. P. 1359.

20. Hirst D.M., Guest M.F. // Mol. Phys. 1982. 46. P. 427.

21. Vivie R., Marian C.M., Peyerimhoff S.D. // Molec. Phys. 1988. 63. P. 3.

22. Yamanishi M., Hirao K., Yamashita K. // J. Chem. Phys. 1998. 108. P. 1514.

23. Hurley M.M., Pacios L.F., Christiansen P.A., Ross R.B., Ermler W.C. // J. Chem. Phys. 1986. 84. P. 6840.

24. Nemukhin A. V., Ermilov A.Yu., Petrukhina M.A., Klocucher W.E., Smets J. // Spectrochimica Acta PartA. 1997. 53. P. 1803.

25. Klotzbächer W.E., Petrukhina M.A., Nemukhin A.V., Ermilov A.Yu., Grigorenko B.L. // Spectrochimica Acta Part A. 2001. 57. P. 1093.

26. DolgM, Stoll H, Savin A., Preuss H. // Theor. Chim. Acta. 1989.

75. P. 173.

27. Stevens W.L., Basch H, Krauss M. //J.Chem.Phys. 1984. 81. P. 6026.

28. Немухин А.В. // ЖНХ. 1994. 39. С. 1714.

29. Schleyer P.v.R., Wurthwein E.-U., Pople J.A. // J. Am. Chem. Soc. 1982. 104. P. 5839.

30. Клименко Н.М., Крамаренко С.С., Чаркин О.П. //Координац. химия. 1984. 10. С. 3250.

31. Reed A.E., Curtiss L.A., Weinhold F. // Chem. Rev. 1988. 88. P. 89.

32. Немухин А.В., Вейнхольд Ф. // Рос. хим. журнал. 1994. 38. С. 5.

33. Nemukhin A.V. // J. Molec. Struct. (Theochem). 1994. 315. P. 225.

34. Nemukhin A.V., Weinhold F. // J. Chem. Phys. 1992. 97. P. 3420.

35. Nemukhin A.V., Almlöf J. // J. Molec. Struct. (Theochem). 1992. 253. P. 101.

36. Немухин А.В., Серебренников Л.В. //Усп. хим. 1993. 62. С. 566.

37. Imuzu Y, Klabunde K.J. // Inorg. Chem. 1984. 23. P. 3602.

38. Klabunde K.J., Whetten A. // J. Am. Chem. Soc. 1986. 108. P. 6529.

39. Nemukhin A.V., Topol I.A., Weinhold F. // Inorg. Chem. 1995. 34. P. 2980.

40. Nemukhin A.V., Solov'ev V.N., Sergeev G.B., Topol I.A. // Mendeleev Comm. 1996. №1. С. 5.

41. Solov 'ev V.N., Sergeev G.B., Nemukhin A. V., Burt S.K., Topol I.A. // J. Phys. Chem. A. 1997. 101. P. 8625.

42.Сергеев Г.Б., Смирнов В.В., Загорский В.В., Косолапое А.М. // ДАН СССР. 1981. 256. С. 1169.

43. Грановский А.А., Ванюшин А.В., Поликарпов Е.В., Ковба

B.М., Немухин А.В. // Вестн. Моск. ун-та. Сер. 2. Химия. 2001. 42. С. 371.

44. Молекулярные постоянные неорганических соединений. Справочник / Под ред. К.С. Краснова. Л., 1979.

45. Solov'ev V.N., Polikarpov E.V., Nemukhin A.V., Sergeev G.B. // J. Phys. Chem. A. 1999. 103. P. 6721.

46. Демьянов П.И., Мышакин Е.М., Крутько Д.П., Петро-сян В.С. // Изв. АН, серия химич. 1997. №12. С. 2111.

47. Хайкин Л.С., Грикина О.Е., Сипачев В.А., Грановский А.А., Никитин В.С. // Известия АН, серия химич. 2000. №4.

C. 616.

48. Соломоник В.Г., Слизнев В.В., Балабанов Н.Б. //ЖНХ. 1997. 42. С. 496.

49. Соломоник В.Г., Марочко О.Ю. // Журн. структурн. химии. 2000. 41. С. 885.

50. Соломоник В.Г., Марочко О.Ю. // ЖФХ. 2000. 74. С. 94.

51. Соломоник В.Г., Маренич А.В.// ЖФХ. 2000. 74. С. 94.

52. Маренич А.В., Соломоник В.Г. // ЖНХ. 2001. 46. С. 462.

53. Bataev V.A., Dolgov E.K., Abramenkov A.V., Pupyshev V.I., Godunov I.A. // Intern. J. Quant. Chem. 2002. 88. P. 414.

54. Годунов И.А., Батаев В.А., Яковлев Н.Н., Пупышев В.И. // ЖФХ. 2002. 76. С. 1912.

55. Долгов Е.К., Батаев В.А., Годунов И.А. // ЖФХ. 2000. 74. С. 94.

56. Немухин А.В. // Рос. хим. журн. 1996. 40. С. 48.

57. Nemukhin A.V., Grigorenko B.L. //Chem. Phys. Lett. 1995. 233. P. 627.

58. Nemukhin A.V., Grigorenko B.L. // Intern. J. Quant. Chem. 1997. 62. P. 55.

59. Klabunde K.J. Free Atoms, Clusters and Nanoscale Particles. N.Y., 1994.

60. Немухин А.В. // Изв. АН, серия физич. 1997. 61. №9. С. 1784.

61. Ермилов А.Ю., Немухин А.В., Сергеев Г.Б. // Изв. АН, серия физич. 1998. 62. №6. С. 1169.

62. Московский А.А., Немухин А.В. // Вестн. Моск. ун-та. Сер. 2. Химия. 1998. 39. С. 83.

63. Ermilov A.Yu., Nemukhin A.V., Kovba V.M. // Mendeleev Comm.

1999. №3. P. 88.

64. Kazankov G.M., Dyachenko O.G., Nemukhin A.V., Ryabov A.D. // Mendeleev Comm. 1997. №12. P. 159.

65. Cundary T.R., Stevens W.J. // J. Chem. Phys. 1993. 98. P. 555.

66. Юхневич Г.В., Тараканова Е.Г., Немухин А.В. //Изв. АН, серия химич. 1997. №3. С. 435.

67. Grigorenko B.L., Moskovsky A.A., Nemukhin A.V. // J. Chem. Phys. 1999. 111. P. 4442.

68. Grigorenko B.L., Moskovsky A.A., Nemukhin A.V. // J. Molec. Struct. (Theochem). 2000. 498. P. 47.

69. Богдан Т.В., Грановский А.А., Немухин А.В. // Вестн. Моск. ун-та. Сер. 2. Химия. 2000. 41. С. 98.

70. Moskovsky A.A., Grigorenko B.L., Granovsky A.A., Nemukhin A.V. // Russ. J. Phys. Chem. 2000. 74. P. 228.

71. Moskovsky A.A., Nemukhin A.V. // J. Chem. Inf. Comput. Sci. 1999. 39. P. 370.

72. Шабатина А.В., Немухин А.В. // Вестн. Моск. ун-та. Сер. 2. Химия. 2000. 41. С. 160.

73. Nemukhin A.V., Granovsky A.A., Firsov D.A. // Mendeleev Comm. 1999. №6. P. 215.

74. Немухин А.В., Тогонидзе В.В., Ковба В.М., Орлов Р.Ю. // Журн. структурн. химии. 1998. 39. С. 460.

75. Шабатина А.В., Немухин А.В. // Вестн. Моск. ун-та. Сер. 2. Химия. 1999. 40. С. 147.

76. Новаковская Ю.В., Степанов Н.Ф. // ЖФХ. 2000. 74. С. 71.

77. Ozhegova N.V., Nemukhin A.V., Shabatina T.I., Sergeev G.B. // Mendeleev Comm. 1998. №6. P. 218.

78. Поликарпов Е.В., Шабатина Т.И., Сергеев Г.Б., Немухин А.В. // Вестн. Моск. ун-та. Сер. 2. Химия. 2000. 41. С. 283.

79. Shabatina T.I., Vovk E.V., Ozhegova N.V., Morosov Y.N., Nemukhin A.V., Sergeev G.B. // Materials Science and Engineering C. 1999. 8-9. P. 53.

80. Sergeev G.B., Shabatina T.I., Solov'ev V.N., Nemukhin A.V. // Spectrochimica Acta, Part A. 2000. 56. P. 2527.

81. Шайтан К.В., Немухин А.В., Фирсов Д.А., Богдан Т.В., Тополь И.А. // Молек. биол. 1997. 31. С. 108.

82. Немухин А.В. // ЖФХ. 1995. 69. С. 889.

83. Nemukhin A.V., Weinhold F. // J. Chem. Phys. 1992. 97. P. 1095.

84. Nemukhin A.V., Weinhold F. // J. Chem. Phys. 1993. 98. P. 1329.

85. Немухин А.В., Колесников И.М., Винокуров В.А. // ЖФХ. 1996. 70. С. 285.

86. Nemukhin A.V., Kolesnikov I.M., Vinokurov V.A. // Mendeleev Comm. 1994. №4. P. 144.

87. Beregovaya I.V., Schchegoleva L.N. // Chem. Phys. Lett. 2001. 348. P. 501.

88. Beregovaya I.V., Schchegoleva L.N. // Intern. J. Quant. Chem. 2002. 88. P. 481.

89. Polikarpov E.V., Granovsky A.A., Nemukhin A.V. //Mendeleev Comm. 2001. №4. P. 150.

90. Nemukhin A.V., Grigorenko B.L., Khriachtchev L., Tanska-nen H., Pettersson M., Rдsдnen M. // J. Am. Chem. Soc. 2002. 124. P. 10706.

91. Khriachtchev L., Pettersson M., Runeberg N., Lundell J., Rдsдnen M. // Nature (London). 2000. 406. P. 874.

92. Pettersson M., Khriachtchev L., Lundell J., Rдsдnen M. //J. Am. Chem. Soc. 1999. 121. P. 11904.

93. Nagel J.F., MorowitzH. // Proc. Natl. Acad. Sci. USA. 1978. 75. P. 298.

94. Григоренко Б.Л., Князева М.А., Немухин А.В. // Вестн. Моск. ун-та. Сер. 2. Химия. 2001. 42. С. 92.

95. Nemukhin A.V., Topol I.A., Grigorenko B.L., Burt S.K. //J. Phys. Chem. B. 2002. 106. P. 1734.

96. Конюхов С.С., Московский А.А., Фирсов Д.А., Немухин А.В. // Вестн. Моск. ун-та. Сер. 2. Химия. 2002. 43. С. 4.

97. Gordon M.S., Freitag M.A., Bandyopadhyay P., Jensen J.H., Kairys V., Stevens W.J. // J. Phys. Chem. A. 2001. 105. P. 293.

98. Nemukhin A.V., Grigorenko B.L., Skasyrskaya E.Ya., Topol I.A., BurtS.K. // J. Chem. Phys. 2000. 112. P. 513.

99. Topol I.A., Nemukhin A.V., Dobrogorskaya Ya.I., Burt S.K. //J. Phys. Chem. B. 2001. 105. P. 11341.

100. Nemukhin A.V., Topol I.A., Burt S.K. // Intern. J. Quant. Chem. 2002. 88. P. 34.

101. Доброгорская Я.И., Немухин А.В. // Вестн. Моск. ун-та. Сер. 2. Химия. 2003. 44. С. 103.

102. Warshel A., LevittM. // J. Mol. Biol. 1976. 103. P. 227.

103. Nemukhin A.V. Grigorenko B.L.,Bochenkova A.V., Topol I.A., Burt S.K. // J. Molec. Struct (Theochem). 2002. 581. P. 167.

104. Grigorenko B.L., Nemukhin A.V.,.Topol I.A., Burt S.K. // J. Phys. Chem. A. 2002. 106. P. 10663.

105. Nemukhin A.V. Grigorenko B.L.,Topol I.A., Burt S.K. // J. Phys. Chem. B. 2003. 107. P. 2958.

106. Nemukhin A.V., Grigorenko B.L., Topol I.A., Burt S,K. //J. Comp. Chem. 2003. 24. P. 1410.

Поступила в редакцию 03.01.04.

MOLECULAR MODELING BY USING THE PC GAMESS PROGRAM: FROM DIATOMIC MOLECULES TO ENZYMES

A.V.Nemukhin, B.L.Grigorenko, A.A.Granovsky

(Division of Physical Chemistry)

Review of the main results in molecular modeling in various fields of Chemistry, obtained by using the computer program PC GAMESS, which is oriented on widespread personal computers and workstations.