高斯計算中分子結(jié)構(gòu)獨特性的判斷

王海燕, 崔文超, 于 爽, 丁益宏

(1. 吉林師范大學(xué) 計算機學(xué)院, 吉林 四平 136000； 2. 吉林大學(xué) 理論化學(xué)研究所, 長春 130021)

化學(xué)反應(yīng)過程中涉及大量的中間體和過渡態(tài), 在確定這些中間體和過渡態(tài)過程中需構(gòu)建數(shù)量更多的中間體和過渡態(tài)作為化學(xué)結(jié)構(gòu)計算程序(如高斯軟件包)的輸入結(jié)構(gòu). 通過高斯計算獲得有效結(jié)構(gòu), 進一步確定低能結(jié)構(gòu)或有效的化學(xué)反應(yīng)通道. 由于存在重復(fù)結(jié)構(gòu), 因此需要將其去除, 留下獨特的結(jié)構(gòu). 分子結(jié)構(gòu)產(chǎn)生有多種途徑, 如通過假設(shè)或借助計算機算法產(chǎn)生初始結(jié)構(gòu), 并判斷產(chǎn)生的結(jié)構(gòu)是否合理, 以及是否存在結(jié)構(gòu)的獨特性等. 在判斷結(jié)構(gòu)獨特性的過程中, 可通過以下過程實現(xiàn)： 將初始結(jié)構(gòu)輸入高斯計算程序中進行深層次計算、 優(yōu)化, 最終得到一個終止結(jié)構(gòu). 但相同的終止結(jié)構(gòu)會存在于大量的高斯計算輸出文件中. 在去重過程中, 簡單結(jié)構(gòu)可通過肉眼分辨, 但若分子中原子數(shù)較多則會導(dǎo)致分子的空間構(gòu)型復(fù)雜, 去重操作困難.

計算化學(xué)的主要目標是利用有效的數(shù)學(xué)近似以及電腦程序計算分子的性質(zhì)(如總能量、 偶極矩、 四極矩、 振動頻率、 反應(yīng)活性等), 并解釋一些具體的化學(xué)問題[1-2]. 在計算化學(xué)應(yīng)用領(lǐng)域, 常將已有的電腦程序和方法針對特定的化學(xué)問題進行研究, 如在構(gòu)建分子結(jié)構(gòu)過程中, 可利用計算機輔助方式; 在結(jié)構(gòu)搜索中, 需判斷結(jié)構(gòu)的獨特性, 如C2Si2Xq(五原子高周期14族平面四配位化合物)[3-4]、 C2Al4H4[5]、 鋁氫團簇和氮雜環(huán)卡賓等[6], 存在大量結(jié)構(gòu)相同的分子, 針對計算機程序產(chǎn)生的大量結(jié)構(gòu)文件, 若不去重, 則會導(dǎo)致大量冗余工作. 而僅憑能量相等判斷分子結(jié)構(gòu)的獨特性易出現(xiàn)失誤. 針對上述問題, 本文將計算機輔助方式應(yīng)用于結(jié)構(gòu)搜索子方向, 提出一種自動化結(jié)構(gòu)獨特性判斷STRUauto算法, 該算法借助能量和鍵長的迭代比較, 依次完成分子結(jié)構(gòu)的批量去重操作, 可避免常規(guī)分子結(jié)構(gòu)獨特性判斷的誤差, 準確得到分子的異構(gòu)體. 為驗證STRUauto算法的優(yōu)勢, 先對C2Si2Xq和C2Al4H4兩類典型的分子結(jié)構(gòu)進行實驗測試, 再對C2Si2Xq中的六類典型體系進一步實驗. 分別考慮常規(guī)狀態(tài)下和應(yīng)用STRUauto狀態(tài)下的多項典型技術(shù)指標比對, 實驗結(jié)果表明, STRUauto算法能有效提升理論化學(xué)中結(jié)構(gòu)搜索問題的效率, 其準確率和可信度良好.

1 自動化結(jié)構(gòu)獨特性判斷算法

