Файл: Методическое пособие Расчет возбужденных электронных состояний многоконфигурационными методами квантовой химии.pdf

16
ORBITALS». Конфигурации детерминантов представлены в виде нулей (орбиталь свободна) и единиц (орбиталь занята). Слова ALPHA и BETA обозначают спин- орбитали с противоположным спином. Справа указан коэффициент разложения волновой функции электронного состояния по соответствующему детерминанту, который показывает с каким весом он входит в волновую функцию. Детерминант, коэффициент которого в разложении является наибольшим, является ведущим детерминантом.
После расчета методом CASSCF необходимо проанализировать орбитали активного пространства. Выбранные на предыдущем этапе в активное пространство молекулярные орбитали могли сильно измениться в ходе CASSCF расчета и перестать быть π орбиталями. Если это произошло необходимо скорректировать активное пространство и перезапустить расчет заново используя в качестве начального приближения только что получившиеся натуральные орбитали.
Известно, что процедура CASSCF быстрее сходится, если в качестве начального приближения использовать натуральные орбитали
(полученные из оптимизированных орбиталей путем ортогонального преобразования и диагонализующие матрицу плотности), а не орбитали полученные в рамках метода
Хартри-Фока. Поэтому оптимальная схема действий, следующая: сначала запускаются несколько расчетов с малым числом итераций (3-5) до тех пор, пока орбитали в активном пространстве не перестанут менять свой характер и состав активного пространства будет соответствовать первоначальным ожиданиям.
Наличие изменений в орбиталях часто возможно контролировать визуально, а также можно следить за строкой «SOSCF IS SCALING ROTATION ANGLE MATRIX» в out файле сразу перед секцией «-MCCI- BASED ON OPTIMIZED ORBITALS». Eе наличие означает, что происходит преобразование исходных молекулярных орбиталей в ходе расчета путем поворотов, причем чем больше эти повороты, тем большие изменения претерпевают орбитали. После того, как орбитали перестанут сильно меняться (повороты будут небольшим или будут отсутствовать вовсе), можно запустить окончательный расчет с большим числом итераций (например,
100). Не забудьте проверить, что расчет завершился без ошибок! Убедиться в этом можно по наличию строки «EXECUTION OF FIREFLY TERMINATED
NORMALLY», а также «DENSITY CONVERGED» или «LAGRANGIAN
CONVERGED». Аналогичный порядок действий выполните для второго изомера.
После успешного выполнения этого этапа у Вас будут орбитали CASSCF для каждого из изомеров. Орбитали находятся в файле *.dat.
Этап 4. Расчет возбужденных состояний методом MRMP2.
После получения молекулярных орбиталей методом CASSCF можно приступить к расчёту энергий вертикальных переходов бензохинонов. Энергии вертикальных переходов будут вычисляться комбинированным методом MRMP2. В качестве примера можно использовать файл bqp_mrmp2.inp, его структура приведена в

17 приложении 5. Необходимо запустить два расчета: один для основного состояния, другой для первого возбужденного. В первой случае в карточке $xmcqdpt надо указать nstate=1 wstate(1)=1 avecoe(1)=1, а во втором nstate=2 state(1)=2 kstate(1)=0,1.
Вставьте полученные на третьем этапе орбитали ($VEC) в inp файл. Понять, что расчет завершился успешно, можно по фразе «EXECUTION OF FIREFLY
TERMINATED NORMALLY». После чего найти в файле фразу *** CAS-CI
ENERGIES *** в файле для основного и в файле для возбужденного состояния и найти вертикальную энергию перехода. Сравните полученный результат с методом
MCQDPT2, в котором энергии состояний получаются в рамках одного расчета
(пример файла в приложении 5). Выполните аналогичные расчеты для обоих изомеров, дополнив активное пространство n орбиталью локализованной на кислороде. Прокомментируйте полученные результаты.
Дополнительное задание. Проверьте как поменяются результаты расчетов, если в качестве орбитальных энергий для теории возмущений использовать энергии, полученные методом DFT. Для этого в карточку $contrl добавьте фразу dfttyp=pbe0, а в карточку $xmcqdpt2 фразы ifitd=2 dftoe=0x0f.
Этап 5. Расчет возбужденных состояний методом MRCI.
Пример файла для расчета находится в приложении 6. Проведите несколько расчетов с различным максимальным числом кратности возбуждения, варьируя параметр IEXCIT от 2 до 4. В файле найдите фразу THE WAVEFUNCTION
CONTAINS запишите число конфигураций в каждом случае. Чтобы узнать рассчитанное значение энергии состояния найдите фразу STATE # 1 ENERGY.
Выполните аналогичные расчеты для обоих изомеров, дополнив активное пространство n орбиталью локализованной на кислороде. Прокомментируйте полученные результаты. Сравните полученные результаты между собой и с теорией возмущений.
Дополнительное задание. Проверьте как поменяются результаты расчетов, если в качестве референсного состояния использовать не S
0
, а S
1
? Прокомментируйте полученные результаты. Для этого в карточке $CIDRT замените соответствующие команды на NFZC=24 NDOC=3 NVAL=3 NAOS=1 NBOS=1 NEXT=-1.

