談談NBO中用bndidx關鍵詞輸出的三種鍵級
注：完全沒必要用NBO算鍵級！Multiwfn（可在http://www.shanxitv.org/multiwfn免費下載）是所有程序中支持鍵級種類最全面的、最好用的。怎么用Multiwfn算各種鍵級看Multiwfn手冊4.9節的例子，各種鍵級原理介紹見Multiwfn手冊3.11節，對于Gaussian用戶只要提供fch文件即可，敲兩下鍵盤鍵級就輸出出來了。更多和鍵級有關的博文看《談談原子間是否成鍵的判斷問題》（http://www.shanxitv.org/414）、《通過鍵級曲線和ELF/LOL/RDG等值面動畫研究化學反應過程》（http://www.shanxitv.org/200）。
1. 前言
用了bndidx關鍵詞后，NBO程序會進行鍵級分析，輸出三種鍵級矩陣，包括：
1 Wiberg bond index matrix in the NAO basis
2 Atom-atom overlap-weighted NAO bond order
3 MO bond order
矩陣中(i,j)元素對應i、j原子間的鍵級。對這三種不同鍵級，有些人對應當選取哪個感到疑惑，此文將談談這三種鍵級的含義及其合理性。值得一提的是，經常有人說“NBO鍵級”，這種說法是很不嚴格的，在NBO程序中，上面已經看到有三種鍵級了，NBO還能產生LMO鍵級、NRT鍵級，另外，NBO程序搜索出的NBO們的占據數本身也可以用于衡量鍵級。所以，不應當使用含糊不清的“NBO鍵級”這個詞。
為了簡明且便于驗證，本文例子使用HF/STO-3G的基態水分子波函數。下面的是用DMNAO關鍵詞輸出的密度矩陣，將在后文用到。后文直接將此矩陣稱為DMNAO。
NAO 1 2 3 4 5 6 7
---------- ------- ------- ------- ------- ------- ------- -------
1. O 1 (S) 2.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
2. O 1 (S) 0.0000 1.7325 0.0000 0.0000 0.3538 0.4113 0.4113
3. O 1 (px) 0.0000 0.0000 2.0000 0.0000 0.0000 0.0000 0.0000
4. O 1 (py) 0.0000 0.0000 0.0000 1.1399 0.0000 0.7001 -0.7001
5. O 1 (pz) 0.0000 0.3538 0.0000 0.0000 1.5321 -0.5439 -0.5439
6. H 2 (S) 0.0000 0.4113 0.0000 0.7001 -0.5439 0.7977 -0.0623
7. H 3 (S) 0.0000 0.4113 0.0000 -0.7001 -0.5439 -0.0623 0.7977
2. 第一種鍵級(Wiberg bond index matrix in the NAO basis)
說起第一種鍵級，先要回顧一下Wiberg鍵級和Mayer鍵級的原始定義。Wiberg鍵級最早來自Tetrahedron 24,1083，鍵級W_AB=∑[i∈A]∑[j∈B]P(i,j)^2，其中P為密度矩陣，這個定義只適合正交基。Mayer在CPL 97,270里將其擴展到非正交基，A、B間的鍵級定義為∑[i∈A]∑[j∈B]PS(i,j)*PS(j,i)，這里PS(i,j)是指P和重疊矩陣S做矩陣乘所得矩陣的(i,j)元素。當使用正交基后，變為∑[i∈A]∑[j∈B]P(i,j)*P(j,i)，由于P是對稱的，所以又還原到了Wiberg鍵級。Mayer鍵級數值比較符合一般化學觀念，比如單重鍵約為1，三重鍵約為3，苯的C-C鍵約為1.5。
NBO的各種分析主要基于由原始非正交基函數轉化成的自然原子軌道(NAO)，NAO是正交基，所以可以直接使用Wiberg鍵級進行分析，也就是NBO輸出的Wiberg bond index matrix in the NAO basis。這是三種鍵級中唯一推薦的鍵級。NBO給出的結果：
Atom 1 2 3
---- ------ ------ ------
1. O 0.0000 0.9552 0.9552
2. H 0.9552 0.0000 0.0039
3. H 0.9552 0.0039 0.0000
可以很容易地通過NAO的密度矩陣DMNAO來計算此矩陣。此矩陣的(2,3)代表H-H的鍵級，由于6、7號NAO分別是兩個氫上的唯一基函數，鍵級也就是DMNAO(6,7)^2=(-0.0623)^2=0.00388。再比如(1,3)代表的O-H鍵級，由于1~5號基函數對應O，也就是∑[i=1,2,3,4,5]DMNAO(i,7)^2=0.95513。
在NBO輸出這個矩陣之后還會輸出Wiberg bond index, Totals by atom，指的就是每個原子形成鍵級的總和，也稱為原子價。對于第二種、第三種鍵級也是如此。
3. 第二種鍵級(Atom-atom overlap-weighted NAO bond order)
這種鍵級號稱是考慮了原子軌道間重疊，引入距離因素，并使得結果能看出成鍵反鍵，但實際結果不好，而且其定義我認為沒有任何根據。A、B間鍵級定義為∑[i∈A]∑[j∈B]P(i,j)*S'(i,j)，其中S'是PNAO重疊矩陣（注：PNAO是原始基函數向NAO轉換的中間基函數，PNAO在原子內正交而在原子間未必正交）。這種定義很莫名其妙，密度矩陣用NAO的，但重疊矩陣卻用PNAO的。從其形式上，可以想起Mulliken鍵級。Mulliken鍵級就是Mulliken布居分析得到的兩原子間交叉項的總和，∑[imo∈占據]∑[i∈A]∑[j∈B]2*占據數(imo)*C(i,imo)*C(j,imo)*S(i,j)，寫成密度矩陣的形式就是∑[i∈A]∑[j∈B]2*P(i,j)*S(i,j)。所以，第二種鍵級實際上是一種既不完全是NAO、也不完全是PNAO下的不倫不類的Mulliken鍵級，而且還少乘了2。這種古怪的鍵級缺乏依據，實際結果又差，所以不要使用。NBO輸出的結果為：
Atom 1 2 3
---- ------ ------ ------
1. O 0.0000 0.5399 0.5399
2. H 0.5399 0.0000 -0.0145
3. H 0.5399 -0.0145 0.0000
這也容易驗證，使用SPNAO關鍵詞可以輸出PNAO重疊矩陣S'
NAO 1 2 3 4 5 6 7
---------- ------- ------- ------- ------- ------- ------- -------
1. O 1 (S) 1.0000 0.0000 0.0000 0.0000 0.0000 0.0935 0.0935
2. O 1 (S) 0.0000 1.0000 0.0000 0.0000 0.0000 0.4675 0.4675
3. O 1 (px) 0.0000 0.0000 1.0000 0.0000 0.0000 0.0000 0.0000
4. O 1 (py) 0.0000 0.0000 0.0000 1.0000 0.0000 0.3205 -0.3205
5. O 1 (pz) 0.0000 0.0000 0.0000 0.0000 1.0000 -0.2265 -0.2265
6. H 2 (S) 0.0935 0.4675 0.0000 0.3205 -0.2265 1.0000 0.2328
7. H 3 (S) 0.0935 0.4675 0.0000 -0.3205 -0.2265 0.2328 1.0000
設T(i,j)=DMNAO(i,j)*S'(i,j)，則T為：
1 2 3 4 5 6 7
1 2.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
2 0.00000 1.73250 0.00000 0.00000 0.00000 0.19228 0.19228
3 0.00000 0.00000 2.00000 0.00000 0.00000 0.00000 0.00000
4 0.00000 0.00000 0.00000 1.13990 0.00000 0.22438 0.22438
5 0.00000 0.00000 0.00000 0.00000 1.53210 0.12319 0.12319
6 0.00000 0.19228 0.00000 0.22438 0.12319 0.79770 -0.01450
7 0.00000 0.19228 0.00000 0.22438 0.12319 -0.01450 0.79770
求O-H鍵的這種鍵級也就是∑[i=1,2,3,4,5]T(i,7)=0+0.19228+0+0.22438+0.12319=0.53985，和NBO輸出的一致。
另外值得一提的是，所謂NPA電荷（無數人謬稱為NBO電荷），就是以NAO為基獲得的Mulliken電荷。A原子布居數=∑[i∈A]∑[j]P(i,j)S(i,j)=∑[i∈A]∑[j]P(i,j)δ(i,j)=∑[i∈A]P(i,i)，所以氧的布居數∑[i=1,2,3,4,5]DMNAO(i,i)=2.0000+1.7325+2.0000+1.1399+1.5321=8.4045，凈電荷為8-8.4045=-0.4045，和NBO程序輸出的一致。
4. 第三種鍵級(MO bond order)
這個鍵級的名字十分令人疑惑，而且手冊上也沒說明。實際上它的定義極其簡單，A、B間的鍵級為∑[i∈A]∑[j∈B]P(i,j)。NBO輸出為
Atom 1 2 3
---- ------ ------ ------
1. O 0.0000 0.5675 -0.8328
2. H 0.5675 0.0000 -0.0623
3. H -0.8328 -0.0623 0.0000
求O-H2鍵級就是∑[i=1,2,3,4,5]DMNAO(i,7)=0.4113-0.7001-0.5439=-0.8327。結果十分荒誕，明明是典型sigma鍵，鍵級卻為負。而另一個等價的O-H3鍵鍵級卻是正值0.4113+0.7001-0.5439=0.5675！可見，這種鍵級從實際應用上看毫無意義。從理論上同樣沒有根據，從形式上看和Mulliken鍵級相仿佛，卻少乘了重疊矩陣元，另外NAO又是正交基，不存在交叉項，也不可能求Mulliken鍵級。而與Wiberg鍵級相比，又少求了平方。我認為NBO程序中應該去掉這種毫無意義的鍵級。
5. 用PNAO計算Mayer鍵級的嘗試
原理上，以非正交基PNAO做Mayer鍵級分析也是可以的，可以代替正交基NAO下的Wiberg鍵級分析，另外，做Mulliken分析通過交叉項還可以研究成鍵反鍵和各分子軌道貢獻，但是NBO程序里沒有提供。筆者嘗試利用PNAO做Mayer鍵級分析，這首先需要獲得PNAO密度矩陣，NBO3.1里不支持輸出它的DMPNAO關鍵詞，故通過.47文件在GENNBO5.0下輸出，所得矩陣如下：
NAO 1 2 3 4 5 6 7
---------- ------- ------- ------- ------- ------- ------- -------
1. O 1( 1s ) 2.0000 0.0000 0.0000 0.0000 0.0000 0.1870 0.1870
2. O 1( 2s ) 0.0000 1.9767 0.0000 0.0000 0.0759 1.0236 1.0236
3. O 1( 2px) 0.0000 0.0000 2.0000 0.0000 0.0000 0.0000 0.0000
4. O 1( 2py) 0.0000 0.0000 0.0000 1.5602 0.0000 0.9391 -0.9391
5. O 1( 2pz) 0.0000 0.0759 0.0000 0.0000 1.7530 -0.7415 -0.7415
6. H 2( 1s ) 0.1870 1.0236 0.0000 0.9391 -0.7415 1.4612 0.3307
7. H 3( 1s ) 0.1870 1.0236 0.0000 -0.9391 -0.7415 0.3307 1.4612
然而，不幸的是，所得Mayer鍵級矩陣很糟
1 2 3
1 0.00000000 7.23697119 7.23697119
2 7.23697119 0.00000000 1.06885215
3 7.23697119 1.06885215 0.00000000
鍵級大得出奇。NBO5.0輸出的PNAO密度矩陣估計錯誤的，理論上Tr(PS)應當等于分子總電子數，然而用PNAO的密度矩陣和重疊矩陣竟得到16.22606912。雖然測試這個方法由于NBO程序的bug失敗了，但原理上是合理的。