24小时热门版块排行榜

返回列表

本帖产生 1 个 QC强帖，点击这里进行查看

zhou2009

版主 (著名写手)

QC强帖: 15
应助: 337 (大学生)
贵宾: 4.053
金币: 42654.6
红花: 153
帖子: 1729
在线: 1274.7小时
虫号: 676800
注册: 2008-12-18
性别: GG
专业: 理论和计算化学
管辖: 量子化学

[交流] 【zhou2009个人文集】分子片ρ计算的一个值得注意的问题已有2人参与

分子片ρ计算的一个值得注意的问题

Zhou2009

前帖“关于卟吩-镁的量化计算及图形解析” （https://muchong.com/bbs/viewthread.php?tid=2934808&page=1#pid362935）曾经提到，作电子密度差Δρ，从电荷的净变化可以考察化学体系相互作用。
对于小分子，可以作分子与组成分子的各孤立原子的Δρ。但是这里常常有一个孤立原子参与成键波函数的取向难以与分子中成键取向一致的问题，特别是当分子中的原子数目很多时，与各孤立原子作的Δρ将变得极其复杂、繁琐，不便观察解析，从而失去了Δρ的直观实用意义。
实际上，分子或化学体系总是可以划分成基团或分子片的，考察它们之间的作用常是关键、才是抓住了主要矛盾。象涉及外界小分子a对基体大分子b的影响作用，用总体分子体系ab与各分子片a、b作密度差即可，这时可以看到孤立的分子片a与孤立的分子片b在发生化学作用后，电荷在a、b上发生了什么样的变化，而相互作用会完全表现在这电荷变化上。
由于是分子片，就不存在孤立原子参与成键波函数的取向难以与分子中成键取向一致的问题了。即使是拉开一个键所形成的两个分子片，分子片也会将拉开处的波函数取向与分子中成键取向形成一致。而且分子片内部的电荷变化可以基本抵消，只表现分子片之间的作用。
分子片的Δρ，适合于氢键、吸附、催化、表面作用等，也可用于断开一个键看两个分子片之间的作用，甚至可以考察同一个MO在受影响的情况下能级发生变化前后电荷分布的净变化。

当涉及到很大的主体分子、分子系综表面时，考察较小分子与之作用，当作Δρ时，需要取能包含主体大分子与小分子的分子体系的cube空间，求总体的ρab。也需要取同样大小的cube空间来分别算这作为分子片的小分子的ρa与大分子的ρb。Δρ=ρab-ρa-ρb。
计算采用b97d/6-31+g(d)。但因为总体的cube空间很大，如果真要计算电荷量，会出现所取格点可能精度不够的问题。这个问题又分两个层次：
首先、对于ρ的计算，它是对cube空间中每一个格点进行的电子密度计算，它是够精确的，在作Δρ时，因为是对应的点对点进行的相减，也是够精确的。并不存在精度不够的问题。
其次、当我们真要用这个ρ来进行全空间或部分空间的电荷加和时，却会出现所取格点可能精度不够的问题。

以优化后的甲基咪唑-镁卟吩（ab）划分为小分子片甲基咪唑a与大分子片镁卟吩b为例，当进入multiwfn计算ρ时，计算ρab时，指定网格点数，令参数由程序自动选择，而计算ρa、ρb则沿用这套数据手动输入，这样就有统一的网格点数（100,100,100）、空间大小（cube的起始点为：-15.306760,-15.602686,-8.55286）、点的间隔（0.310913,0.312060,0.248017）。这里的长度为原子单位，即bohr。

程序会自动对全空间的电子进行加和：甲基咪唑本应有44个电子，但计算得到的加和却为46.32，镁卟吩本应有172个电子，但计算得到的加和却为170.25。计算的正负误差竟达到约2个电子！这种误差在对全空间或部分空间的进行电荷加和、进一步对电荷定量分析时显然是不允许的。

进一步查找误差的原因，不是程序的问题，是cube参数选择时出的问题。Cube主要参数通常表示成为一个矩阵：
X0, Y0, Z0    （网格原点）
N1, X1, Y1, Z1    （X-方向的点数和点的间隔）
N2, X2, Y2, Z2    （Y-方向的点数和点的间隔）
N3, X3, Y3, Z3    （Z-方向的点数和点的间隔）

我们通常会先给出空间格点参数（N1,N2,N3），如上面的100,100,100。然后根据分子整体的空间大小选择或者由程序自动选择cube的起始点（X0,Y0,Z0），如上面列出的数据。这两组数据其实就决定了网格点的间隔，即实际的xyz方向上的长度（参数只列出了网格的起始点）分别除以相应的网格点数就得到相应的网点的间隔(X1,Y2,Z3)。从上面例子看，网点的间隔取得大了点是引起误差的原因。

在multiwfn计算中，对全空间的电荷计算，是这样进行的：ρ是几率密度，ρ的量纲是1/(bohr3) ，ρdv是几率，其意义是在dv这个空间体积微元中找到电子的几率。空间各个体积微元的几率加和就是电子数，∫ρdv，归一或归N（分子或分子片或MO的电子数）。
由于cube的格点已是离散化的，ρ仅是一个格点处的值，具体的体积微元dv却是这样一个空间X1×Y2×Z3，这就要求在这样一个空间中各处都仍具有ρ这个恒定的值而没有变化，然而在原子核附近和原子之间，ρ的变化还是很密集的，只有dv足够小，才能使得其中ρ的变化是可以忽略的。

可见，在在网格点数、空间大小、格点间隔这三组参数中，首先要确定的是格点的间隔，在甲基咪唑-镁卟吩这个例子中，格点间隔取0.2 bohr即可。其次空间大小是由分子体系的大小已经决定了的，取一个足够装得下所算分子体系即可。此例取cube的起始点为-15（终点15，总长30）、-15（终点15，总长30）、-10（终点20，总长30）。由于xyz三边的总长各为30，分别除以隔点间隔0.2即得到有足够精度的网格点数，此例为150,150,150。（当然分子体系不会是空间对等的，为充分利用空间，各数字不一定等同。）

经过这样的选择，甲基咪唑-镁卟吩算得电子数为146.038（实应为146），甲基咪唑为32.038（实应为32），卟吩镁为113.9999（实应为114）。

但这里还是要提请注意的是，上述分子片的取法与使用multiwfn的“自定义运算”功能实际上是一样的，而multiwfn步骤简化了许多。但这种取法与在multiwfn中从分子整体中取的分子片则有原则的不同，后者取的分子片是在分子整体中存在的，不是孤立分子片的，这时对分子片作的密度差，分子片本身的电荷被完全抵消了，求得的Δρ，仅是分子片之间波函数所形成的重叠布居电荷，这是另外一种信息、另外的考察方式。

[ Last edited by 375642546 on 2011-4-19 at 15:34 ]

回复此楼