1.1 問題描述 兩個具有相同原子數(shù)目的分子其結(jié)構(gòu)可能存在差異, 如C2Al4H4可排列成多種分子結(jié)構(gòu). 研究表明, 當(dāng)原子數(shù)目增加時, 分子結(jié)構(gòu)數(shù)目會成倍增加. 除通過肉眼和經(jīng)驗值判斷C2Al4H4分子可能存在獨特的結(jié)構(gòu), 目前還有許多輔助性的識別方法, 如借助總能量值判斷. 計算出分子結(jié)構(gòu)后, 可計算出一個能量. 總能量值相近的分子, 易誤認為結(jié)構(gòu)相同. 一些能量相近但結(jié)構(gòu)不同的實例列于表1, 表1中實例的分子結(jié)構(gòu)構(gòu)型如圖1所示.

表1 能量相近但結(jié)構(gòu)不同實例Table 1 Examples of similar energy but different structures

圖1 表1中實例的分子結(jié)構(gòu)構(gòu)型Fig.1 Molecular structural configuration of examples in tabel 1

由表1可見: 在C2Si2Ge-體系中, 第一組中C2Si2Ge-_046和C2Si2Ge-_010的能量接近, 其能量差僅為0.002 646 39 a.u., 但其結(jié)構(gòu)構(gòu)型完全不同； C2Al4_0025+H_mt_2+H_md_2+H_ms_3(簡寫)和C2Al4_0108+H_mt_3+H_ms_2+H_md_2(簡寫)的能量更接近, 能量差僅為0.000 068 a.u., 但其結(jié)構(gòu)構(gòu)型完全不同； C2Al4_0079+H_mt_4+H_ms_3+H_mt_1(簡寫)和C2Al4_0079+H_mt_4+H_mt_3+H_mt_3(簡寫)的能量差僅為0.000 01 a.u., 結(jié)構(gòu)構(gòu)型相似, 但甲基的轉(zhuǎn)動位置不同, 因此二者屬于兩種不同結(jié)構(gòu). 所以僅憑能量相等或能量差小于某個閾值判斷分子結(jié)構(gòu)的獨特性會得到錯誤結(jié)論.

1.2 算 法 通?？捎煞肿又性觽€數(shù)及各原子的空間位置判斷一種分子的結(jié)構(gòu)獨特性. 但對于原子個數(shù)相同且空間位置相似的分子, 即使加入經(jīng)驗值, 判斷也很困難. 因此， 需要一種自動化的結(jié)構(gòu)獨特性判斷算法， 以提高判別的準確度和效率. 研究表明, 分子總能量是判別結(jié)構(gòu)獨特性的一個關(guān)鍵因素. 若最初按分子結(jié)構(gòu)的總能量進行判斷, 則可先區(qū)分出一些異構(gòu)體(中間體), 即能量不同, 其結(jié)構(gòu)不同, 進而確定部分分子的獨特性. 其次, 分子中各原子間的鍵長是另一個關(guān)鍵因素. 總能量相同的分子, 結(jié)構(gòu)內(nèi)部原子間鍵長不同也會導(dǎo)致兩種分子的異構(gòu). 基于此, 本文提出一種自動化判斷兩種分子結(jié)構(gòu)獨特性的STRUauto算法, 算法描述為： 1) 預(yù)處理過程; 2) 判斷總能量的一致性； 3) 判斷鍵長的一致性; 4) 去重或確定分子結(jié)構(gòu)獨特性. 其工作流程如圖2所示.

在判斷分子總能量一致性的過程中, 先通過化學(xué)結(jié)構(gòu)計算程序(如高斯軟件包等)確定所有中間體總能量, 再計算中間體之間的能量差, 這是決定異構(gòu)體獨特性的首要因素. 考察能量差需預(yù)先設(shè)定一個閾值, 實驗結(jié)果表明, 閾值應(yīng)設(shè)為9×10-5a.u., 若2種分子結(jié)構(gòu)的能量差大于閾值, 則2個結(jié)構(gòu)必不同, 進而確定為獨特結(jié)構(gòu)； 若能量差小于9×10-5a.u., 則需進一步判斷. 在進一步判斷過程中, 主要決定因素是原子與原子間的鍵長以及2個原子的類型. 鍵長決定2個結(jié)構(gòu)的內(nèi)部距離, 內(nèi)部距離構(gòu)成分子的空間結(jié)構(gòu). 由于相同的鍵長值可能對應(yīng)不同的原子組合, 但不需要判斷鍵長的一致性, 因此, 需要判斷是否對鍵長升序排序. 對于相同的原子組合, 才將這些距離按升序規(guī)則排序.

