Метод+базис для растворителя

[tabaki321]

Привет всем.
Есть необходимость рассчитывать ионы (для начала бензоат) в воде. На данном этапе рассматривается применение для этого ONIOM (NWChem-6.1.1). Кластер 1*бензоат ион + около 100 молекул воды. Для аниона и нескольких молекул воды у карбоксильной группы B3LYP/6-311++g(2d,2p). Вопрос в том, какой метод+базис применить при расчете следующих слоев воды? Достаточно ли, что-то вроде B3LYP/6-31G(d) или надо подключать диф. функции?
Всем замечаниям и предложениям буду рад.

[VTur]

Какое свойство Вы хотите получить?

[ Post made via Android ]

[sanya1024]

Давайте уточним: кластер из бензоата+100 "живых" H2O погружен в ONIOM-ную воду? или эти 100 молекул предполагается распределить по степени "продвинутости" учета?
И действительно, какие свойства моделируете: надеюсь, не возбуждения с переносом заряда на растворитель? Если нет, то вполне достаточно внутреннюю координационную сферу воды сделать "живой" с полноценным базисом (таким же, как для самого аниона), а более далекие -- либо в виде ММ, либо в виде какой-то жесткой поляризуемой модели (что для этого есть в NWChem, я не знаю, а в GAMESS-US и Q-Chem есть эффективные фрагментные потенциалы, EFP).

[tabaki321]

Спасибо за то что откликнулись.
Планируется сравнить теплоты реакций с участием бензоат (реакция 1) и 4-нитробензоат (реакция 2) ионов. Экспериментально измеренная разница составляет порядка 30кДж/моль
Второй исходник бензоил хлорид, а продукты - соответствующие ангидриды. Оба находятся в органической фазе (дихлорметан).
Уточняю структуру. 1 бензоат ион + 6 молекул воды (расположенных около карбоксильных кислородов) - это планируемый "high level" B3LYP/6-311++g(2d,2p). Остальная вода (порядка еще 100 молекул) - это "low level" B3LYP/6-31G(d) (здесь вопрос, можно ли еще опустить уровень теории, например до RHF/3-21g).
Причины использования явного учета растворителя: специфическая сольватация аниона и нежелание использовать COSMO.
Про EFP1 читал, но он только для воды. EFP2, как я понял, не дает возможности проводить оптимизацию ab initio фрагмента. А в дальнейшем планируется расчет скоростей реакций, протекающих в растворах (вода, дихлорметан, хлороформ, тетрахлорметан, этанол). Поэтому хотелось бы иметь опробованный алгоритм явного учета растворителя.
Всегда рад предложениям и замечаниям.

[sanya1024]

Это хорошо, что Вы про EFP читали EFP2 (в варианте fragments only) позволяет гонять какое-то подобие мол.динамики. Но для этого придется и бензоат тоже превратить в EFP2. Ну, это так если вдруг на будущее решите моделировать неводные или смешанные растворители. А к тому времени, глядишь, и сделают оптимизацию ab initio в EFP2 Я бы скорее в этом направлении работала.
Дело в том, что 100 молекул воды даже в 6-31G базисе -- это еще куда ни шло, а 100 молекул хлороформа -- это уже ой! С другой стороны, вода тянет сетку водородных связей довольно далеко, а в неводных растворителях 1-2 координационные сферы -- и растворитель не "помнит", что где-то там сидит молекула (даже если это ион). Так что для неводных можно кластеры поменьше брать и комбинировать с континуальным растворителем. А утверждение насчет 1-2 координационных сфер можно как раз проверить с помощью мол.динамики фрагментов, построив радиальную функцию распределения. Кстати, можно и для воды посмотреть, на каком расстоянии она "забывает" о растворенной молекуле, и на этом основании выбрать размер кластера для моделирования.
Отдельный вопрос, имеет ли смысл "дальние" молекулы воды моделировать на уровне RHF/3-21g. Мне кажется, любые параметризованные дискретные модели воды, хоть даже мол. механические (TIP5P, SPC) будут лучше описывать водородные связи, по крайней мере связи воды с водой. Тогда можно в ab initio ограничиться кластером из 1 сольватной оболочки.
И еще: в неводных растворителях особенно, но и в воде тоже не забывайте про дисперсионную поправку Гримме к функционалу. В неводных она дает бОльшую часть межмол. взаимодействия, но и в воде тоже влияет.

[tabaki321]

Спасибо за совет. На выходных обязательно изучу вопрос.

[VTur]

Планируется сравнить теплоты реакций с участием бензоат (реакция 1) и 4-нитробензоат (реакция 2) ионов. Экспериментально измеренная разница составляет порядка 30кДж/моль.

Вот Вы получите полную энергию Вашей системы. Одну единую цифру для молекулы и растворителя. (точнее 2 цифры для реакции - до и после) Как Вы собираетесь превратить ее в теплоту реакции индивидуального вещества, да еще с учетом растворителя?

[tabaki321]

Предположительно следующим образом:
В ONIOM энергия системы: E = E(низкий уровень, вся система) - Е(низ. ур., внутрю часть) + Е(выс. ур., внутр. часть)
В ходе оптимизации снижается конечная E. А для сравнения двух реакций будут браться Е(выс. ур., внутр. часть).
Система в целом: анион (выс. ур., внутр. часть), 1 сольватационный слой (выс. ур., промеж. часть) и остальной растворитель (низкий уровень, вся система).
Как я (от Sanya1024) понял низким уровнем при системе таких размеров (суммарно около 120-150 молекул, особенно если растворитель уже не вода) целесообразно делать ММ (полуэмпирика в NWChem отсутствует). Отсюда вопрос про силовые поля для органики (хлороформ, дихлорметан, растворенные в них пиридин и его пара замещенные производные)

