计算化学/几何优化

几何优化

需要注意的重要特点是，势能函数在内侧更陡峭，在外侧更平缓，尤其是在向解离方向移动时。这个概念可以扩展到多原子分子，在多原子分子中，我们拥有 $n$ 维度的势能面。 ( $n$ 是 3N-6，其中 N 是原子数。-6 来自分子的平移和旋转，即 MOPAC 打印输出中提到的平凡振动。）势能井的底部位于平衡键长处。在这个区域，势能函数类似于抛物线 $E={\frac {1}{2}}fr^{2}+c$

如果我们可以计算总能量，那么双原子分子的几何优化就变得微不足道了，我们可以像牛顿方法那样，通过改变 $r$ 的值来向下跟踪势能函数，直到我们找到最小值两侧的点。从这个位置开始，我们可以不断二分间隔，并使用二次插值，直到达到所需的精度。

这项程序有了很大的改进，因为我们能够对哈密顿量相对于 3N 个核位置坐标 ( $x,y,$ 和 $z$ 在每个原子中心上) 进行微分。所有现代量子化学程序现在都使用这种技术。类似牛顿-拉弗森的方法被用来沿着势能面下降，并且记住旧的梯度来改进优化。一旦进入二次区域，就能快速收敛到任何精度。

笛卡尔坐标和Z矩阵

分子几何结构的枯燥问题。

到目前为止，我们一直在讨论几何结构，假设它们用各个原子的 $x,y,z$ 坐标矩阵表示。有时，更实用且更高效的方法是用分子的键长和键角来表示相同的独特几何排列。这种方法已经消除了对应于分子旋转和平移的 6 个自由度。

创建分子几何结构的实用方法

通过观察来写出中等大小分子的 $z$ -矩阵是一个可能的但艰巨的任务。在没有计算机辅助的情况下，创建笛卡尔坐标是不可能的。在实践中，人们使用分子建模系统，例如在 UNIX 系统上可用的 MACROMODEL（哥伦比亚大学，纽约）。分子建模系统使用几种标准力场之一来优化构建的几何结构。（注意，对于经典计算，如 MM2，术语为力场；对于量子方法，如 AM1（Austin Method-1），术语为哈密顿量。）

力场计算建立了分子模型系统中所有原子的连接性，因此生成矩阵涉及以最具化学意义的方式遍历连接性的树，选取叶节点（氢原子），以及它们的二面角作为 $z$ -矩阵条目。在计算术语中，树向下生长而不是向上生长。树的节点通常是原子。终端节点称为叶节点。

n-戊烷的树状形式

                                        C
                                      . . . .
                                    .   .   .   .
                                  C     H    H   H
                                . . .
                              .   .   .
                             C    H    H
                           . . .
                         .   .   .
                        C     H    H
                      . . .
                    .   .   .
                   C    H    H
                 . . .
               .   .   .
              H    H    H

化学连接性的理论研究称为化学图论，它在计算机辅助化学合成和文献检索等领域非常重要。

极坐标

                          z
                          |
                          |
                          |               * (atom)i
                          | theta(i)    . .
                          |           .   .
                          |))       .     .
                          |  )    . r(i)  .
                          |   ) .         .
                          |   .           .
                          | .             .
                          0------------------------y
                         /   .            .
                        /     ( .         .
                       /(   ((     .      .
                      /  (((          .   .
                     /   phi(i)          ..
                    /                     .
                   /
                  /
                 /
                /
              x

请注意，轴系具有手性，遵循右手定则。 $z$ -轴（主轴）是拇指，食指指向 $x$ ，中指指向 $y$ 。强烈建议遵守此惯例，因为虽然左手坐标系在保持一致的情况下可以使用，但立体化学可能导致巨大的问题（蛋白质仅由L型氨基酸构成）。

二氧化碳

                      :z
                      :
                      :
                      :
       00             CC           00
      0  0 ----------C  ----------0  0  ----------x
      0  0 ----------C  ----------0  0
       00 1           CC1          00 2
                    /
                   /
                  /
                 /
                y

通过观察，可以写出笛卡尔坐标

	x	y	z
C	0	0	0
O1	r	0	0
O2	-r	0	0

对称点群为 $D_{\infty h}$ ，但在量子化学计算中，算法通常要求仅使用阿贝尔群（~非简并~），因此 $D_{2}h$ 用于二氧化碳。（ $D_{2}h$ 是从头算量子化学家最喜欢的群，因为它有 8 个不可约表示，使得计算比更常见的 $C_{s}$ 或 $C_{2}v$ 群小得多，例如只有两个唯一重原子的醌。旋转主轴几乎总是定义为 $z$ -轴。

氯甲烷

在氯甲烷中，C-Cl 键长为 1.8 埃，C-H 键长为 1.094 埃，CL-C-H 角为 110.6667 度。在图中，Z 轴沿 C-CL 键。

                       CL
                        :
                        :
                        :
                        :
                        :
                        C --------------------y
                       /=.
                      / = .
                    //  =  ..            (x out of page)
                   //   =   . .
                 / /    =    . .
                / /    H(3)   . .
               /-/             .-.
              H(1)             H(2)

将 C 原子作为原点是显而易见的。Cl 将位于 $z$ -轴上方的 1.8 处。所有氢的 $z$ -坐标为 $1.094~*~\cos ~110.67$ ，即 $-~(1.094~*~\cos ~69.33)$ 。

H(3) 的 $y$ 坐标为零， $x$ 坐标为负（即在右手笛卡尔框架中，参见极坐标图）。

最好始终使用传统的右手坐标系。记忆技巧是想象一张坐标纸， $x$ 像往常一样横向，而 $y$ 纵向。正的 z 则从纸面向外。有趣的是，如果弄错这一点，不会影响除了旋光性等之外的任何计算结果，但更严重的问题是当坐标用作构建块时。例如，如果您将来自计算的氨基酸坐标集与来自晶体结构或分子模型构建器的坐标集组合，您可能会无意中在同一个肽中混合 D 和 L 型。

因此，上面负的 $x$ 值将是 $1.094~*~\sin ~69.33$ 。现在记住这个值作为 ' $xc$ '，我们可以轻松地得到 H2 和 H3 的 $x$ 和 $y$ 坐标。它们通过 120 度旋转相关联，并且除了 $y$ 方向的符号变化之外，具有相同的坐标。现在， $y$ 坐标将是 + 或 - $xc~*~\sin ~120$ 。 $x$ 坐标是 cos~120。

记住常用三角函数值的实用表格


0	0	90
30	${\frac {1}{2}}$	60
45	${\frac {1}{\sqrt {2}}}$	45
60	${\frac {\sqrt {3}}{2}}$	30
90	1	0

正弦值从表格向下读取，余弦值从表格向上读取。

使用 MACROMODEL 建模系统，人们可能会押注 CH4 的图标，并使用 'Cl' 图标推动其中一个氢原子。然后可以使用 MM2 力场优化几何结构。

苯

通过将原点放置在对称中心，苯的几何结构计算得到了极大的简化。碳原子距离对称中心 1.39 埃，氢原子距离对称中心 2.47 埃。两个相邻碳原子和对称中心形成的三角形是等边三角形。可以使用 60 度的余弦和正弦计算 H2 和 C2 的坐标。所有其他坐标都通过符号变化相关联。

        Y                                      H1
        :                                      *
        :                                      *
        :                         H6           *            H2
        :                          *          1*           *
        *                            *       *   *1.39   *
      * : *                            *   *       *   *
    *   :   *                            *6         2*
  *     :     *                          *           *
  *     :.....*.........X                *           *
  *           *                          *           *
  *           *                          *5         3*
    *       *                          *   *       *   *
      *   *                          *       * 4 *       *
        *                          *           *           *
                                  H5           *            H3
                                               *
                                               *
                                               H4

                 (Z-axis vertically out of page)
 Benzene
 Centre          0.000    0.000    0.000
 C               0.000    1.390    0.000
 H               0.000    2.470    0.000
 C               1.204    0.695    0.000
 H               2.139    1.235    0.000
 C               0.000   -1.390    0.000
 H               0.000   -2.470    0.000
 C              -1.204   -0.695    0.000
 H              -2.139   -1.235    0.000
 C               1.204   -0.695    0.000
 H               2.139   -1.235    0.000
 C              -1.204    0.695    0.000
 H              -2.139    1.235    0.000
 **

我用计算器计算出上述坐标的一个问题是精度。轴上的原子具有无限精度，（~零 | ~），而计算出的原子位置只有 4 位数字，（原子 2、3、5 和 6）。这表现为从完整的 $D_{6}h$ 对称性中轻微的滑动，这会导致简并轨道和相等的原子性质不完全相等。

尽管有限精度对输入数据和计算的影响通常并不严重，但在计算结果中经常会发生这种意想不到的情况。在量子化学中，我们需要比平时更高的精度，这可以通过双精度算术来获得，因为我们关心的是包含大量核心电子能量的实体之间的能量差异，这些能量会“消耗”总能量的左端有效数字。（与高能现象相比，大多数化学都与能量的微小差异有关。）如果您正在计算包含苯基的分子坐标，您可能会使用分子建模系统，但在某些情况下，您会从上述高度对称的坐标开始，因为存在高度对称性。

投影

Newman Projection                             Sawhorse Projection
*****************                             *******************
                                                                   .H6
                                                                  .
               *)                                                .
               *  )                                             * C2
               *    )                                         .. .
               * 27  )   @                         H1      . .   .
  @           .* .    )@                           .    .  .      .
    @      .   *    .@                             .   H4.        .
      @        *)                                  .   .           H5
        @.     *  )   .                            . .
               *120                             C1 *
        .      *   )   .                         .   .
             *   *)                            .       .
         . *       *  .              Z       .           .
         *           *              /       H2            H3
      *    .        .  *           /
    *         .  .       *        /
  *          @             *     /
                                :------X      Z       Y
            @                   :        (    ^      /         )
                                :        (    ^    /           )
           @                    :        (    ^ /              )
                                Y        (    *---------X      )
                                         (                     )