2 常規(guī)方法與自動化方法比較

為證明STRUauto自動化結(jié)構(gòu)獨特性判斷算法的優(yōu)越性, 將STRUauto算法與常規(guī)異構(gòu)體判斷方法進行對比. 實驗環(huán)境為： Intel(R) Xeon(R), X5570, 2.93 GHz處理器, 32 GB內(nèi)存[7].

為驗證STRUauto算法的效率[8], 選取典型的分子組合C2Si2Xq和C2Al4H4為研究對象, 分別研究具有相同原子組合的分子結(jié)構(gòu)獨特性, 其技術(shù)指標包括初始結(jié)構(gòu)數(shù)目、 異構(gòu)體數(shù)目、 異構(gòu)體時間、 去重比例和準確率. 通過上述指標比較常規(guī)方法與自動化判斷STRUauto算法的差異, 并對上述指標進行統(tǒng)計分析.

2.1 結(jié)構(gòu)去重結(jié)果比較 C2Si2Xq和C2Al4H4的初始結(jié)構(gòu)分別為24 892,4 372個， 用STRUauto算法找出異構(gòu)體的數(shù)目分別為3 442,552個. C2Si2Xq體系內(nèi)的結(jié)構(gòu)去重結(jié)果列于表2. 由表2可見, 由于STRUauto算法考慮了結(jié)構(gòu)總能量、 原子與原子間的鍵長以及2個原子類型等因素, 因此去除了初始結(jié)構(gòu)中的重復(fù)結(jié)構(gòu), 有效得到了相應(yīng)數(shù)目的異構(gòu)體.

表2 C2Si2Xq體系內(nèi)的結(jié)構(gòu)去重結(jié)果Table 2 Deduplication results of structure in C2Si2Xq system

以C2Si2Ge-為例對數(shù)據(jù)進行統(tǒng)計分析. C2Si2Ge-體系的初始結(jié)構(gòu)有396個, 經(jīng)程序去重, STRUauto找到same文件夾(所有存在相同結(jié)構(gòu)的文件集合)中的結(jié)構(gòu)有316個, 最終確定80個異構(gòu)體, 去重率為79.8%. 經(jīng)確認, 找出的80個結(jié)構(gòu)均為獨特結(jié)構(gòu), 即發(fā)現(xiàn)異構(gòu)體的準確率為100%. 表2中所有分子體系找出的結(jié)構(gòu)均為獨特結(jié)構(gòu), 發(fā)現(xiàn)異構(gòu)體的準確率均為100%.

為驗證STRUauto算法在異構(gòu)體去重方面的優(yōu)勢, 將自動化結(jié)構(gòu)獨特性判斷STRUauto算法與常規(guī)人工分子結(jié)構(gòu)獨特性判斷方法在去重率方面進行比較. 實驗結(jié)果表明, 常規(guī)方法的去重率均小于50%, STRUauto算法的去重率均大于70%.

綜上所述, 本文利用計算機輔助分子異構(gòu)體構(gòu)建了STRUauto算法， 結(jié)果表明： 在結(jié)構(gòu)搜索過程中, 即使2個結(jié)構(gòu)的能量相同, 其分子結(jié)構(gòu)也可能不同, 因此需要進一步借助其他輔助因素判斷, 最合適的輔助因素是原子間的鍵長以及2個原子的類型； 當(dāng)結(jié)構(gòu)能量差值較小時, 通過2組n(n-1)/2個鍵長的比較可確定2個結(jié)構(gòu)是否相同. STRUauto算法針對結(jié)構(gòu)搜索體系中原子個數(shù)增多與篩選準確率之間的矛盾, 較精準地完成了化學(xué)結(jié)構(gòu)搜索的去重問題. 當(dāng)分子總能量、 原子數(shù)和鍵長均相同時, 可排除存在異構(gòu)體的可能. 通過STRUauto算法與常規(guī)方法的數(shù)據(jù)比較, 證明其準確率和效率均有較大提高.