Расчет и оптимизация структуры нафталена и азулена с помощью MOPAC
Так же, как и в предыдущем практикуме, надо создать пространственные структуры нафталена и азулена, используя SMILES-аннотации. Затем проводим оптимизацию с помощью MOPAC:
export PATH=${PATH}:/home/preps/golovin/progs/bin
export MOPAC_LICENSE=/home/preps/golovin/progs/bin
echo "C1=CC=C2C=CC=C2C=C1 azulene" > az.smi
obgen az.smi > az.mol
babel -imol az.mol -omop az.mop -xk "PM6"
MOPAC2009.exe az.mop
babel -imopout az.out -opdb az.pdb
echo "c1ccc2ccccc2c1 naphthalene" > napht.smi
obgen napht.smi > napht.mol
babel -imol napht.mol -omop napht.mop -xk "PM6"
MOPAC2009.exe napht.mop
babel -imopout napht.out -opdb napht.pdb
Ниже представлены изображения полученных структур (Рисунки 1 и 2)
Рисунок 1. Структура азулена, полученная из SMILES-аннотации и оптимизированная с помощью MOPAC.
Рисунок 2. Структура нафталена, полученная из SMILES-аннотации и оптимизированная с помощью MOPAC.
Полученная структура нафталена лежит в одной плоскости (как и положено ароматической молекуле). Однако, ароматическая структура азулена оказалась не плоской, поэтому продолжаем оптиизацию (на этот раз попробуем использовать силовое поле MMFF94 в obgen):
obgen az.smi -ff MMFF94 > az_mmff94.mol babel -imol az_mmff94.mol -omop az_mmff94.mop -xk "PM6" MOPAC2009.exe az_mmff94.mop babel -imopout az_mmff94.out -opdb az_mmff94.pdb
Изображение структуры азулена, полученной с использованием силового поля MMFF94 представлено на Рисунке 3.
Рисунок 3. Структура азулена, полученная с использованием силового поля MMFF94.
Полученная структура также не плоская. Поэтому пробуем еще раз получить данную структуру, используя силовое поле MMFF94s:
obgen az.smi -ff MMFF94s > az_mmff94s.mol babel -imol az_mmff94s.mol -omop az_mmff94s.mop -xk "PM6" MOPAC2009.exe az_mmff94s.mop babel -imopout az_mmff94s.out -opdb az_mmff94s.pdb
Полученную структуру азулена с использованием силового поля MMFF94s можно увидеть на рисунке 4.
Рисунок 4. Структура азулена, полученная с использованием силового поля MMFF94s.
И опять структура неплоская. Пробуем использовать силовое поле UFF:
obgen az.smi -ff UFF > az_uff.mol babel -imol az_uff.mol -omop az_uff.mop -xk "PM6" MOPAC2009.exe az_uff.mop babel -imopout az_uff.out -opdb az_uff.pdb
Изображение структуры - на Рисунке 5.
Рисунок 5. Структура азулена, полученная с использованием силового поля UFF.
Ура, плоская структура азулена получена!
Оптимизация геометрии средствами GAMESS
В предыдущем задании было получено 2 файла: azu.out и nap.out (azu.out был переименован из azu_uff.out, nap.out - из napht.out). Сначала с помощью babel переформатируем координаты в gamin формат:
babel -imopout azu.out -ogamin azu_opt.inp babel -imopout nap.out -ogamin nap_opt.inp
Теперь поменяем заголовок в полученных файлах, чтобы он выглядел так:
$CONTRL COORD=CART UNITS=ANGS SCFTYP=RHF RUNTYP=OPTIMIZE $END $BASIS GBASIS=N31 NGAUSS=6 $end $system mwords=2 $end $DATA
Таким образом получили входные файлы для оптимизации геометрии средствами GAMESS. Проводим оптимизацию для обеих молекул и получаем .pdb-файл с их структурой:
gms azu_opt.inp 1 >& azu_opt.log babel -igamout azu_opt.log -opdb azu_gamess.pdb gms nap_opt.inp 1 >& nap_opt.log babel -igamout nap_opt.log -opdb nap_gamess.pdb
В результате оптимизации в структуре азулена изменилась длина двух связей с 1.4 до 1.5 (рисунок 6). В структуре нафталена в результате оптимизации изменений не произошло (рисунок 7).
Рисунок 6. Структура азулена до оптимизации GAMESS (слева) и после (справа).
Рисунок 7. Структура нафталена до оптимизации GAMESS (слева) и после (справа).
Расчет энергии методом Хартри-Фока и по теории функционала плотности (DFT)
На основе полученных координат составляем новые входные файлы для расчёта энергии. Сначала переформатируем .log файл gamout в gamin для исползования babel:
babel -igamout azu_opt.log -ogamin azu_gamin.inp babel -igamout nap_opt.log -ogamin nap_gamin.inp
Tеперь получаем файлы azuhf.inp и naphf.inp для расчета по Хартри-Фоку и файлы azudf.inp и napdf.inp для расчета по теории функционала плотности. Файлы для расчета методом Хартри-Фока содержат следующий заголовок:
$CONTRL COORD=CART UNITS=ANGS SCFTYP=RHF RUNTYP=ENERGY $END $BASIS GBASIS=N31 NGAUSS=6 POLAR=POPN31 NDFUNC=1 $END $GUESS GUESS=HUCKEL $END $system mwords=2 $end $DATA
Файлы для расчета по теории функционала плотности содержат следующий заголовок:
$CONTRL COORD=CART UNITS=ANGS dfttyp=b3lyp RUNTYP=ENERGY $END $BASIS GBASIS=N31 NGAUSS=6 POLAR=POPN31 NDFUNC=1 $END $GUESS GUESS=HUCKEL $END $system mwords=2 $end $DATA
Запустим GAMESS для расчета энергии:
gms azuhf.inp 1 >& azuhf.log gms azudf.inp 1 >& azudf.log gms naphf.inp 1 >& naphf.log gms napdf.inp 1 >& napdf.log
Сравнение методов
| Вещество | ENaphthalene | EAzulene | Δ E, Hartree* | Δ E, kCal/mol* |
|---|---|---|---|---|
| Хартри-Фок | -383.3523603875 | -383.2825020374 | 0,06985835 | 43,83681327 |
| DFT | -385.6417198702 | -385.5857377223 | 0,055982148 | 35,12935763 |
Из эксперимента известно, что энергия изомеризации нафталина в азулен составляет 35.3±2.2 kCal/mol, что численно ближе к результатам подсчета по DFT.