Newman 投影对于定量处理二面角最为有用。存在一些程序可以从晶体结构或 MACROMODEL 坐标中打印出二面角。

环系中的局部原点

通常希望描述一个子系统并将其转移到几何位置。这个旋转/平移过程可以用一个向量和三个欧拉角来描述。模型片段非常合适的取向是使任何芳香部分尽可能地平面化（可能通过最小二乘拟合）。笛卡尔坐标系的原点可以设置为由质心定义的环的中心（假设所有环原子都是碳，并且不考虑其他原子）。如果考虑到 C₆H₄I- 基团，原因就很明显了。质心将位于 C-I 键中，而环系统的中心是一个更实用的原点。

有很多方法可以定义欧拉角。QM 使用三个角度 $\alpha ,\beta$ 和 $\gamma$ ，它们是绕着 $z$ ， $y^{'}$ ，然后是 $z^{'}$ 的旋转（即 $z$ ，新的 $y$ ，然后是新的 $z$ 轴。

$\cos(\alpha )$	$\sin(\alpha )$	0
$-\sin(\alpha )$	$\cos(\alpha )$	0
0	0	1

那么

$\cos(\beta )$	0	$\sin(\beta )$
0	1	0
$-\sin(\beta )$	0	$\cos(\beta )$

那么

