向 垚 鄭斯瑞 賴向京

（南京郵電大學(xué)先進技術(shù)研究院 南京 210023）

1 引言

全局優(yōu)化問題廣泛存在于各個領(lǐng)域，其中通過搜索原子團簇勢能函數(shù)的全局最小值，是計算化學(xué)領(lǐng)域中一個著名的全局優(yōu)化問題。在計算化學(xué)領(lǐng)域，原子團簇的許多物理和化學(xué)性質(zhì)是由其幾何結(jié)構(gòu)決定的，因此確定它們的基態(tài)結(jié)構(gòu)對于理解它們的各種性質(zhì)具有重要意義［1］。由于團簇的局部最優(yōu)結(jié)構(gòu)的數(shù)量隨團簇規(guī)模呈指數(shù)型增長，確定原子團簇基態(tài)結(jié)構(gòu)（即團簇的幾何優(yōu)化）是一項非常具有挑戰(zhàn)性的任務(wù)。事實上，原子團簇優(yōu)化問題也已經(jīng)被證明為NP-hard問題［2］。

原子團簇的結(jié)構(gòu)優(yōu)化問題由于其NP-hard 特性和應(yīng)用價值，受到了許多研究者的廣泛研究。在文獻中已有大量相關(guān)算法被提出，這些方法可以分為兩種，包括基于單解的算法和基于種群的算法。對于基于種群算法，具有代表性的算法包括隨機隧道優(yōu)化算法（Random Tunneling Algorithm，RTA）［3］，自適應(yīng)免疫優(yōu)化算法（Adaptive Immune Optimization Algorithm，AIOA）［4］，構(gòu)型空間退火算法（Conformational Space Annealing，CSA）［5］，基于種群的盆地跳躍算法（Population Basin Hopping，PBH）［6］等。對于基于單解的算法，具有代表性的算法包括basin-hopping 算法（BH）［2］，帶有表面算子和內(nèi)部算子的啟發(fā)式算法（Heuristic Algorithm with the Surface and Interior Operators，HA-SIO）［7～8］，動態(tài)格子搜索算法（Dynamic Lattice Searching，DLS）及其變種［9～11］，和模糊全局優(yōu)化算法（Fuzzy Global Optimization，F(xiàn)GO）［12］等。

PBH算法是一種混合進化算法，在許多全局優(yōu)化問題上表現(xiàn)出了很高的性能，例如圓形packing問題［13］。本文通過整合PBH 算法的框架和基于連通表（Connectivity Table，CT）［4］的距離函數(shù)，為Gupta 勢能函數(shù)建模的原子團簇的結(jié)構(gòu)優(yōu)化提出了PBH算法的一個變種，并將該算法成功應(yīng)用于銀團簇的結(jié)構(gòu)優(yōu)化。計算結(jié)果很好地表明了該算法的有效性。

2 Gupta勢能函數(shù)模型

Gupta 勢能［14］常用于給金屬團簇建模，其勢能函數(shù)具有如下形式：

其中N為原子數(shù)目，rij為原子i與原子j之間的歐式距離，r0為團簇兩原子之間的平衡距離，r0=1.0，A，p，q，UN為勢能函數(shù)的參數(shù)，這些參數(shù)的不同組合能表示不同類型的原子團簇。對于銀團簇，其相應(yīng)的參數(shù)取值如下：

3 基于連通表的種群盆地跳躍算法（PBH-CT）

3.1 PBH-CT算法的主框架

本文提出的PBH-CT 算法是一種混合進化算法，它由一個種群初始化方法、一個局部優(yōu)化方法，一種擾動算子和一個種群更新策略組成，其偽代碼如算法1所述。

從一個隨機生成的m個個體構(gòu)成的初始種群開始，算法執(zhí)行若干次迭代，直到連續(xù)MaxNoImp次迭代無法改進當前最好解為止，其中MaxNoImp 是一個參數(shù)。在每個迭代，首先找到種群X中與當前解Yk距離最近的構(gòu)型Xr，然后通過將兩個構(gòu)型之間的距離與接收罰值dcut進行比較，來判斷當前構(gòu)型Yk是否與種群X的解Xr相似。

若d(Yk，Xr)＞dcut，則表示兩個構(gòu)型Yk與Xr是不相似的。此時，我們比較Yk與種群X中的最差構(gòu)型Xw的目標函數(shù)值；若E(Yk)＜E(Xw)，則用Yk替換種群X中的Xw；若d(Yk，Xr)≤dcut，則認為兩個構(gòu)型是相似的，此時，我們比較這兩個構(gòu)型的目標函數(shù)值；若E(Yk)＜E(Xr)，則用Yk替換種群X中的Xr。

3.2 初始種群

在基于種群的算法中，種群的多樣性對算法性能起著關(guān)鍵作用。因此，為了保證初始種群X的多樣性，本文采用了如下的初始種群策略：隨機生成m個構(gòu)型來構(gòu)成初始種群，每個個體都是在一個半徑為R=（3N/（4Π））1/3·r0的球體中隨機撒入N個原子中心點，從而得到一個隨機構(gòu)型，并通過局部優(yōu)化方法對該構(gòu)型進行局部優(yōu)化，從而得到初始種群X。在本文算法中，局部優(yōu)化方法采用的是L-BFGS［16］算法，其中L-BFGS 算法是一種高效的局部優(yōu)化算法，已被廣泛應(yīng)用于可導(dǎo)函數(shù)的局部優(yōu)化。