18
Список литературы
Основная
1. Барановский В.И. Квантовая механика и квантовая химия. Москва:
Издательство Лань, 2017. 428 с.
2. Ермаков А.И. Квантовая механика и квантовая химия. Москва: Издательство
Юрайт, 2015. 556 с.
3. Уилсон С. Электронные корреляции в молекулах. Пер. с англ. Издательство
Мир, 1987. 304 c.
Дополнительная
1. Szabo A., Ostlund N. Modern Quantum Chemistry: Introduction to Advanced
Electronic Structure Theory. New York: Dover Publications, 1989. 466 p
2.
Björn O. Roos, Roland Lindh, Per Åke Malmqvist, Valera Veryazov, Per-
Olof Widmark
. Multiconfigurational Quantum Chemistry. Wiley, 2016. 224 p.
3. Степанов Н.Ф., Пупышев В.И. Квантовая механика молекул и квантовая
химия. Москва. Издательство МГУ, 1991. 384 c.

19
Приложения
Приложение 1. Пример входного файла для оптимизации геометрии.
Входной (*.inp) файл для оптимизации геометрии пара-бензохинона методом MP2 в базисе (aug)-cc-pVDZ (с диффузными функциями, локализованными только на атомах кислорода).
$contrl icharg=0 scftyp=rhf d5=.t. runtyp=optimize nzvar=1 mult=1 maxit=100 icut=20 inttyp=hondo itol=30 gencon=.f. wide=1 mplevl=2 fstint=.t. $end
$moorth syms=1 symden=1 symvec=1 tole=0 tolz=0 nostf=1 nozero=1 $end
$mp2 method=1 $end
$system mwords=100 timlim=100000 nojac=1 splitf=.f. ioflgs(30)=0 $end
$p2p p2p=.t. dlb=.t. tracebs=.t.
$end
$smp csmtx=.t. $end
$trans altpar=1 mptran=2 dirtrf=.t. mode=112 cuttrf=1d-30 $end
$SCF DIIS=.t. SOSCF=.f. DIRSCF=.T. fdiff=.f. nconv=7 $END
$mp2grd dblbf=0 fused=0 asysnc=1 tol1=1d-12 tol2=1d-12
$end
$zmat dlc=1 auto=1 orttol=1d-8 $end
$BASIS GBASIS=avdz extfil=.t. $END
$statpt nstep=1000 opttol=1d-4 hssend=.t. $end
$FORCE NVIB=2 VIBSIZ=0.005 $END
$data para bq opt
Cs
CARBON 6.0 0.340927711 0.000000000 -1.558223531
CARBON 6.0 -0.864191680 0.000000000 -0.835714967