$\cos(\gamma )$	$\sin(\gamma )$	0
$-\sin(\gamma )$	$\cos(\gamma )$	0
0	0	1

以下是来自这三个旋转的乘积所得到的变换矩阵的 FORTRAN 代码片段。

  tran (1,1) = (cos (alpha) * cos (beta) * cos (gamma)) _
  - sin (alpha) * sin (gamma)
  tran (2,1) = (cos (alpha) * cos (beta) * sin (gamma)) _
  + (sin (alpha) * cos (gamma))
  tran (3,1) = - cos (alpha) * sin (beta)
  tran (1,2) = - ( (sin (alpha) * cos (beta) * cos (gamma)) _
  +  ( cos (alpha) * sin (gamma))  )
  tran (2,2) = - ( sin (alpha) * cos (beta) * sin (gamma) ) _
  + ( cos (alpha) * cos (gamma) )
  tran (3,2) = sin (alpha) * sin (beta)
  tran (1,3) = sin (beta) * cos (gamma)
  tran (2,3) = sin (beta) * sin (gamma)
  tran (3,3) = cos (beta)

优化

一旦你有了你的试验几何结构，通过图形程序或三角学获得，正如之前描述的那样，强大的量子或分子力学程序之一将通过最小化 (3N-6) 个自由坐标来解决几何结构。

在许多情况下，这不会给出唯一的几何结构，而仅仅是许多可能的构象之一。

下一章 - 分子量子力学的应用