量子化学计算方法-Gaussian的使用

量子化学计算方法 章永凡福州大学化学系2011年2月 课程主要内容 一 有限尺度体系 分子 团簇等 电子结构计算方法 G03程序的使用二 无限周期体系 一维链状化合物 二维层状化合物或固体表面 三维固体体相 电子结构计算方法 VASP CASTEP程序的使用主要参考资料 G03用户手册或G03的帮助文件相应网址 VASP程序用户手册 简介 量子化学软件目的在于将量子化学复杂计算过程程序化 从而便于人们的使用 提高计算效率并具有较强的普适性 绝多数量子化学程序是采用Fortran语言编写的 Fortran77或Fortran90 通常由上万行语句组成 软件分类 计算原理 基于从头算或第一性原理方法 abinitio firstprinciples Gaussian ADF Dalton Gamess Crystal VASP Wien Dmol Castep等 基于半经验或分子力学方法MOPAC EHMO NNEW3等 研究对象 有限尺度体系 分子 簇合物等 Gaussian ADF Dalton Gamess MOPAC EHMO等 无限周期重复体系 晶体 固体表面 链状聚合物等 Crystal NNEW3 VASP Wien等 本研究室目前常用的量化软件 Gaussian03 由Pople等人编写 经过几十年的发展和完善 该软件已成为国际上公认的 计算结果具有较高可靠性的量子化学软件 它包含从头算 半经验以及分子力学等多种方法 可适用于不同尺度的有限体系 除了部分稀土和放射性元素外 它可处理周期表中其它元素形成的各种化合物 Crystal09 该软件由意大利都灵大学理论化学研究所开发 采用基于原子轨道线性组合的从头算方法来研究固体及表面的电子结构 VASP 该软件由奥地利维也纳大学开发 采用基于平面波基组的密度泛函理论来研究固体及表面的构型以及动力学过程 CASTEP MS软件模块之一 与VASP程序类似 Dmol MS软件模块之一 主要用于有限尺度体系电子结构研究 采用理论方法要解决的问题 当前的研究状况 包括实验和理论研究现状 已解决和尚未解决的问题 一般计算过程 化合物构型的确定 具体途径包括 利用实验测定结果 或者采用软件进行构造等 根据现有的计算条件 模型的大小以及所要解决的问题 选择可行的计算方法和相应程序 对计算结果进行加工和提取有用的信息 一般包括构型描述 能量分析 轨道组成 电荷和成键分析等 并与实验结果比较 计算模型和方法的选取是保证计算结果可靠性的关键 理想的情况是 1 所选取的计算模型与实际情形一致 2 采用高级别的计算方法 但是 由于受到计算软硬件的限制 在多数情况下 很难同时做到上述两点要求 实际操作中 当计算模型较大时 只能选择精确度较低的计算方法 只有对较小的模型才能选取高级的计算方法 因此 当确定了一种计算模型和方法后 最好对其进行验证 以保证计算结果的可靠性 假设当前的研究对象是化合物A 可通过下列途径进行验证 1 与A化合物现有实验结果之间的比较 2 若无实验方面的报道 可对与A类似的化合物B进行研究 此时以B的实验结果作为参照 3 当上述方法行不通时 可以采用较大模型和较为高级的计算方法得到的计算结果作为参照 该方法主要用于系列化合物的研究 如对A1 A2 A3 先用大模型和基组对A1进行研究 然后以该结果为参照 确定计算量适中的模型和方法并应用于A1 A2 A3 Gaussian03程序的使用 G03的安装和运行 G03的功能和程序结构 输入文件的编写与主要功能的使用 补充说明 G03程序的安装和运行 1 G03程序的安装 1 确定运行平台 Windows或Linux 2 对Windows平台 直接运行setup exe 其余步骤按提示操作即可 也可将其它机器上将已安装好的G03直接拷贝到本机 但需设置运行环境 对Linux平台 a 若G03是经过压缩过的 文件结尾为gz 用gunzip命令解压 例如 gunzipg03 linux tar gzb 若G03是打包的 文件结尾为tar 用tar命令将其释放 例如 tarxvfg03 linux tar c 设置环境变量 以cshell为例 在用户根目录下的 cshrc文件添加下列内容 也可在执行g03前逐条运行 setenvg03root home USER 设置g03所在目录 根据实际情况修改 source g03root g03 bsd g03 login 激活g03运行时所需环境变量 setenvGAUSS SCRDIR home USER g03 tmp 设置临时目录 d 运行bsd install 自动配置并行计算环境注 对Linux平台 运行g03时 需注意权限问题 可用chmod命令更改权限 将所安装的g03对所有用户开放 对于多节点并行计算还需保证各节点之间能无密码互相访问 2 G03程序的运行 1 对Windows平台 a 对于刚安装好的g03 先检查环境设置情况 左侧至上而下依次为 默认的文本编辑器 g03可执行文件所在目录 计算中间结果存放目录 缺省的计算结果存储目录 缺省的输入文件所在目录 PDB分子构型浏览器 右侧至上而下依次为 设置显示属性 如背景色等 设置文本编辑器属性 计算过程控制属性 尤其是批作业过程 Default Rou文件的编辑 该文件内容为默认情况下 计算所花费的内存及硬盘大小 需设置正确 否则运行将出错 b 编写或打开g03输入文件 点击RUN 并给定输出文件名后开始运行 c g03运行过程的控制 最上行按钮的功能从左至右依次为 开始运行g03 暂停进程 运行至下一模块 link 时暂停进程 重新启动进程 清除进程 停止运算 编辑批作业 运行完当前任务后 暂停批作业 停止批作业的运算 观看计算结果 打开文本编辑器 不要随意点击 交换机 计算节点 计算节点 计算节点 网关 用户终端 基于Linux系统的计算拓扑结构 内部网 高速 外部网 普通 Window系统 Linux系统 2 Linux平台 说明 网关作用类似于防火墙 用于保证内部网的安全和稳定 作为网关的计算机通常配有2个网卡 分别用于外部网和内部网的连接 本实验室网关机子IP地址 219 229 135 210 非固定IP 可能发生改变 计算作业提交过程 a 用户登录网关 通过SSH远程登录软件实现SSH软件 SSHSecureShellClient 3 2 9 exe 可从网络上免费下载 安装过程与通常软件安装类似 安装完毕后 设置网关外部网的IP地址以及账号名即可使用 点击Profiles设置IP地址及用户名 b 从网关登录到计算节点 采用telnet命令实现例如 telnet134 14 83 5 c Linux常用命令 1 ls 显示文件清单 相当于DOS下的dir命令 文件属性 所属用户 大小 创建时间 注 Linux系统下字符是大小写区分的 2 cp 复制文件命令 相当于DOS下的copy命令 cp rf 3 mkdir 创建目录 相当于DOS下的md命令 4 rm 删除文件或目录 相当于DOS下的del命令 删除文件 删除目录 5 top 显示当前进程和CPU以及内存使用情况 6 kill 终止某个进程 格式为 killPID号 PID号由top命令可得 受权限限制 7 pwd 显示用户当前所在目录 格式为 pwd 8 cat 显示文件内容 格式为 cat文件名 9 grep 一般用于从某个或多个文件中搜索某串字符 格式为 grep 字符串 文件名例 grep F vasp out 10 scp 用于网关与内部网内各计算节点或外部网络之间的文件传输格式为 从其它到网关 scp文件用户名 网关IP 目录例 scpvasp outzyf zyf 2400 trans从网关到其它 scp用户名 网关IP 目录 文件名目标目录例 scpzyf zyf 2400 trans vasp out该命令也可用于同一节点不同用户之间的文件传输 11 vi 文本编辑命令该命令常用但较为复杂 它有2种模式 命令模式和插入模式 二者之间关系为 iEsccommandmode insertmode commandmode在命令模式下 可实现以下功能及其对应按键 deleteacharacter xdeletealine ddsearchastring 向后 向前 savethechange wsavethechangeandquit wqquitewithoutsaving q pagedown Ctrl dpageup Ctrl ugotofileend shift ggotonline n 12 tar 文件打包命令 适用文件扩展名为tar 该命令用于多个文件 目录的打包或解包 常用格式有 文件打包 tarcvf要打包成的文件名要打包的文件例 tarcvfmodel tar 生成model tar文件解包 tarxvf要解包的文件名例 tarxvfmodel tar 13 gzip gunzip 文件压缩 解压命令 文件扩展名为gz 例 gzipmodel tar生成model tar gz文件gunzipmodel tar gz生产model tar文件 14 du 察看当前目录所占硬盘空间大小 类似命令df 例 du hdf h 察看硬盘各分区大小 15 rsh或ssh 用于从某个节点登录到其它节点例 rshc0102登录到c0102节点上 为当前用户名rshzyf c0102与上相同 但用户名为zyfsshc0102sshzyf c0102根据rsh或ssh服务的具体设置来确定是否需要提供密码 16 su 从当前用户转变为超级用户或其它用户例 su转变为超级用户suzyf将用户转变为zyf用户 17 ifconfig 察看网络设置 18 dmesg 察看系统日志 19 adduser passwd d 运行g03过程 1 编写输入文件 用vi命令编写或在Windows下编写完毕后ftp至Linux系统 vitest gjf 2 运行g03 g03输出文件名 例 g03test out 说明 1 末尾的 符号表示将作业提交到后台计算 否则在用户退出Linux时 作业将终止 2 若运行g03出错 请检查环境变量是否设置正确 尤其是用户权限上的问题 3 观看计算结果 使用vi命令 或采用tail命令跟踪计算输出 tail f输出文件名 4 运行过程的控制 采用top命令观察g03运行到那个模块 通过renice命令改变进程的优先级来调整g03的运行速度 课堂练习安装G03Linux版本采用vi命令编辑Gaussian输入文件 具体内容如下 mem 32mb pb3lyp 6 311 G optGeomOptimizationofC2H40 1CC11 5H11 02120 0H11 02120 03180 0H21 01120 030 0H21 01120 03180 03 运行g03 采用top tail命令察看进程以及用vi和grep命令参看计算输出等 G03的主要功能和程序结构1 主要功能 分子构型的优化 基态 Groundstate 激发态 Excitedstate 反应过渡态 Transitionstate 能量计算 基态和激发态能量 化学键的键能 电子亲合能和电离能 化学反应途径和势能面 光谱计算 IR光谱 Raman光谱 吸收 发射光谱以及二阶或三阶非线性光学性质 NMR 其它功能 电荷分布和电荷密度 偶极矩和超极矩 热力学参数 适用体系 气相和溶液 2 程序结构 a 由主引导模块 g03 exe 和各分模块 l exe 组成 b 常用模块的功能 L0 初始化模块 L1 读入输入文件 根据所给关键词确定将要使用的模块 L101 102 与构型优化和反应过渡态相关的模块 L202 输出距离矩阵 判断化合物点群及确定新的坐标系 L301 302 309 与基组和赝势有关模块 L310 319 计算单电子及双电子积分模块 L401 402 SCF初始猜测模块 L502 503 508 SCF模块 L601 608 Mulliken布居以及自然键轨道 NBO 分析模块 L701 702 计算能量一阶和二阶导数模块 L8 9 10 11 与Post SCF方法有关模块 L9999 进程结束模块 说明 1 根据不同的任务 某些模块需重复调用多次 2 通常耗时较多的模块有 L5 L7 L8 L9 L10 L11等 此外 L8 L11这些模块的执行对内存和硬盘的需求较大 3 若L9999未能正常执行完毕 则表明计算过程存在问题 需检查之 4 可根据各个模块的功能 对g03程序进行简化 例如如果用户通常只用g03进行能量计算 则可只保留L1 6和L9999模块其它模块可以删除去 c g03运行过程所使用的文件 在scratch目录 或工作目录下有下列文件 gxx 打头的文件为临时文件 计算结束后将自动删除 其中对于结尾为inp的文件 记录了当前g03所执行的输入文件内容 有时可通过该文件确定当前运行作业 chk文件 该文件记录了g03运行的结果 包括分子结构 基组 分子轨道 电荷密度以及偶极矩等 通常该文件在计算结束后要保留 便于以后作补充计算或计算结果处理 rwf文件 该文件记录了计算的中间结果 以便在计算过程非正常中断后用于续算 该文件通常较大 当作业正常结束后 可删除之 Default Rou文件 该文件设置一些系统默认参数 例如g03运行时内存和硬盘的大小 其内容如下 M 256MB 内存大小 MaxDisk 2000MB 硬盘大小 说明 1 对于chk文件并不自动产生 需用户自行指定 在大多数情况下 最好给定chk文件 并在计算结束后保留 以便后续处理 例如计算结果的图像化等 2 对于rwf文件 在运行过程中 g03会自行产生gxx打头的rwf文件 当计算非正常中断后 可通过更改该文件来续算 但用户在编写输入文件时 最好还是指定rwf文件 在计算结束后 再删除 G03输入文件的编写与使用 1 G03输入文件的组成 L0命令部分 可无 关键词部分 标题部分 体系电荷和自旋多重度 分子构型 1 L0命令部分 该部分内容均以 打头 主要用于指定计算过程所需内存 chk以及rwf文件名 其典型内容是 mem 100mb 计算所需内存 若无 由Default Rou指定 chk h2o chk文件的名称为h2o chk rwf h2o 中间结果文件名称为h2o rwf NProcShared 4 多核系统中的cpu数目 LindaWorkers g0101 2 Linda并行计算节点及节点数目 2 关键词部分 该部分内容由一个或多个关键词组成 用于指定计算类型 方法和计算输出的控制等 3 标题部分 由一行文本组成 该内容是必需的 4 电荷及自旋多重度 体系所带电荷以及自旋多重度S 2s 1 成单电子数 1例如H2O S 1 Ti原子的基态S 3 5 分子构型描述部分 有三种描述分子构型的方法 a 直角坐标系方法 适用于全自由度构型优化情况 格式为 元素符号xyz例如 O 0 4640 1770 0H 0 4641 1370 0H0 441 0 1430 0说明 1 元素符号大小写均可 也可直接采用原子序数 2 有时为了便于区别 可在元素符号后加一整数 如 O 0 4640 1770 0H1 0 4641 1370 0H20 441 0 1430 03 x y z数值必须以小数格式输入 O 0 4640 1770 O 0 4640 1770 4 g03的数据输入均为自由格式 即除了用空格来分隔数据外 也可用逗号或混合使用 b 内坐标 z matrix 方法 适用于构型的局部优化 内坐标与直角坐标之间的区别在于 它侧重于从原子之间的键连角度来描述原子间的相对位置 具体参数包括 1 键长 需用两个原子描述 即两个原子间的距离 注 该两个原子并非要具有化学直观意义上的成键 此外 在默认情况下 键长单位为埃 2 键角 需用三个原子描述 确定了二根键之间的夹角 默认单位为度 范围为 180 180deg之间 3 二面角 需用四个原子描述 二面角加上键长和键角就确定了四个原子的位置 其默认单位为deg 范围为 360 360deg 当二面角等于0 180和 360deg时四个原子共面 键长 键角和二面角数目的总和 3N 6 内坐标的输入格式为 原子1 原子2 键长 原子3 键角 原子4 二面角 1 2 3 4 键长 键角 二面角 例1 O H H 表示一 OH 1 1 0O 1 1 2 2 104 0H 3 1 0 1 104 0 2 170 0表示二 O1H1 O1 1 0O2 O1 1 2 H1 104 0H2 O2 1 0 O1 104 0 H1 170 0 O 1 2 3 4 表示三 OH 1 r1O 1 r2 2 a1H 3 r1 1 a1 2 d1Variables 本行内容可省 r1 1 0r2 1 2a1 104 0d1 170 0 当对分子的构型进行局部优化时 需采用该表示方法 例2 乙烯 C C 1 2 3 4 5 6 CC11 3H11 02120 0H11 02120 03180 0H21 01120 030 0H21 01120 03180 0 在同侧共面 共面但不同侧 对同一构型 内坐标的表示并不唯一 C C 2 4 1 3 5 6 HC 1 1 0H 2 1 0 1 120 0C 2 1 3 1 120 0 3 180 0H 4 1 0 2 120 0 1 0 0H 4 1 0 2 120 0 1 180 0 虚原子的使用 有时为了保证所描述的构型符合特定的点群 利用虚原子便于做到这一点 虚原子的符号为X 例1 CO2 在g98版本 键角不能为180 0 C O O X 1 2 3 4 XC 1 1 0O 2 1 1 1 90 0O 2 1 1 1 90 0 3 180 0 该键长值可任意 C O O 2 1 3 例2 NH3 N H H H X 1 2 3 4 5 要使输入的构型满足C3v点群 需要准确提供H N H键角以及四面体相邻两个平面间的二面角 为此 在3个H所在三角形中心引入一个虚原子X 则 XN12 0H11 0290 0H11 0290 03120 0H11 0290 03 120 0 注 在本例中N X和H X不能任给说明 1 根据需要 有时可同时用到多个虚原子 2 在大多数场合 虚原子通常取在对称元素所处位置或它们相交处 c 直角坐标和内坐标混合输入方法 对于该方法 只需在采用直角坐标方法输入的原子的元素符号后加一个整数0即可 例如 X01 01 01 0N12 0H11 0290 0H11 0290 03120 0H11 0290 03 120 0 d 分子构型的输入准确性是保证计算结果可靠性的前提 对于复杂体系 在计算前均需对所输构型进行检查 具体包括 构型的可视化处理 即采用一些分子构型软件 例如Gaussview和Chem3D 观察所给构型是否合理 在g03运行到L2模块 会给出所输入分子所属点群 此时 可检查点群是否合理 6 g03输入文件编辑时的注意事项 除了可采用g03所提供的输入文件编辑器来编写输入文件外 在更多场合下 是采用其他文本编辑器来编写 此时应注意到 在标题部分的前后各有一空行 例如 chk h2o HF 6 31G d 此处为空行 waterenergy 此处为空行 01O 0 4640 1770 0H 0 4641 1370 0H0 441 0 1430 0 练习 采用内坐标方法输入苯 甲烷和立方烷的构型 并通过运行g03来验证对称性是否正确 CH11 H11 2aH11 2a3120 H11 2a3 120 a 109 47如果输入的键角小数点后没有三位有效数字 则程序判断为c3v群 此时 可以结合采用symm loose关键词来降低对精度的要求 甲烷内坐标 XC1aC1a260 C1a360 2180 0C1a460 3180 0C1a4120 2180 0C1a5120 3180 0H2b3120 7180 0H3b2120 4180 0H4b3120 5180 0H5b4120 6180 0H6b5120 7180 0H7b6120 2180 0a 1 42b 1 0 苯的内坐标 2 g03主要功能的使用 g03功能的使用主要由用户所给的关键词 keyword 内容而定 在输入关键词时注意以下事项 1 关键词的输入是自由格式 且不区分大小写 2 当存在多个关键词时 可用空格 逗号来隔开 3 通常每个关键词有多个选项 option 若要选择单个或多个选项时 书写方式有如下几种 keyword optionkeyword option keyword option1 option2 keyword option1 option2 例如 opt z matrix opt z matrix opt z matrix 三者是等价的 opt z matrix maxcycle 20 与opt z matrix maxcycle 20 是等价的 建议统一采用第四种表示方式 4 最简单的关键词输入是 或 p 其含义是采用HF方法和STO 3G基组计算体系的能量 a 能量的计算 如何计算一个体系的能量是获取分子各种性质的基础 因此首先来看如何计算体系的能量 即进行单点能计算 1 计算方法的选择 g03提供的常用计算方法有 1 半经验方法 关键词 AM1 PM3 CNDO INDO MINDO它们主要用于大的有机分子体系 由上百个原子组成 一般对于含金属体系不适用 这些方法只有在特殊场合适用 2 从头算 abinitio 方法 HF方法 即基于Hartree Fock原理的方法关键词 HF RHF UHF ROHF说明 i 当关键词为HF时 会自动根据自旋多重度选择RHF还是UHF ii ROHF为限制性开壳层HF方法 与UHF区别在此时除了成单电子外 其余的 和 电子仍配对 通常该方法得到的能量要较UHF略高 iii HF方法可以看作是最低级的从头算方法 该方法除了在构型优化时有使用外 不适合计算能量 密度泛函方法 DFT 基于电荷密度自洽的方法关键词 B3LYP等根据所采用的相关和交换泛函 可以选择不同的DFT方法 具体参见g03的帮助文件 其中B3LYP方法是使用最为广泛的DFT方法 由于DFT方法考虑了电子之间的相关作用 因此得到的能量要较HF来得精确 它是目前最常用的量子化学计算方法 交换泛函 Slater S X XA Becke88 B Perdew Wang91 PW91 PBE PBE 等相关泛函 VWN LYP PW91 P86 PBE等组合形式如 BP86 SVWN PW91PW91 PBEPBE等 MPn方法 关键词 MP2 MP3 MP4 MP5说明 i 这些方法在HF基础上 进一步根据MP微扰理论考虑电子相关作用 微扰项截至到二阶则为MP2 截至到三阶则为MP3 其它类推 理论上考虑的微扰项越多 得到的能量越精确 但将大大增加计算量 而且通常也无此必要 多数场合选取MP2即可 ii 对于该类方法 硬盘和内存通常开销较大 应考虑具体的硬件考虑之 其中对于硬盘空间的设置见文件Default Rou内容 另外 必须注意到由于受到操作系统的限制 中间文件不能超过2GB 32位系统 此时需设置多个中间文件 具体见g03说明 耦合簇 CoupledCluster 方法 关键词 CCD CCSD T 说明 i 该类方法与MPn方法一样 也是属于较高精度的计算方法 其中CCD方法 只考虑了双取代 CCSD则在CCD基础上进一步考虑了单取代 ii 与MPn方法类似 该类方法计算量较大 通常只适用小体系 其它后自洽场方法 组态相互作用 Configurationinteraction 方法 关键词 CID CISD QCISD CASSCF等 此外 g03还提供一些高精度的组合计算方法 如G1 G2等 练习 分别用HF B3LYP MP2 MP4 CCSD方法计算H2基态能量 假使H H键长为0 7A 采用STO 3G基组 基态的电子态为 1 g HF能量 单位a u 2 基组的选择 全电子基组 赝势基组 高斯型函数 gp p r Npexp pr2 原子轨道可表示为高斯型函数的线性组合 1S p 1 L dpgp p r 系数 指数 最终体系分子轨道为这些原子轨道的线性组合 1 全电子基组 关键词 sto 3g 3 21g 4 31g 6 21g 6 31g 6 311g d95 d95v cc pV UGBS等说明 I 不同的基组适用范围是不同的 如 STO 3G H Xe 3 21G H Xe 6 21G H Cl 4 31G H Ne 6 31G H Kr 6 311G H Kr D95 H Cl除了Na Mg 说明 ii 基组的大小决定了基函数的数目 即体系的原子轨道数目 因此可从所选择的基组来推断MO数目 sto 3g 为最小基组 每个原子轨道用三个高斯函数 GF 来描述 原子轨道数即为基函数数目 如O 1s2s2p 原子轨道数为1 1 3 5GF数目为3 5 153 21g 为劈裂 split 型基组 其含义是 内层的每个AO用3个GF描述 价层的AO劈裂为两组 分别用2个和1个GF描述 显然 3 21g的GF数与sto 3g是相同的 如O 内层为1s AO数为1 GF数为3价层2s的AO数为2 1 2 GF数为2 1 3价层2p的AO数为2 3 6 GF数为3 2 3 1 9共1 2 6 9个AO和3 3 9 15个GFMg 1s2s2p3s3p内层1s 2s和2p共有1 1 3 5个AO和3 5 15个GF价层3s有2个AO和3个GF 价层3p有6个AO和9个GF共5 2 6 13个AO和15 3 9 27个GF 6 31g 为劈裂型 split 基组 其含义与3 21g类似 内层的每个AO用6个GF描述 价层的AO劈裂为两组 分别用3个和1个GF描述 如O 内层为1s AO数为1 GF数为6价层2s的AO数为2 1 2 GF数为3 1 4价层2p的AO数为2 3 6 GF数为3 3 3 1 12共1 2 6 9个AO和6 4 12 22个GFMg 1s2s2p3s3p内层1s 2s和2p共有1 1 3 5个AO和6 5 30个GF价层3s有2个AO和4个GF 价层3p有6个AO和12个GF共5 2 6 13个AO和30 4 12 46个GF 6 311g 也为劈裂基组 只是价层的AO劈裂为3组 分别用3个 1个和1个GF描述 对于4 31g和6 21g类似 iii 极化 polarization 函数的使用 在实际计算中 有时需在上述标准基组的基础上 添加一个或多个极化函数 极化函数是指具有比原子价轨道更高角量子数的高斯函数 例如H的价轨道为1s 则其极化函数为p型GF 同样对C O等价层为p轨道的原子 它们的极化函数应为d型或f型轨道 类似地 对于过渡金属原子的极化函数为f型轨道 极化函数的使用目的在于使原子价轨道在空间取向上变得更 柔软 从而使之易于与其它原子的轨道成键 例如对于羰基基团中的C O原子 它们极化轨道 即d轨道 之间形成的d 轨道 使得C原子的p 轨道朝O方向极化 O的p 轨道向C方向极化 从而增强了C与O之间的 作用 由于极化轨道的使用增强价轨道的空间柔软程度 因此其使用场合主要在于环状化合物 例如有机环状化合物以及含有桥联结构的金属簇合物等 使用方法 在标准基组后加 d d p 3df 3pd 等例 6 31G 6 31G df pd 说明 当只添加一个极化函数时 可使用下述表示 等价于 d 含义是对非氢原子添加一个极化轨道 等价于 d p 含义是对非氢原子添加一个极化轨道的同时对H原子添加一个p极化函数 当添加多个极化函数时 使用下述表示 d p 表示对非H原子和H原子分别添加 个极化轨道如 6 31G 3d 2p df pd 表示对非H原子添加 个极化函数和1个具有更高角动量的极化函数 对H原子分别添加 个p极化函数和1个d极化函数 例如 对H2O计算时采用6 31G 2df 3pd 基组 其含义是 对其中的H原子 在6 31G基组上在添加3个p轨道和1个d轨道 对O原子 则在6 31G基组上添加2个d轨道和1个f轨道 对于不同的基组 可添加的极化函数数目是不同的 具体可参照程序手册 iv 弥散函数 diffuse 的使用 对于带有较多电荷的体系 采用标准的基组来描述是不够的 此时需添加弥散函数 以增加价轨道在空间上的分布范围 即 极化函数用于改进价轨道的角度分布 弥散函数则用于改进价轨道的径向分布 弥散函数是指具有较小轨道指数的高斯函数 其表示方法是在标准基组后加上 或 如6 31 G 6 31 G 等 其中第一个 表示对非H原子添加弥散函数 第二个 则对H原子添加弥散函数 弥散函数主要用于带有电荷的体系 包括离子 以及弱作用体系 2 赝势基组 关键词 lanl1mb lanl2mb lanl1dz lanl2dz cep 4g cep 31g cep 121g sdd等说明 通常对于重元素采使用赝势基组 其中 对于过渡金属原子一般使用lanl系列的基组 对于主族元素采用cep系列的基组 sdd基组相对较少使用 赝势基组的适用范围要较全电子基组来得广 一般除了少数稀土和放射性元素外 均可使用赝势基组 赝势基组是将内层原子的影响用势函数来代替 因此 在使用时必须注意根据所采用的赝势基组 检查体系的总电子数是否正确 对lanl1系列的基组 只考虑价层电子 如V原子为5 对lanl2系列的基组 除了价层电子外 还考虑了次外层的电子 对部分主族元素例外 如V原子为5 8 13 对于使用较多的lanl系列基组 需注意 对于H Ne范围的原子 lanl1mb和lan2mb基组等价于sto 3g最小基组 lanl1dz和lanl2dz基组等价于d95基组 对于含金属 金属键体系一般使用lanl2系列 而不用lanl1系列mb结尾含义为 minimalbasis dz结尾含义为 double zeta 故后者基函数数目要较前者多 3 混合基组的使用 有时 因需要 对同一化合物中的不同原子采用不同的基组 受基组限制 例如 对于Mo CO 2 化合物 C和O原子可采用6 31 G基组 但Mo原子不在6 31G基组的使用范围之内 受体系大小限制 当体系较大时 不可能对所有原子均采用较大的基组 此时可对其中的局部原子采用精度较高的基组来描述 而对其余原子采用小基组来描述 混合基组输入方法 关键词 GEN和PSEUDO READ 对赝势基组 mem 32mb chk hf pb3lyp genpop fullscfcyc 500Mo CO 2 Line 1 6MoX11 0C1r1290 0O1r2290 030 0C1r1290 03180 0O1r2290 03180 0r1 2 228r2 3 329Mo03 21G CO06 31 G 整数0为结束标志 4个 用于分隔不同基组 该行也可用数字代替 即23450 编号时不包含虚原子 mem 32mb chk hf pb3lyp genpseudo readpop fullscfcyc 500Mo CO 2 Line 1 6MoX11 0C1r1290 0O1r2290 030 0C1r1290 03180 0O1r2290 03180 0r1 2 228r2 3 329Mo0Lanl2mb CO06 31G Mo0Lanl2mb 使用赝势基组时 应注意不能漏 基组定义部分 赝势定义部分 注意不能漏 另外其后没有 号 赝势部分和基组部分之间有一空行 4 补充说明 可采用关键词GFINPUT或GFPRINT 对于赝势基组 输出赝势参数 输出基组的内容 可用于对基组的优化 如O的6 31G 基组如下 S61 000 000000000000 为s型AO 6个GF 0 5484671660D 040 1831074430D 02 分别为GF的指数和系数 0 8252349460D 030 1395017220D 010 1880469580D 030 6844507810D 010 5296450000D 020 2327143360D 000 1689757040D 020 4701928980D 000 5799635340D 010 3585208530D 00SP31 000 000000000000 为sp型 GF共享同一个指数 0 1553961625D 02 0 1107775495D 000 7087426823D 010 3599933586D 01 0 1480262627D 000 3397528391D 000 1013761750D 010 1130767015D 010 7271585773D 00SP11 000 0000000000000 2700058226D 000 1000000000D 010 1000000000D 01D11 000 0000000000000 8000000000D 000 1000000000D 01 练习 输出W原子Lanl2MB和Lanl2DZ基组的高斯轨道指数和系数 观察指数和系数的变化规律 并对两种基组进行比较 3 G03输出结果 H2O能量计算的输出结果为例 L1输出G03版权信息 部分内容 关键词部分 计算所需调用的模块 Title内容 电荷 S以及分子构型 各原子的直角坐标 距离矩阵 需核查 注意 计算结果是以该坐标系为准 需核查 MO初始猜测 自洽场迭代求解部分 收敛指标 迭代次数 通常接近2 0 布居分析结果 MO对称性及能级 a u Mulliken键级和电荷 偶极矩及多极矩 对主要计算结果进行总结 4 能量计算中使用到的其它常用关键词 SCF收敛指标的设置 关键词 scf conver 含义 当SCF前后电荷密度的差值小于10 时 SCF收敛 例 scf conver 5 表示收敛指标为10 5 说明 对于能量计算 缺省情况下收敛指标为5x10 5 对于构型优化 则为10 8 SCF迭代轮数的设置 关键词 scf maxcyc 含义 SCF最大的允许迭代次数为 例 scf maxcyc 100 表示最大迭代次数为100次 说明 缺省情况下 SCF迭代次数为128次 通常对于含有过渡金属的体系 SCF较难收敛 迭代次数会多些 分子轨道系数的输出 关键词 pop full说明 有时为了分析MO成分 则需利用该关键词输出各MO的成分 以H2O为例 轨道对称性 O表示占据轨道 V为空轨道 轨道能级 AO轨道 轨道系数 HOMO LUMO 注 输出的MO系数是未经归一化过的 1s 2s等符号并非真正意义上的AO 它与计算所采用的基组有关 对于上例 采用的是6 31G基组 可知 对于价层AO是分裂为两组的 故对于O原子 输出的2s和3s实际上均为价轨道 其它类推 当基函数较多时 若只考察前线附近的轨道成分 可用关键词 pop regular 此时只给出5个占据轨道和5个空轨道组成 练习 计算O2的基态能量 并分析MO轨道组成假设基组为6 31G O O键长为1 2145A 5 能量计算的几个问题 基组的选取 理论上而言 计算中所选取的基组越大 计算结果越准确 但由于受到硬件上的限制 需根据实际情况选择基组 此外 对于环状分子或存在 作用的体系 通常需要考虑使用极化基组 对于带有较多电荷的体系 或考察弱相互作用的体系 通常需考虑使用弥散基组 以O2为例 O O为1 208A 6 31G 150 26755a u 6 31G 150 31998a u dE 0 05243a u 6 31 G 150 27678a u dE 0 00923a u 由于O2中存在 作用 故极化函数的影响要比弥散函数来得大 以O22 为例 6 31G 149 78721a u 6 31G 149 82396a u dE 0 03675a u 6 31 G 149 98117a u dE 0 19396a u 6 31 G 150 02156a u dE 0 23435a u 此时 因带较多负电荷 故弥散函数的影响较极化函数来得大 体系多重度的选取 对于复杂体系 特别是含有过渡金属原子的体系 多重度的选取凭经验不易得到 此时 可通过比较不同多重度时体系的能量来确定基态能量 以O2为例 基组为6 31G S 1 150 25734a u S 3 150 31998a u S 5 149 79200a u S 7 148 80784a u 由上结果可见 O2的基态应为三重态 自旋污染判断 对开壳层体系的计算 需判断自旋污染情况 以便保证计算结果的可靠性 s s 1 例如对O2基态 有2个单电子 s 2 1 2 1 则 1 1 1 2 对上例子s 3 2 体系电子态的选取 多数情况下 体系的电子态 组态 与体系的点群和电子占据情况有关 由程序自动确定 但在某些情况 特别是在团簇体系 采用这种默认方法得到的电子态可能出现问题 其中包括 无法确定轨道的对称性或者电子填充顺序出现问题 此时可人为对电子的填充轨道进行调整 所采用的方法为Guess关键词中的Alter选项 此外 采用该方法也可实现对特定激发态进行研究 pb3lyp6 31GoptNH20 2NH1rH1r2ar1 03a120 0 NH2自由基的基态 NH2基态的电子填充情况为 a12a12b22a12b1 即相应的电子态为2B1 相应能量为 55 85291a u 如果我们想考察NH2的2A1激发态 则可对电子填充顺序进行人为调整 具体的办法是将 态中的4号和5号轨道进行对调 使得最终电子填充顺序为 a12a12b22b12a1 即体系电子态为2A1 pb3lyp6 31GoptGuess AlterNH20 2NH1rH1r2ar1 03a120 045 该行前面有一空行 对应于 态电子填充的改变情况 最终计算结果 键能的计算 以O2为例 计算O O键能 步骤如下 首先确定O2基态的能量 假设基组为6 31G 则总能量为 150 32005a u 基态为 3 g O O 1 2145A 其次 计算O原子基态能量 结果为 75 06062a u 则O2中O O键的键能为 2 75 06062 150 32005 0 19881a u 522 0KJ mol该结果略高于实验值 494KJ mol 在有关键能的计算中 除了计算方法以及基组对最终的结果有较大影响外 若要得到精确的键能 还需考虑基组叠加误差即BSSE BSSE是指源于能量计算中基组带来的误差 它可以粗略理解为化学实验中的背景误差 要消除BSSE 可采用counterpoise关键词来实现 G98是采用message关键词 以O2为例具体使用过程如下 mem 32mb pb3lyp 6 31g counterpoise 2BSSEofO20 3 0 3 0 3O0 00 00 01O1 21450 00 02 表示体系由2个碎片组成 整个体系的电荷和自旋多重度碎片1的电荷和自旋多重度碎片2的电荷和自旋多重度 该原子属于那个碎片 由此得到考虑了BSSE时的键能为 0 19881 0 00436 0 19445a u 510 4KJ mol 此外 进一步的修正是在计算O2总能量时 需考虑零点能 ZPE 校正 上机练习 采用B3LYP方法和6 31G 基组 计算CO中的C O键能 C O键长为1 138 分子轨道成分的分析 以O2为例 mem 32mb pROB3LYP STO 3GPOP FULLThesinglepointcalc OfO20 3O0 0 00 0O1 2080 00 0 为了便于计算结果分析 采用限制性开壳层计算方法 及将 和 电子强制配对 输出所有分子轨道 首先确定计算结果中各原子所处的位置 2 1 o z Xory z 1 2 因1s为内层电子 故O与O之间成键情况可不考虑 No 1 No 2与上一MO类似 也以1s成分为主 处在能量较低的位置 z z 1 2 No 3 4分别为O2s轨道之间形成的 成键和 反键轨道 z No 5为O2pz轨道形成的 成键轨道 z z No 6和7轨道分别为O2px和2py形成的 成键轨道 No 8和9轨道分别为O2px和2py形成的 反键轨道 z No 10轨道为O2pz形成的 反键轨道 z 根据上述分析结果以及各轨道的能级 可以得到O2的MO分布为 2s 2s 2pz 2px 2py 2px 2py 2pz 电荷分布 在能量计算结果中 除了能量以及上述的MO成分外 其它有用的信息包括原子电荷和键级 对于原子电荷 Mulliken电荷分布 说明 电荷的绝对值是没有意义的 其数值受到所用方法 尤其是基组的影响较大 以H2O为例 方法 基组OHB3LYP STO 3G 0 3290 165B3LYP 3 21G 0 6370 318B3LYP 6 31G 0 7050 353HF 6 31G 0 7950 398因此 计算得到的电荷只具有相对意义 即只能考察其相对数值或变化趋势 当然 当考察变化趋势时 如考察H2X X O S Se 系列化合物中 X为不同原子时 对H原子电荷的影响时 必须使用相同的计算方法和基组 键级 键级的大小反映了某根键的强度 同样地 键级的数值也受到计算方法和基组的影响 O H键的Mulliken键级 以H2O为例 方法 基组O H键级B3LYP STO 3G0 243B3LYP 3 21G0 239B3LYP 6 31G0 231B3LYP 6 31 G0 214HF 6 31G0 249由上结果可见 键级数值与所用的基组和计算方法是有联系的 但其影响不象电荷那么显著 此外 必须注意到 弥散基组对键级的数值影响较大 因此 在实际应用中 对于电荷和键级 我们更关心的是它们的变化情况而不是它们的绝对数值 此外 必须强调的是 上述结果是采用Mulliken布居分析的方法 若采用其它布居分析方法 结果将有所不同 自然键轨道 NBO 分析方法 我们前面所讲分子轨道均为正则分子轨道 它是未经定域化处理 这导致计算结果与我们通常的成键概念有所不同 例如在乙烯分子中 碳碳之间为双键 但在正则MO中 反映C与C之间成键作用的MO可能有多个 因此根据正则MO的结果 我们无法断定C C是单键还是双键 此时 通过对正则MO的定域化处理 可以得到通常意义上的成键图像 正则MO的定域化处理方法较多 其中较为常用的是NBO方法 其使用方法是在输入文件中添加关键词 POP NBO mem 32mb pb3lyp 3 21gpop nboTheNBOanalysisofethylene0 1CC 1 CCH 1 CH 2 HCCH 1 CH 2 HCC 3 180 0H 2 CH 1 HCC 3 180 0H 2 CH 1 HCC 4 180 0CC 1 31477CH 1 07363HCC 121 8867 进行NBO成键分析 NBO计算结果的主要内容 自然布居分析结果 自然电荷 自然电子