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

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

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

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

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

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

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

圖1 表1中實(shí)例的分子結(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). 所以僅憑能量相等或能量差小于某個(gè)閾值判斷分子結(jié)構(gòu)的獨(dú)特性會得到錯(cuò)誤結(jié)論.

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

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

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

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

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

2.1 結(jié)構(gòu)去重結(jié)果比較 C2Si2Xq和C2Al4H4的初始結(jié)構(gòu)分別為24 892,4 372個(gè)， 用STRUauto算法找出異構(gòu)體的數(shù)目分別為3 442,552個(gè). C2Si2Xq體系內(nèi)的結(jié)構(gòu)去重結(jié)果列于表2. 由表2可見, 由于STRUauto算法考慮了結(jié)構(gòu)總能量、 原子與原子間的鍵長以及2個(gè)原子類型等因素, 因此去除了初始結(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ù)進(jìn)行統(tǒng)計(jì)分析. C2Si2Ge-體系的初始結(jié)構(gòu)有396個(gè), 經(jīng)程序去重, STRUauto找到same文件夾(所有存在相同結(jié)構(gòu)的文件集合)中的結(jié)構(gòu)有316個(gè), 最終確定80個(gè)異構(gòu)體, 去重率為79.8%. 經(jīng)確認(rèn), 找出的80個(gè)結(jié)構(gòu)均為獨(dú)特結(jié)構(gòu), 即發(fā)現(xiàn)異構(gòu)體的準(zhǔn)確率為100%. 表2中所有分子體系找出的結(jié)構(gòu)均為獨(dú)特結(jié)構(gòu), 發(fā)現(xiàn)異構(gòu)體的準(zhǔn)確率均為100%.

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

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