20
CARBON 6.0 1.546047101 0.000000000 -0.835714967
CARBON 6.0 -0.864515536 0.000000000 0.558137192
CARBON 6.0 1.546370957 0.000000000 0.558137192
CARBON 6.0 0.340927711 0.000000000 1.261784767
HYDROGEN 1.0 -1.807408195 0.000000000 -1.373122642
HYDROGEN 1.0 2.489263617 0.000000000 -1.373122642
HYDROGEN 1.0 -1.808618909 0.000000000 1.095422867
HYDROGEN 1.0 2.490474331 0.000000000 1.095422867
OXYGEN 8.0 0.340927711 0.000000000 -2.998223531
OXYGEN 8.0 0.340927711 0.000000000 2.701784767
$end
В данном входном файле в первой строке в карточке $contrl задается тип расчета – в нашем случае задается оптимизация геометрии (RUNTYP=optimize) в методе MP2
(mplevl=2). Параметр icharg устанавливает заряд молекулярной системы, в данном случае равный нулю. Параметры scftyp, mult, nzvar, maxit задают тип волновой функции, мультиплетность электронного состояния, тип координат, используемых при заданном типе расчета (в данном случае оптимизации), количество итераций процедуры самосогласования, соответственно. В данном расчете мы используем волновую функцию типа ограниченного метода Хартри-Фока, мультиплетность электронного состояния равна 1 (синглет), расчёт проходит во внутренних координатах, а для итеративной процедуры самосогласования выбрано 100 шагов.
Раздел с ключевым словом $basis позволяет задать базисный набор атомных орбиталей, используемый в расчете. Параметр extfil указывает, считывать ли базис из внешней библиотеки (файл *.lib) или базис указан внутри основного файла, при значении параметра .t. считывание происходит из внешнего файла. Раздел $statpt контролирует поиск стационарных точек, в данном случае поиск минимума энергии основного электронного состояния при оптимизации геометрии. Параметр nstep указывает максимальное количество шагов, которое будет сделано программой при оптимизации. Параметр opttol определяет значение наибольшей компоненты градиента, при котором расчет градиента останавливается, и геометрия молекулы считается оптимизированной. Наконец, команда
HSSEND=.TRUE. запускает расчет гессиана после того, как оптимизация геометрии успешно завершится. Раздел
$SYSTEM отвечает за ограничения по времени в минутах (TIMLIM=100000) и по выделяемой памяти в миллионах машинных слов (MWORDS=100). Остальные параметры разделов описывают технические детали расчета (точность расчетов интегралов и т.п). Далее идет раздел $DATA отвечающий за исходную геометрию и

21 имеющий четко определенную структуру. Следующей строкой после названия раздела идет строка комментария, который будет вставляться в заголовки расчетов.
Ее можно оставлять пустой, но пропускать нельзя. Затем идет строка, отвечающая за симметрию молекулы. Необходимо помнить, что при задании симметрии программа сама размножает атомы по заданной симметрии, поэтому эквивалентные атомы указывать несколько раз не нужно. Если же симметрия С
1
, то геометрия молекулы начинает задаваться со следующей строки без пропусков. Геометрия молекулы задается следующим образом: в каждой строке указывается атом в формате название, заряд ядра, декартовы координаты. После окончания задания геометрии без пропуска строк идет окончание раздела $end. Более подробную информацию смотрите в руководстве пользователя Firefly version 8.2.0. http://classic.chem.msu.su/gran/firefly.

Приложение 2. Пример скрипта для запуска оптимизации.
#!/bin/sh
#SBATCH -n 36 -N 1
#SBATCH --exclusive
#SBATCH --no-requeue
#SBATCH -t 3-00:00:00 ulimit -s 4096 export I_MPI_DEBUG=10 export I_MPI_PIN=0 export I_MPI_WAIT_MODE=0 export I_MPI_FALLBACK=1 export MALLOC_CHECK_=1
/opt/software/intel/impi/oldmpi/ia32/bin/mpirun -bootstrap slurm -perhost 36
/opt/software/firefly/ff9/pfirefly -prealloc:350 -rf -ff -p -stdext -V -prof -t
/scratch/antonibo/ttt."$SLURM_JOB_ID" -ex /opt/software/firefly/ff9/ -outbuf:80 -i
/home/antonibo/EFP/bqp_opt.inp -o /home/antonibo/EFP/bqp_opt.out -b
/home/antonibo/EFP/avdz.lib
-eoi
- stderr=/scratch/antonibo/xlogs."$SLURM_JOB_ID" -synchro:nodebug #&
Строка «#SBATCH -n 36 -N 1» указывает, что для расчета будет использоваться один узел, на котором расположено 36 ядер. Строка «#SBATCH -t 3-00:00:00» указывает максимальное время расчета в днях. Команда «-i» задает путь до входного файла, команда «-o» - путь до выходного файла и его название. Команда «-b» задает путь до базисного набора. Команда «-t» задает путь до папки, в которую будут записываться временные файлы. Остальные параметры относятся к техническим характеристикам расчета.