3.3 擾動

對于團簇優(yōu)化問題，其勢能函數(shù)具有大量的局部最優(yōu)解，且局部最優(yōu)值的個數(shù)隨著原子個數(shù)N呈指數(shù)關(guān)系增長。因此，對于規(guī)模較大的算例，如何跳出這些局部陷阱變成了解決團簇優(yōu)化問題的主要難點。為了避免陷入局部最優(yōu)陷阱，本文采用經(jīng)典的隨機擾動策略。具體而言，給定一個候選解x=(x1，x2，…，x3N)，通過隨機小擾動來生成后代解，即xi←xi+rand(-Δ，Δ) (1 ≤i≤3N)，其中rand(-Δ，Δ)表示[-Δ，Δ]中的一個隨機數(shù)，Δ 是一個參數(shù)，在本文中設(shè)置成0.7。

3.4 構(gòu)型間的距離函數(shù)

在該算法中，采用團簇的最近鄰接數(shù)（記為nnn）和連通表（CT）［4］來定義兩個團簇構(gòu)型的距離函數(shù)d（A，B），具體而言，對于一個團簇來說，nnn的定義如下：

其中rij表示原子i和原子j之間的歐式距離，r0為原子間的平衡距離。

另外，CT 是一個基于團簇構(gòu)型的拓撲信息的向量，由團簇中原子的最近鄰居個數(shù)生成。CT 為14 維向量，其中CT［i］（1 ≤i≤14）表示在團簇中有i個鄰居原子的原子個數(shù)。值得注意的是，在所研究的原子團簇中，原子的最近鄰居數(shù)均不超過14。

借助連通表（CT）和最近鄰接數(shù)（nnn），可以如下定義團簇構(gòu)型A和團簇構(gòu)型B之間的距離函數(shù)d（A，B）：

其中，（或）是構(gòu)型A（或B）的最近鄰接數(shù)。CT（Ai）（或CT（Bi））是構(gòu)型A（或B）中有i個鄰居原子的原子個數(shù)。

4 實驗結(jié)果與比較

4.1 實驗條件

本文算法通過C 語言實現(xiàn)，并在Windows 10操作系統(tǒng)、主頻為2.5GHz、內(nèi)存為8GB 的PC 機上進行實驗。算法參數(shù)設(shè)置如下：m=20，MaxNoImp=400。另外，由于算法的隨機性，對于n∈［10，61］范圍內(nèi)的每個銀團簇，本文算法獨立地運行了10次。

4.2 實驗結(jié)果

實驗結(jié)果總結(jié)在表1 中，其中第1 列給出了團簇中的原子個數(shù)（N），第2 列給出了文獻中已知的最優(yōu)解，本文所提出的算法的計算結(jié)果在隨后3 列中給出，包括算法的最優(yōu)解，達到最優(yōu)解所需要的平均計算時間和達到最好解的成功率。

表1 計算結(jié)果與當前文獻中最好的結(jié)果比較，其中最好的結(jié)果來自參考文獻［3］和［9］

從表1 中可以看出，所提出的算法在大多數(shù)測試的算例中都可以找到當前文獻中最好的結(jié)果，這表明所提出的算法在10 ≤N≤61 的范圍內(nèi)具有強的搜索能力。然而，該算法在4 個難例（即Ag53、Ag54、Ag58、Ag59）上未找到當前已知的最優(yōu)結(jié)果，其計算結(jié)果在表1中用斜體表示。

對于10 ≤N≤29 的小團簇來說，所提出的算法在平均計算時間和算法成功率上都具有較高的性能。但是隨著原子個數(shù)N的增加，算法的計算時間會顯著增加。

另外，為了直觀地表示計算結(jié)果，我們在圖1中分別給出4 個代表性團簇構(gòu)型的最優(yōu)構(gòu)型的圖形。如圖1所示，Ag33的構(gòu)型是基于十面體的結(jié)構(gòu)，Ag39的構(gòu)型是基于面心立方體（FCC）的結(jié)構(gòu)，Ag46的構(gòu)型是無序形態(tài)的，而Ag49的構(gòu)型是基于二十面體的結(jié)構(gòu)。

圖1 4個代表性算例的最優(yōu)構(gòu)型

5 結(jié)語

本文針對Gupta 勢能描述的銀團簇的結(jié)構(gòu)優(yōu)化問題，提出了一個PBH 算法的新變種，并對原子數(shù)在10 ≤N≤61范圍內(nèi)的銀團簇進行了優(yōu)化。與文獻中已有的PBH 算法不同的是，本文提出的PBH 算法的變種利用了連通表來定義團簇構(gòu)型間的距離函數(shù)。計算結(jié)果表明，該算法對于Gupta 勢能描述的銀團簇的結(jié)構(gòu)優(yōu)化是有效的。具體而言，在10 ≤N≤61范圍內(nèi)，提出的算法能夠獲得除4個難例以外的已知全局最優(yōu)構(gòu)型。實驗結(jié)果再次表明，連通表是定義團簇構(gòu)型間距離函數(shù)的一種很好的工具。在將來，我們打算在其他原子團簇上驗證它的有效性。