1、结构优化 Wbfeng 势能面的方程 分子的完全Schr dinger方程 Born Oppenheimer近似后方程分解为核运动和电子运动两个方程 几何结构优化问题的数学描述 势能面的势能函数 势能函数求极值问题 分子几何结构优化的数学过程 1 早期优化方法 逐点优化法 基于能量本身 计算量大 收敛慢 不利于程序化2 现代优化方法 能量梯度法 基于能量的一阶 二阶导数 更准确快速 易于程序化 多维优化方法 一阶导数法 最陡下降法共轭梯度 方向 法二阶导数法 Newton Raphson方法准Newton方法 最陡下降法 基本原理 从指定点出发 循梯度的负方向搜寻到极值点后作为新的起点 进行下
2、一步搜寻例 函数从 9 9 出发 在 18 36 方向找到极值点 4 1 后 a 求该点的负梯度方向 8 4 得到下一点为 4 3 b 得到该方向的方程y 0 5x 1 c 继续使用Lagrange乘因子法 求得该方向极值点 2 3 2 3 重复上述步骤 得到 0 296 0 074 优点 对于远离驻点的结构 优化效率非常高 能很快释放分子内的力缺点 每一步都要进行直角转向 收敛慢 校正过度 振荡 共轭梯度 方向 法 基本原理 做完一次线性搜索后 后一次优化的方向取该点的梯度方向与前一次优化的方向的组合 例 使用共轭梯度法求的极值点1 起始点 9 9 的负梯度为 18 36 此方向极值点 4
3、1 2 点 4 1 处的负梯度为 8 4 搜索方向表达式为y 1 4x 可以找到极值点 0 0 优点 对于有M个变量的函数 可以通过M步优化找到极值 两种共轭梯度法 1 纯的二次函数Fletcher Reeves算法2 非纯二次函数Polak Ribiere算法 Newton Raphson方法 将势能函数展开成Taylor级数 如果势能函数是纯二次函数 那么存在条件在势能极小点处 对函数求导 可以得到 例 求的极值点一阶导数Hessian矩阵 优点 对于纯二次函数 可以一步找到极值点 缺点 要求Hessian必须正定 否则将得到能量更高的坐标 对于非纯二次函数 需要多步计算 Hessian矩
4、阵的计算量和存储量都非常大 主要适用于小分子体系 准Newton方法 基本原理 初猜一个Hessian矩阵 开始优化后 每步更新一次Hessian矩阵 每次更新Hessian矩阵都只使用上一步的Hessian矩阵和当点的一阶导数 DFP法 BFGS法 优化算法的选择 算法的选择由多种因素决定 1 大分子体系多使用最陡下降法或者共轭梯度法 2 小分子多用Newton Raphson法 3 对于远离驻点的结构 结合最陡下降法和Newton Raphson法 收敛判据 1 力最大力 均方根力2 位移最大位移 均方根位移 PRHF 6 31G d OptTestEthyleneGeometryOpti
5、mization01CC1CCH1CH2HCCH1CH2HCC3180 H2CH1HCC3180 H2CH1HCC4180 Variables CC 1 31CH 1 07HCC 121 5 输入文件 优化作业执行过程 1 18 20 38 1 1 3 2 9 110 17 6 18 5 40 1 2 3 5 1 6 6 7 1 11 1 16 1 25 1 30 1 1 2 3 4 7 1 1 5 5 2 38 5 2 6 7 2 8 2 9 2 10 2 28 1 1 7 1 2 3 16 1 18 20 3 1 99 99 2 9 110 2 3 5 1 6 6 7 1 11 1 16
6、1 25 1 30 1 1 2 3 4 5 5 7 1 16 3 1 5 5 2 38 5 2 7 1 2 3 16 1 18 20 3 2 9 110 2 6 7 2 8 2 9 2 10 2 19 2 28 1 1 99 9 1 99 L101L103L202L301L302L303L401L502L601L701L702L703L716L9999 优化开始的标志 第一次做L103初猜Hessian矩阵 优化结束的标志 Hessian矩阵 冗余内坐标 程序自动生成跟据共价半径确认成键 包括氢键和分子间键 所有键原子之间的键角所有键原子间的二面角 Gaussiane3 08 优化时的坐标选项
7、 OPT RedundantOPT AddRedundantOPT Z matrixOPT Cartesian部分冻结优化 HF 6 31G Opt ModRedundantTestPartialoptimization11lCH1R1H1R12A1O1R22A23120 0H4R33A32180 0 A1 120 0 R3 1 1451 3F5432F 优化收敛问题 a 默认的优化次数不够 增加优化次数opt RestartMaxcycle N b 初猜力常数与实际不符 解决办法 1 从振动分析计算的chk文件中读入力常数 Opt ReadFc2 使用给定的方法和基组计算力常数 Opt Ca
8、lcFc3 优化中的每一步都计算力常数 Opt CalcAllc 从能量最低的最优结构重起优化作业Opt CalcFc Geom Check Step n 其他优化选项 NoEigenTest 在使用Berny优化计算过渡态时 不做曲率测试NoFreeze 激活所有冻结的变量 Geom Check Tight VeryTight LooseIOP 1 8 m 改变步长0 01 m 优化结果的确认 计算完整的Hessian矩阵检查Hessian矩阵的负本征值极小点0个负本征值过渡态有且仅有1个负本征值 寻找极小点如果有负本征值 循负本征值方向得到能量更低的结构 寻找过渡态如果没有负本征值 沿最小
9、本征值 上山 通过振动分析计算来确认优化结果 复合作业 例 HF 6 31goptfreq振动分析与优化计算必须使用同样的方法和基组其他类型的复合作业 HF 6 31G HF 6 31G d p 振动分析计算 内坐标与简正坐标下Hessian矩阵的转换 内坐标下的势能表示 平衡点的势能规定为零 核运动的哈密顿 交叉项存在 无法分离变量求解 简正坐标和简正振动模式 核运动的哈密顿质权坐标再找到一个U阵将K对角化 可以得到核运动的哈密顿为 通解为据此能求算相关的热化学性质 例 1595 1826 cm 13756 4188 cm 13652 4070 cm 11 1451 1151 114HF 6 31G 实验 练习 H3O 分子结构的优化NH4 分子的结构优化