22
Приложение 3. Пример входного файла для расчета энергии методом ХФ.
$contrl icharg=0 mult=1 scftyp=rhf d5=.t. runtyp=energy maxit=100 icut=20 inttyp=hondo itol=30 gencon=.f. wide=1 fstint=.t. exetyp=run $end
$moorth syms=1 symden=1 symvec=1 tole=0 tolz=0 nostf=1 nozero=1 $end
$system mwords=360 timlim=100000 nojac=1 $end
$system mxbcst=-1 $end
$p2p p2p=.t. dlb=.t. $end
$basis gbasis=avdz extfil=.t. $end
$guess guess=huckel
$end
$trans mptran=2 dirtrf=.t. mode=112 cuttrf=1d-20 $end distci=256 maxv=4
$END
$data bqp rhf
Cs
CARBON 6.0 -1.447291768 0.000000322 0.0000000
CARBON 6.0 -0.679455541 -1.272277194 0.0000000
CARBON 6.0 -0.679454700 1.272276609 0.0000000
CARBON 6.0 0.679461301 -1.272276854 0.0000000
CARBON 6.0 0.679461775 1.272275816 0.0000000
CARBON 6.0 1.447289099 -0.000000524 0.0000000
HYDROGEN 1.0 -1.267419253 -2.196119151 0.0000000
HYDROGEN 1.0 -1.267417800 2.196118933 0.0000000
HYDROGEN 1.0 1.267440498 -2.196107501 0.0000000

23
HYDROGEN 1.0 1.267441310 2.196106291 0.0000000
OXYGEN 8.0 -2.688290088 0.000000793 0.0000000
OXYGEN 8.0 2.688279641 0.000000666 0.0000000
$end
Структура файла почти полностью совпадает со структурой файла для оптимизации геометрии на уровне RHF. Здесь типом расчета является вычисление энергии в точке при заданной геометрии (runtyp=energy), а раздел $statpt, посвященный расчёту градиента, который был в файле для оптимизации, отсутствует.
Обратите внимание, в данном и последующих расчетах мы используем симметрию
Cs, что и указано в начале карточки $data. Для корректного выполнения расчета молекулу необходимо расположить в плоскости XY.
Приложение 4. Пример входного файла для расчета энергии методом
CASSCF.
$contrl icharg=0 mult=1 scftyp=mcscf d5=.t. runtyp=energy maxit=100 icut=20 inttyp=hondo itol=30 gencon=.f. wide=1 fstint=.t. exetyp=run $end
$moorth syms=1 symden=1 symvec=1 tole=0 tolz=0 nostf=1 nozero=1 $end
$system mwords=360 timlim=100000 nojac=1 $end
$system mxbcst=-1 $end
$p2p p2p=.t. dlb=.t. $end
$basis gbasis=avdz extfil=.t. $end
$guess guess=moread norb=160 norder=1
Iorder(22)=23,24,27,28,29,30,31
Iorder(32)=-45
$end
$trans mptran=2 dirtrf=.t. mode=112 cuttrf=1d-20 $end
$MCSCF CISTEP=ALDET soscf=.t. focas=.f. MAXIT=100 fors=.t. acurcy=1d-7 $END
$DET NCORE=24 NACT=8 NELS=8 NSTATE=2 WSTATE(1)=1,1