使用FCClasses3程序计算系间窜越速率
Calculating Intersystem Crossing Rate via FCClasses3 Program
Jiawei Xu
Released: 2023-06-12 / Updated: 2023-06-12
计算系间窜越(Intersystem Crossing, ISC)速率是准确考察多态反应机理的关键之一，其方法并不唯一。本文将介绍由意大利金属有机化学研究所Fabrizio Santoro等开发的FCClasses程序的安装以及计算ISC速率的方法。
1. 程序安装
FCClasses的最新版本为3.0.1，可在http://www.iccom.cnr.it/en/fcclasses登记并下载程序源代码和手册。FCClasses3是Fortran写的，需要gfortran或ifort编译器，支持LaPack或MKL数学库，默认使用LaPack，同时可选FFTW快速傅立叶变换库。
tar -xf fcclasses3-3.0.1.tar.gz
mkdir fcclasses3
cd fcclasses3-3.0.1
./configure --prefix=/home/jwxu/programs/fcclasses3 --with-mkl --with-fftw (如不使用MKL或FFTW，不写即可)
make -j && make install -j
make test
cd ..
rm -rf fcclasses3-3.0.1*make test后会运行一系列任务检查是否安装正确，一般不需要全部运行完，如前几个任务都显示成功运行，就可以结束测试。安装完成后，所有的可执行文件都在/home/jwxu/programs/fcclasses3/bin下，并需要在.bashrc中添加环境变量：export PATH=$PATH:/home/jwxu/programs/fcclasses3/bin。FCClasses的绝大多数任务都可以直接在登陆节点上运行。
2. 计算系间窜越速率
以S1-T1之间的系间窜越过程为例，首先通过Gaussian16进行构型优化和频率计算，得到.fchk文件。FCClasses需要提供结构、梯度和Hessian矩阵信息，由其自带的脚本可以生成对应格式的输入文件：
gen_fcc_state -i S1.fchk
gen_fcc_state -i T1.fchk以上命令产生了S1.fcc和T1.fcc两个文件。准备以下输入文件fcc.inp：
$$$
PROPERTY = NR0 ; OPA/EMI/ECD/CPL/RR/TPA/MCD/IC
MODEL = AH ; Adiabatic Hessian; AS/ASF/AH/VG/VGF/VH
DE = 0.150 ； Energy difference in eV
NR0_COUPL = 5.00E-6 ; Spin-orbit coupling integral in a.u.
TEMP = 298.15 ; Temperature in K
BROADFUN = GAU ; GAU/LOR/VOI
HWHM = 0.06 ; Broadening width in eV
METHOD = TD ; Time-dependent; TI/TD
ROT = 1 ; Rotate state 1
;VIBRATIONAL ANALYSIS
NORMALMODES = COMPUTE ; COMPUTE/READ/IMPLICIT
COORDS = CARTESIAN ; CARTESIAN/INTERNAL
;INPUT DATA FILES
STATE1_FILE = S1.fcc
STATE2_FILE = T1.fcc注意修改输入文件中的两态能量差（末态减去初态）以及旋轨耦合积分的数值。然后运行：
fcclasses3 fcc.inp运行成功后得到一系列输出文件，可从fcc.out读取系间窜越速率，同时目录下的其他输出文件还提供了Huang-Rhys因子、Duschinsky矩阵以及位移矢量的相关信息。
========================================================
Non-radiative rate constant (s-1) 1.425E+08
(for Ead = -0.177 eV)
========================================================3. 便捷脚本
以下是一个全自动计算ISC速率的脚本，可供参考：
#!/bin/bash
# Usage: ./runfcc mol1 S1 T1
mol=$1
state1=$2
state2=$3
exc="S1 T2"
mkdir -p $mol/fcc-$state1-$state2; cd $mol/fcc-$state1-$state2
### Generate .fcc inputs
cp ../$state1/$mol-$state1.chk ../$state2/$mol-$state2.chk .
for i in *.chk
do
name=`basename $i .chk`
fchk16 $i $name.fchk
gen_fcc_state -i $name.fchk
done
rm *.chk
### Get energy difference
if [[ "$exc" =~ "$state1" ]]
then
set `grep "E(TD-HF/TD-DFT)" ../$state1/$mol-$state1.log | tail -1`
else
set `grep "SCF Done" ../$state1/$mol-$state1.log | tail -1`
fi
E1=$5
if [[ "$exc" =~ "$state2" ]]
then
set `grep "E(TD-HF/TD-DFT)" ../$state2/$mol-$state2.log | tail -1`
else
set `grep "SCF Done" ../$state2/$mol-$state2.log | tail -1`
fi
E2=$5
dE=`echo "scale=10; ($E2)-($E1)" | bc`
dE=`echo "$dE*27.2114" | bc`
echo "The energy difference is $dE eV."
### Get spin-orbit coupling integral from DALTON output
grep "B excited state no., symmetry, spin: ${state1:0-1}" ../../dalton/$state1*dalton/$mol/*.out -A 3 -B 2 > null1
grep "C excited state no., symmetry, spin: ${state2:0-1}" null1 -A 3 -B 2 > null2
set `grep X1SPNORB null2 -A 5`; socx=${39}
set `grep Y1SPNORB null2 -A 5`; socy=${39}
set `grep Z1SPNORB null2 -A 5`; socz=${39}
soc2=`echo "scale=20; ($socx^2)+($socy^2)+($socz^2)" | bc`
soc=`echo "sqrt($soc2)" | bc`
rm null1 null2
echo "The spin-orbit coupling integral is $soc a.u."
### Generate FCClasses main input
rm -f $mol-fcc-$state1-$state2.inp
cat >> $mol-fcc-$state1-$state2.inp << EOF
\$\$\$
PROPERTY = NR0 ; OPA/EMI/ECD/CPL/RR/TPA/MCD/IC
MODEL = AH ; AS/ASF/AH/VG/VGF/VH
DE = $dE
NR0_COUPL = $soc ; SOC
TEMP = 298.15 ; (temperature in K)
BROADFUN = GAU ; GAU/LOR/VOI
HWHM = 0.06 ; (broadening width in eV)
METHOD = TD ; TI/TD
ROT = 1 ; Rotate state 1
;VIBRATIONAL ANALYSIS
NORMALMODES = COMPUTE ; COMPUTE/READ/IMPLICIT
COORDS = CARTESIAN ; CARTESIAN/INTERNAL
;INPUT DATA FILES
STATE1_FILE = $mol-$state1.fcc
STATE2_FILE = $mol-$state2.fcc
EOF
### Run FCClasses
fcclasses3 $mol-fcc-$state1-$state2.inp