电子结构计算方法和计算软件-刘晓辉.pdf
电子结构计算方法和计算软件-刘晓辉.pdf
电子结构计算方法和计算软件 软件系列一：Gaussian和GaussView使用介绍 刘晓辉 (liuxiaohui@ustc.edu.cn) 中国科学技术大学 超级计算中心 2023年5月31日 高性能计算、超级计算、科学计算 高性能计算 ◼ HPC: High Performance Computing • 通过互联技术将多个计算机系 统连接在一起，用来处理大型 计算问题 • 更快的计算速度、更大的负载 能力、更高的可靠性 • 主要组件：计算、网络和存储 随着高性能计算技术的快速发展，高性能计算（计算 科学）成为继理论科学和实验科学之后，人类科学研 究的第三大范式。 超级计算 ◼ Super Computer • 具有很强的计算和处理数据能力 的计算机。一般来说（2021年 ）， 5.2P 科 总 大 超 算 – 运算速度平均每秒1万亿次以上 1.6/1.1E 新一代神威？ 125/93P – TFLOPS：每秒万亿次，1012 – PFLOPS：每秒千万亿次，1015 – EFLOPS：每秒百亿亿次，1018 神 威 Frontier ➢ 计算能力：每秒浮点运算次数 （FLoating-point Operations Per Second, FLOPS ） 科学计算 ◼ 科学计算 • 科学研究，计算模拟 • 是指使用计算机数学建模和数值 分析技术解决科学与工程领域中 的数学计算问题 • 典型应用领域：物理、化学、材 料科学、地球与环境等，人工智 能与大数据处理 NERSC 应 用 分 布 应 用 领 域 计 算 需 求 电子结构计算软件 占比>60% 2022年TC4600系统 应用软件使用机时分布 电子结构计算方法简介 基于量子力学基本原理 使用模拟计算的方法 研究体系中电子的运动状态 电子结构计算方法 • 量子力学 1933年 诺贝尔物理学奖 微观粒子的状态 描述和运动规律 Erwin Schrödinger Paul Dirac 电子结构计算方法 量子力学基本概念：波函数 • 用来描述量子力学体系的状 态，是体系内所有粒子的位 置以及时间的函数 • 波函数的模方表示在空间某 处粒子出现的几率 • 通常用符号𝜓 或Ψ 表示 • 包含了实验上所有可观测量 的信息 氢原子的波函数 归一化条件： 𝛹 2 𝑑𝜏 = 1 Ψ 𝑥1 , 𝑦1 , 𝑧1 , 𝑥2 , 𝑦2 , 𝑧2 , … , 𝑥𝑁 , 𝑦𝑁 , 𝑧𝑁 , 𝑡 = Ψ 𝒓1 , 𝒓2 , … , 𝒓𝑵 , 𝑡 电子结构计算方法 量子力学基本概念：算符 • 任何一个可观测的物理量都有对应的算符，比如 位置、能量、动量等 • 算符的作用是将一个函数变换成另外一个函数， 𝒅 比如微分算符： 𝒔𝒊𝒏𝒙 = 𝒄𝒐𝒔𝒙 𝒅𝒙 • 常用的量子力学算符有： 位置 势能 ෝ 𝑥⇔𝒙 𝑉 𝑥 ⇔𝑽 动量 ෝ𝒙 𝑝𝑥 ⇔ 𝒑 ෝ𝜓 𝑥 = 𝒙𝜓 𝑥 𝒙 𝜓 𝑥 = 𝑽 𝒙 𝜓 𝑥 𝑽 𝝏 ෝ𝒙 𝜓 𝑥 = −𝒊ℏ 𝒑 𝜓 𝑥 𝝏𝒙 电子结构计算方法 量子力学基本概念：算符 ◼ 物理量𝑄，对应的算符为𝑄，其期望值计算公式为： 𝑄 = Ψ∗ 𝑄Ψ𝑑𝜏 Ψ∗ Ψ𝑑𝜏 算符的期望值 即：测量的平均值 波函数满足归一化条件： Ψ ∗ Ψ𝑑𝜏 = 1 𝑄 = න Ψ ∗ 𝑄Ψ𝑑𝜏 电子结构计算方法 量子力学基本概念：算符 ◼ 一个特殊的算符，本征值方程 𝑄Ψ = 𝑞Ψ 𝑞 是一个实数 本征函数 本征值 𝑸 = න Ψ∗ 𝑄Ψ𝑑𝜏 = න Ψ ∗ 𝑞Ψ𝑑𝜏 = 𝑞 න Ψ∗ Ψ𝑑𝜏 = 𝒒 电子结构计算方法 量子力学基本概念：算符 ◼ 哈密顿量算符：𝐻，体系总能量对应的算符 𝐻 = 𝑇 + 𝑉 ℏ2 𝜕 2 𝜕2 𝜕2 ℏ2 2 𝑇 = − + 2+ 2 =− 𝛻 2 2𝑚 𝜕𝑥 𝜕𝑦 𝜕𝑧 2𝑚 𝑞 1 𝑉 = 4𝜋𝜖0 𝑟 电子结构计算方法 量子力学基本方程 ◼ 薛定谔方程 𝜕 ℏ2 2 𝑖ℏ Ψ 𝒓1 , 𝒓2 , … , 𝒓𝑵 , 𝑡 = − 𝛻 + 𝑉 𝒓, 𝑡 Ψ 𝒓1 , 𝒓2 , … , 𝒓𝑵 , 𝑡 𝜕𝑡 2𝑚 电子结构计算方法 量子力学基本方程 ◼ 定态薛定谔方程 如果势函数不显含时间变量，即𝑉 𝒓 通过分离变量法，可得， 𝐻𝛹 𝒓𝟏 , 𝒓𝟐 , … , 𝒓𝑵 = 𝐸𝛹 𝒓𝟏 , 𝒓𝟐 , … , 𝒓𝑵 即：定态薛定谔方程（概率密度与时间无关） Ψ 𝒓1 , 𝒓2 , … , 𝒓𝑵 , 𝑡 Ψ 𝒓1 , 𝒓2 , … , 𝒓𝑵 , 𝑡 𝒊 − 𝑬𝒕 = Ψ 𝒓1 , 𝒓2 , … , 𝒓𝑵 𝒆 ℏ 2 = Ψ 𝒓1 , 𝒓2 , … , 𝒓𝑵 2 电子结构计算方法 电子结构的基本近似 玻恩-奥本海默近似：原子核运动 & 电子运动 Ψ 𝒓, 𝑹 = Ψ𝑛𝑢𝑐 𝑹 Ψ𝑒𝑙𝑒𝑐 𝒓 • 电子的运动方程： 𝐻𝑒𝑙𝑒𝑐 𝑹 𝛹𝑒𝑙𝑒𝑐 𝒓 = 𝐸𝑒𝑙𝑒𝑐 𝑹 𝛹𝑒𝑙𝑒𝑐 𝒓 其中，电子的哈密顿量算符： 𝐻𝑒𝑙𝑒𝑐 𝑹 = 𝑇𝑒𝑙𝑒𝑐 + 𝑉 𝒓; 𝑹 质子、电子质量比 𝒎𝒑 Τ𝒎𝒆 ≈ 𝟏𝟖𝟑𝟔 电子结构计算方法 电子结构的基本近似 N-电子体系的哈密顿量：复杂、难求解 𝑁 𝐻𝑒𝑙𝑒𝑐 = − 𝑖 ℏ2 2𝑚𝑒 𝑀 𝛻𝑖2 − 𝐼 𝑒 2 𝑍𝐼 4𝜋𝜖0 𝒓𝑖 − 𝑹𝐼 𝑁 + 𝑖,𝑗 𝑒2 4𝜋𝜖0 𝒓𝑖 − 𝒓𝑗 𝒊 𝒉 单电子近似：假设电子之间没有相互作用， 单电子运动方程： ℎ 𝑖 𝜓𝑖 = 𝜀𝑖 𝜓𝑖 多体波函数： Ψ = 𝜓1 𝜓2 … 𝜓𝑛 体系的总能量： 𝐸 = 𝜀𝑖 𝑖 变分原理、自洽场（SCF）方法 𝑵 𝒊 𝑯𝒆𝒍𝒆𝒄 = 𝒉 𝒊 Hartree 积（HP） 𝜳𝑯𝜳 𝑬𝜳 = 𝜳𝜳 𝑬𝟎 = 𝐦𝐢𝐧 𝑬 𝜳 𝜳 电子结构计算方法 电子结构的基本近似 泡利不相容原理：𝑥𝑖 包含空间坐标和自旋指标 Ψ 𝒙𝟏 , 𝒙𝟐 , … , 𝒙𝒊 , … , 𝒙𝒋 , … , 𝒙𝑵 = −Ψ 𝒙𝟏 , 𝒙𝟐 , … , 𝒙𝒋 , … , 𝒙𝒊 , … , 𝒙𝑵 • Slater行列式𝜳𝑺𝑫 ：最简单的满足交换反对称性的N-电子波函数 Ψ 𝑆𝐷 𝒙𝟏 , 𝒙𝟐 , … , 𝒙𝒊 , … , 𝒙𝒋 , … , 𝒙𝑵 𝜓1 𝒙1 1 𝜓1 𝒙2 = ⋮ 𝑁! 𝜓1 𝒙𝑁 𝜓2 𝒙1 𝜓2 𝒙2 ⋮ 𝜓2 𝒙𝑁 ⋯ 𝜓𝑁 𝒙1 ⋯ 𝜓𝑁 𝒙2 ⋱ ⋮ ⋯ 𝜓𝑁 𝒙𝑁 Fock建议将𝜳𝑯𝑷 的自洽场过程推广到𝜳𝑺𝑫 ，即Hartree-Fock 方法 𝐸0 = min 𝐸 Ψ = min Ψ 𝑆𝐷 𝐻 Ψ 𝑆𝐷 Ψ Ψ 电子结构计算方法 电子结构的基本近似 Hartree-Fock 方法 𝐸0 = min 𝐸 𝜓𝑖 Ψ = min Ψ 𝑆𝐷 𝐻 Ψ 𝑆𝐷 𝜓𝑖 𝜓𝑗 = 𝛿𝑖𝑗 Ψ Hartree-Fock 方程 𝑴 𝝍𝒋∗ 𝒓𝒋 𝝍𝒋 𝒓𝒋 𝝍𝒋 𝒙𝒊 𝝍𝒋∗ 𝒙𝒋 𝟏 𝟐 𝒁𝑰 − 𝜵𝒊 − + න 𝒅𝒓𝒋 𝝍𝒊 𝒙𝒊 − න 𝝍𝒊 𝒙𝒋 𝒅𝒓𝒋 = 𝜺𝒊 𝝍𝒊 𝒙𝒊 𝟐 𝒓𝒊 − 𝑹𝑰 𝒓𝒊 − 𝒓𝒋 𝒓𝒊 − 𝒓𝒋 𝑰 动能 外势𝒗𝒆𝒙𝒕 𝒓𝒊 𝒋≠𝒊 𝒋≠𝒊 库仑势𝒗𝑯 𝒓𝒊 库仑积分 𝑱 交换势 交换积分 𝑲 电子结构计算方法 Hartree-Fock 方程求解 自洽场(self-consistent field, SCF)方法 1 − 𝛻𝑖2 + 𝑉𝑒𝑥𝑡 𝒓𝑖 + 𝑉𝐻 𝒓𝑖 2 𝜓𝑖 𝑥𝑖 − න 𝑗≠𝑖 𝜓𝑗 𝑥𝑖 𝜓𝑗∗ 𝑥𝑗 𝑟𝑖 − 𝑟𝑗 𝜓𝑖 𝑥𝑗 𝑑𝑟𝑗 = 𝜀𝑖 𝜓𝑖 𝑥𝑖 1. 首先猜测一套波函数 𝜓0 ，作为输入 2. 由 𝜓0 计算得到库仑积分J 和交换积分K 3. 求解HF方程，得到新的波函数 𝜓 ′ 4. 重复步骤2和3，直到体系的总能量降至最低 电子结构计算方法 Hartree-Fock 方程求解 波函数的构造：LCAO (Linear Combination of Atomic Orbitals) 𝑁𝐵𝐹 𝜓𝑖 𝒓 = 𝑐𝑖𝛼 𝜙𝛼 𝒓 波函数/分子轨道 分子体系的电子波函数 𝛼=1 基函数/原子轨道 单个原子的电子波函数 LCAO-MO方法：分子轨道由原子轨道线性组合而成，表示原子 轨道的函数称为基函数，基函数组成的集合称为基组。适用于 由任意数量原子组成的任意分子。 ✓ 分子轨道的个数等于分子内所有原子轨道的总数； ✓ 每个分子轨道包括每个原子的每个轨道的贡献。 电子结构计算方法 Basis Set：基组 基组的选取：要适合计算机有效地计算 𝑁𝐵𝐹 𝜓𝑖 𝒓 = 𝑐𝑖𝛼 𝜙𝛼 𝒓 𝛼=1 基函数/原子轨道 • 氢原子轨道： 𝜓𝑛𝑙𝑚 𝑟, 𝜃, 𝜑 = 𝑅𝑛𝑙 𝑟 𝑌𝑙𝑚 𝜃, 𝜑 • Slater-type orbital (STO)：类氢原子轨道 𝜁 𝑆𝑛𝑙𝑚 𝑟, 𝜃, 𝜑 = 𝑁 × 𝑟 𝑛−1 𝑒 −𝜁𝑟 𝑌𝑙𝑚 𝜃, 𝜑 • Gaussian-type orbital (GTO)：更适合用来计算双电子积分 𝐺𝑖𝑗𝑘,𝛼 𝑟 2 𝑖 𝑗 𝑘 −𝛼𝑟 =𝑥 𝑦 𝑧 𝑒 𝑖+𝑗+𝑘 =𝑙 电子结构计算方法 Basis Set：基组 STO & GTO • STO 可以很好地描述电子的分布，但计算双电子积分很复杂 • 靠近原子核，cusp • 远离原子核，decay 1个STO 𝝓𝒊 = 𝒂𝟏 𝑮𝟏 + 𝒂𝟐 𝑮𝟐 + 𝒂𝟑 𝑮𝟑 1个CGs • Contracted Gaussians (CGs)： • GTO 计算双电子积分相对比较高效 • 用n 个GTO 的线性组合来逼近STO • 每个GTO的系数和展宽是固定的 电子结构计算方法 Gaussian 基组 one CGs per AO → one CGs per ‘core’ AO + more than one CGs per valence AO 基组：STO-3G、3-21G、6-31G、6-31G(d)、6-31G(d,p)、6-31+G(d)、6-31++G(d)、def2-XVY… 含义 计算特征 基组类型 STO-3G 每个轨道(内层和外层)由1个STO描述； 每个STO由3个GTOs线性组合表示 计算量小，计算精度不高 minimal basis 3-21G 每个内层轨道由1个STO描述，每个STO由3个GTOs表示； 每个价层轨道由2个STO描述，2个STO分别由2个和1个GTOs 表示 兼顾计算效率和计算精度 (使轨道在径向上的分布 更灵活) split-valence basis 6-31G 每个内层轨道由6个GTOs表示， 每个价层轨道劈裂成两个基函数，分别由3个和1个GTOs表示 同上 split-valence basis 6-31G(d)/ 6-31G* 在价层劈裂基组的基础上，为重原子添加1个d轨道 计算精度明显提高 (使轨道在角度上的分布 更灵活) polarization function 6-31G(d,p)/ 6-31G** 在价层劈裂基组的基础上，为重原子添加1个d轨道，也为氢原 子添加1个p轨道 同上 polarization function 6-31+G(d) 在价层劈裂&极化函数的基础上，为重原子添加1个弥散函数 计算量大，计算精度较难 评估 diffuse function 6-31++G(d) 在价层劈裂&极化函数的基础上，为重原子添加1个弥散函数， 也为氢原子添加1个弥散函数 同上 diffuse function 电子结构计算方法 The Roothaan Procedure 将体系整体的波函数表示为单电子波函数的组合（满足交换反对称）， 使用自洽场方法求解Hartree-Fock方程，得出一套单电子波函数的最优解， 将分子轨道表示为分子内所有原子轨道的线性组合（LCAO-MO）， 为了方法计算双电子积分，通常将每个原子轨道用一套Gaussian函数表示。 电子结构计算方法 Hartree-Fock 方法存在的问题 Slater行列式𝜳𝑺𝑫 ：最简单的满足交换反对称性的N-电子波函数 Ψ 𝑆𝐷 𝑥1 , 𝑥2 , … , 𝑥𝑖 , … , 𝑥𝑗 , … , 𝑥𝑁 𝜓1 𝒙1 1 𝜓1 𝒙2 = ⋮ 𝑁! 𝜓1 𝒙𝑁 𝜓2 𝒙1 𝜓2 𝒙2 ⋮ 𝜓2 𝒙𝑁 ⋯ 𝜓𝑁 𝒙1 ⋯ 𝜓𝑁 𝒙2 ⋱ ⋮ ⋯ 𝜓𝑁 𝒙𝑁 Fock建议将𝜳𝑯𝑷 的自洽场过程推广到𝜳𝑺𝑫 ，即Hartree-Fock 方法 • 假定了一个电子在由原子核和其它电子形成的平均势场中运动，没有 考虑电子间的瞬时相关作用，即电子关联（correlation）； • HF理论中采用的Slater行列式波函数满足泡利不相容原理，所以一定 程度上包含了同自旋电子的关联作用，这部分被单独命名为交换相互 作用（exchange）。 电子关联在考虑化学反应和电子激发等问题时十分重要！ 电子结构计算方法 Post-HF方法 组态相互作用（Configuration Interaction, CI） 通过引入激发组态来考虑电子关联 耦合簇方法（Coupled-cluster method, CC） 通过指数形式的耦合算符改进波函数 多体微扰理论（Moller-Plesset perturbation theory, MPx） 以HF波函数为零阶波函数的微扰方法 密度泛函理论（Density functional theory, DFT） … Errors: 𝑬𝑯𝑭 > 𝑬𝑴𝑷𝟐 > 𝑬𝑩𝟑𝑳𝒀𝑷 > 𝑬𝑪𝑰 > 𝑬𝑪𝑪 CPU time: 𝒕𝑪𝑪 > 𝒕𝑪𝑰 > 𝒕𝑴𝑷𝟐 > 𝒕𝑩𝟑𝑳𝒀𝑷 > 𝒕𝑯𝑭 计算软件系列一： Gaussian & GaussView Gaussian Gaussian是一款做半经验计算和从头计算的量子化学软件，使用非常广泛 主要功能： 软件信息： – 分子结构和能量 • 官网：http://www.gaussian.com – 过渡态的能量和结构 • License：商业软件 – 化学键以及反应能量 • 当前版本：Gaussian 16 C.01 – 振动频率 • 运行平台：Windows, Unix/Linux, Mac OS – 红外和拉曼光谱 – 核磁共振谱 关于授权： – 热力学性质 超算中心已统一购买并提供给校内师生使用 – 反应路径 学校主页→信息门户→正版软件→高斯 使用方法：从头计算法、半经验算法、密度泛函理论、分子力学、杂化方法 GaussView GaussView是专门针对Gaussian软件开发的可视化图形设置界面 主要功能： 软件信息： – 建模：导入、构建分子结构 • 官网：http://gaussian.com/gaussview6/ – 设置Gaussian输入文件 • License：商业软件 – 检索、查看计算数据 • 当前版本：GaussView 6 – 可视化计算结果 • 运行平台：Windows, Unix/Linux, Mac OS Gaussian & GaussView GaussView 建模、输入、可视化 输入/输出 文本文件 Gaussian 计算、输出 申请使用权限 1. 下载申请表：学校主页→信息门户→正版软件→高斯 2. 填表并发送至邮箱：sccadmin@ustc.edu.cn Windows版本：将获取许可证号，之后在学校主页下载、安装软件 Linux版本：不提供官方下载，超算系统上已安装，申请之后相应账户将获取软件使用权限 开始计算 一个例子：水分子结构优化计算 Gaussian输入文件：H2O.gjf %nprocshared=1 CPU核数 %mem=500MB 内存 %chk=C:\Users\SCC-Liu\Desktop\H2O.chk checkpoint文件(续算、分析) #p b3lyp/6-31g sp 关键字 空格（必须的） single point energy of water 标题 01 空格（必须的） 电荷和自旋多重态 O -1.33828988 0.39033457 0.00000000 H -0.37828988 0.39033457 0.00000000 H -1.65874447 1.29527040 0.00000000 元素符号和原子坐标 两个空格（必须的） 开始计算 Windows 版本运行 开始计算 以瀚海20系统为例： $ module load gaussian/g16.C01 $ sbatch gauss_job.sh gauss_job.sh：slurm作业脚本文件 Linux 版本运行 #!/bin/sh #SBATCH -p test #SBATCH - -qos=testqos #SBATCH -o job-%j.log #SBATCH -e job-%j.err #SBATCH -N 1 -n 1 echo Time is `date` echo Directory is $PWD echo This job runs on the following nodes: echo $SLURM_JOB_NODELIST echo This job has allocated $SLURM_JOB_CPUS_PER_NODE cpu cores. g16 H2O Gaussian计算 典型的计算流程 1. 准备初始构型 2. 编写输入文件 3. 开始计算 4. 分析计算结果 Gaussian计算 准备初始构型：GaussView 周期表元素：不同成键类型 有机链：… 有机环：… Gaussian计算 准备输入文件：GaussView 作业类型：单点能、结构优化、频率… 计算方法：HF、半经验、DFT… 基组：STO-3G、3-21G、6-31G… Gaussian计算 准备输入文件：GaussView %nprocshared=1 %mem=500MB %chk=C:\Users\SCC-Liu\Desktop\H2O.chk #p b3lyp/6-31g sp single point energy of water Gaussian输入文件：H2O.gjf 01 O -1.33828988 0.39033457 0.00000000 H -0.37828988 0.39033457 0.00000000 H -1.65874447 1.29527040 0.00000000 Gaussian计算 准备输入文件：输入文件格式 %nprocshared=1 1. Link 0命令段：以%开始，内容包括 CPU与内存、chk文件路径和名称等 %mem=500MB %chk=C:\Users\SCC-Liu\Desktop\H2O.chk 2. Route Section：以#开始，内容包括 输出方式、计算方法、基组等 #p b3lyp/6-31g sp 空行：Route section必须以空行结束 single point energy of water 3. Title Section：计算任务的简洁描述 空行：Title section必须以空行结束 01 O -1.33828988 0.39033457 0.00000000 H -0.37828988 0.39033457 0.00000000 H -1.65874447 1.29527040 0.00000000 4. Molecule Specification：分子结构， 内容包括电荷和自旋多重度、原子坐标 空行：Molecule Specification必须以空行结束 5. 其他可选项：通常以空行结束 Gaussian计算 开始计算：Gaussian Windows 版本运行… Linux 版本运行… Gaussian计算 其他：频率、光谱 分析计算结果：GaussView 分 子 轨 道 电 荷 密 度 静 电 势 谢谢！