[sanya1024]

Тут, к сожалению, ничего хорошего не подскажу. Мои коллеги, кто непосредственно этим занимается, тоже в поиске. Обычно выбор -- из тех полей, с к-рыми работает программа. С чем умеет работать NWChem, я не знаю. Но Амбер, с к-рым работают почти все программы, точно не годится -- в нем нет параметризации для хлорорганики. Может, обратить внимание на поля семейства MM2 (Аллинджер) и их современные модификации?

[tabaki321]

Спасибо. В NWChem похоже заложено что-то из Amber. Спрошу у них на форуме про другие силовые поля. Возможно придется уменьшать кластер и задействовать COSMO. Буду изучать вопрос дальше.

[sanya1024]

Это не худший вариант. Правда, COSMO -- это COnductor-like Screening MOdel, а органические растворители -- не очень-то conductor-like. Во всяком случае, обычный C-PCM (отличающийся от COSMO одним коэффициентиком в формуле) в GAMESS и Гауссиане давал неважные результаты, D-PCM (Dielectric PCM) был более адекватным. Вроде бы миннесотская модель SMD (Solvation Model Density), недавно появившаяся в тех же GAMESS и Гауссиане, тоже получше обычного C-PCM будет. Но все от системы зависит, надо разные варианты пробовать.

[VTur]

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

[ Post made via Android ]

[tabaki321]

В свете сказанного VTur, боюсь сейчас буду задавать глупые вопросы. Но для начала уточню схему интересующей реакции
1) Есть двухфазная система вода/дихлорметан
2) В органике сидит бензоил хлорид, в воде - бензоат (4-нитробензоат) ион и катализатор (пиридин-1-оксид или его производное).
3) Продукт - соответствующий ангидрид образуется в водной фазе, но сразу переходит в органику.
Есть экспериментальные величины теплот этих реакций. Задача: сопоставить эксп. с теор.
Если я правильно все понимаю, то с термодинамически интересны только начальные и конечные состояния системы. Получается такое уравнение:
BzСl[орг] + BzO-[водн] = BzOBz[орг] + Cl-[водн]
Отсюда:
dH[реакции] = dH[BzOBz_орг] + dH[Cl_водн] - dH[BzO-_орг] - dH[BzO-_водн]
ddH[реак_BzO - реак_4-O2NBzO] = dH[BzOBz_орг] - dH[BzO-_водн] - (dH[4-O2NBzOBz_орг] - dH[4-O2NBzO-_водн])
Таким образом, мне надо бы получить энергии отдельных молекул в растворителе.
Их геометрии, полные энергии и термодинамические поправки уже были получены для газовой фазы B3LYP/6-311++G(2d,2p). ddH составило, насколько я помню, около 70 кДж. Экспериментальное значение ddH около 30 кДж, т.е. разница эксп. и теор. порядка 10 ккал. dH[реакции] эксп. и теор. отличаются уже очень сильно (около 100 кДж). Но для расчетов проведенных для газовой фазы чего-то другого было ожидать глупо. Вот и встал вопрос о расчетах с учетом растворителя.

[sanya1024]

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

[tabaki321]

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

[sanya1024]

Химическая интуиция А собс-но, почему структура достаточно жесткой молекулы в растворе (изотропном, что характерно) должна принципиально отличаться от газофазной? Я имею в виду: то, что в газовой фазе было минимумом, вряд ли станет седлом в изотропном континууме. Конечно, "газофазная" геометрия уже не будет соответствовать минимуму на новой, "растворной" ППЭ, но уйдет оттуда недалеко -- сдвинется на стенку потенциальной ямы, и несколько итераций оптимизации геометрии приведут ее в новый, почти такой же минимум.

[VTur]

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

[ Post made via Android ]

[VTur]

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

[ Post made via Android ]

[tabaki321]

VTur, спасибо за объяснения. Значит буду пробовать различные эмпирич. модельки (PCM и др.). Этот пример не для диссертации, а для годового отчета коллеги. Сам я пока занимался только оптимизацией в газе, расчетом колебательных частот в гарм. приближении (самый простой вариант) и прочими несложными задачками.
Если получится получить адекватные результаты для систем с бензоатами, тогда я буду пробовать на своих. Причина в том, что они значительно тяжелее ("цветочки": хлордифенилфосфат, триарилфосфаты, тетрафенилпирофосфат и "ягодки": пиридиниевые соли хлорфосфатов и ТС с их участием). Все это дополняется наличием нескольких ротамеров для каждого из объектов.
Понимаю, что во многих случаях "прокатывает" недобросовестное отношение к расчетам. Я даже уверен, что если проведу все расчеты в газовой фазе, то вряд ли кто-то в нашем местном совете будет резко против (максимум сделают рядовое замечание). Но так как хотелось бы сделать работу не для галочки, то буду дальше изучать вопрос про модельки растворителей. Пока наиболее реальным представляется: использование continuum моделей с добавленными при необходимости нескольких молекул растворителя.

[sanya1024]

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