雙重變異遺傳算法及其在臨界滑動(dòng)面搜索中的應(yīng)用

覃 偉

（重慶工程職業(yè)技術(shù)學(xué)院，重慶 402260）

在公路工程、建筑工程、礦山工程、水利工程等領(lǐng)域中，常會(huì)遇到大量的邊坡穩(wěn)定性問題。確定邊坡臨界滑動(dòng)面對(duì)邊坡穩(wěn)定性分析及坍岸寬度預(yù)測(cè)等研究具有重要意義[1]。

邊坡臨界滑動(dòng)面的搜索，本質(zhì)上是一個(gè)優(yōu)化問題，其具有變量多、約束條件復(fù)雜、局部極值點(diǎn)多、非凸性等特點(diǎn)，使得傳統(tǒng)的優(yōu)化算法（如模式搜索法）在搜索臨界滑動(dòng)面時(shí)，常常陷入局部極值點(diǎn)，給邊坡的穩(wěn)定性分析帶來許多困難。

遺傳算法（genetic algorithm，GA）由 Holland 教授于 20世紀(jì) 70 年代提出，是一種模擬達(dá)爾文遺傳選擇和自然淘汰的生物進(jìn)化過程的搜索尋優(yōu)算法，具有思想簡(jiǎn)單、易于編程實(shí)現(xiàn)、算法健壯等優(yōu)點(diǎn)，同時(shí)還具有隱含并行性和全局搜索等顯著特性[2]。因此，一些學(xué)者將遺傳算法運(yùn)用到邊坡臨界滑動(dòng)面搜索研究中，取得了一定成效[3?8]。如，賀子光等[3]將單純形法和回溯機(jī)制引入GEP 中，提出了混合GEP 方法搜索邊坡的臨界滑動(dòng)面；梁冠亭等[4]采用自適應(yīng)遺傳優(yōu)化算法搜索采用抗滑樁支護(hù)邊坡的臨界滑動(dòng)面；朱劍鋒等[5]提出自適應(yīng)禁忌變異遺傳算法搜索土釘墻臨界滑動(dòng)面。但是目前的標(biāo)準(zhǔn)遺傳算法及其改進(jìn)算法仍在不同程度上存在收斂速度較慢、容易進(jìn)入早熟等缺陷，有待進(jìn)一步改進(jìn)和完善。

本文在傳統(tǒng)遺傳算法的基礎(chǔ)上，提出雙重變異策略對(duì)遺傳算法進(jìn)行改進(jìn)，形成雙重變異遺傳算法（DMGA），并結(jié)合簡(jiǎn)化Bishop 法，編寫Python 程序，通過對(duì)土坡算例臨界滑動(dòng)面的計(jì)算，進(jìn)一步為邊坡臨界滑動(dòng)面搜索提供有效方法。

1 臨界滑動(dòng)面搜索問題的目標(biāo)函數(shù)

求解邊坡的臨界滑動(dòng)面，即是搜索使邊坡安全系數(shù)最小的滑動(dòng)面。邊坡穩(wěn)定性計(jì)算方法有很多，主要包括剛體極限平衡方法和數(shù)值模擬方法[9?10]。對(duì)于滑動(dòng)面呈圓弧形的土質(zhì)邊坡，由于簡(jiǎn)化Bishop 法計(jì)算簡(jiǎn)單、物理意義明確，在工程實(shí)踐中應(yīng)用廣泛。因此本文采用簡(jiǎn)化Bishop 法計(jì)算土質(zhì)邊坡的穩(wěn)定性。

1.1 簡(jiǎn)化Bishop 法

簡(jiǎn)化Bishop 法是一種適用于滑動(dòng)面呈圓弧形的邊坡的穩(wěn)定性分析方法。其中邊坡安全系數(shù)F的計(jì)算公式[11]為：

式中：ci—第i條塊滑動(dòng)面上巖土體的黏聚力/kPa；

φi—第i條塊滑動(dòng)面上巖土體的內(nèi)摩擦角/（°）；

li—第i條塊滑動(dòng)面的弧長(zhǎng)/m；

Wi—第i條塊的單寬重量/（kN· m?1）；

αi—第i條塊滑動(dòng)面傾角/（°），當(dāng)滑動(dòng)面傾向與滑動(dòng)方向一致時(shí)，αi取正；當(dāng)滑動(dòng)面傾向與滑動(dòng)方向相反時(shí)，αi取負(fù)。

當(dāng)變量mi很小時(shí)，計(jì)算結(jié)果與實(shí)際情況不符。根據(jù)一些學(xué)者的意見，當(dāng)mi?0.2時(shí)，不能采用簡(jiǎn)化B ishop 法計(jì)算邊坡安全系數(shù)[11]。

1.2 目標(biāo)函數(shù)的確定

邊坡臨界滑動(dòng)面的位置與圓弧形滑動(dòng)面的圓心半徑R、圓心橫坐標(biāo)xo、圓心縱坐標(biāo)yo有關(guān)，因此將這3 個(gè)量作為目標(biāo)函數(shù)的求解變量。

并不是任意半徑及圓心位置的圓弧都能計(jì)算出邊坡的安全系數(shù)，如發(fā)生圓弧與邊坡不存在交點(diǎn)，或所采用的邊坡安全系數(shù)計(jì)算方法不適用等情況。因此，臨界滑動(dòng)面搜索問題可表示為在一定約束條件下求解邊坡安全系數(shù)最小值，即求解下面約束問題：

運(yùn)用外點(diǎn)罰函數(shù)法[12]將該約束問題轉(zhuǎn)化為無約束問題：

其中：

其中，當(dāng)圓弧滑動(dòng)面與坡體無交點(diǎn)時(shí)，F(xiàn)（x）不存在，則令其為一較大值100；σ為很大的正數(shù)，取值1 000。

2 雙重變異遺傳算法的設(shè)計(jì)

2.1 個(gè)體編碼與初始種群產(chǎn)生

本文對(duì)變量采用一種特殊的實(shí)數(shù)編碼，在編碼中設(shè)置符號(hào)基因位與數(shù)字基因位。符號(hào)基因位存放實(shí)數(shù)的符號(hào)信息，當(dāng)實(shí)數(shù)為正或0 時(shí)，該基因位的值取1；實(shí)數(shù)為負(fù)時(shí)，該基因位的值取?1。數(shù)字基因位存放實(shí)數(shù)的數(shù)字部分，實(shí)數(shù)各位上的數(shù)字占一個(gè)基因位，且該位上的數(shù)字即為對(duì)應(yīng)基因位的值。對(duì)每個(gè)變量xi均采用上述實(shí)數(shù)編碼，并將它們的編碼依次連接構(gòu)成一個(gè)染色體。例如，若目標(biāo)函數(shù)的解由3 個(gè)變量構(gòu)成，其中變量x1為36.12，變量x2為?5.099，變量x3為72.562，則對(duì)其進(jìn)行實(shí)數(shù)編碼如圖1所示。

圖1 實(shí)數(shù)編碼Fig.1 Real number coding

初始種群內(nèi)的個(gè)體數(shù)為N，每個(gè)個(gè)體通過隨機(jī)數(shù)產(chǎn)生，個(gè)體染色體中的每個(gè)數(shù)字基因位隨機(jī)生成[0,9]區(qū)間的整數(shù)，每個(gè)符號(hào)基因位隨機(jī)生成1 或?1，若該變量不為負(fù)，則對(duì)應(yīng)符號(hào)基因位僅產(chǎn)生1。

2.2 適應(yīng)度值函數(shù)

遺傳算法在進(jìn)化時(shí)，以個(gè)體的適應(yīng)度值作為搜索依據(jù)。求解優(yōu)化問題時(shí)，需要將目標(biāo)函數(shù)轉(zhuǎn)化為適應(yīng)度函數(shù)。由于目標(biāo)函數(shù)G（x）為最小值問題，構(gòu)造適應(yīng)度值函數(shù)如下：

適應(yīng)度函數(shù)值f（x）越大，邊坡的安全系數(shù)就越小。

2.3 雙重變異遺傳算法的進(jìn)化過程

2.3.1 選擇操作

父代個(gè)體的適應(yīng)度值確定后，按照適應(yīng)度值進(jìn)行個(gè)體的選擇。本文種群的選擇采用精英保留策略[13?14]，每一代種群進(jìn)化完成后，首先從中找出適應(yīng)度值最大的一個(gè)個(gè)體，直接讓其進(jìn)入下一代，不再對(duì)其進(jìn)行交叉、變異操作。剩余的個(gè)體，則采用轉(zhuǎn)盤式選擇方式進(jìn)行選擇操作，轉(zhuǎn)盤式選擇的基本原理是根據(jù)每個(gè)個(gè)體的適應(yīng)度值的比例來確定該個(gè)體的選擇概率[15]。

轉(zhuǎn)盤式選擇的計(jì)算方法：先根據(jù)個(gè)體的適應(yīng)度值fi按式（8）計(jì)算個(gè)體的選擇概率pi，并令p0=0，然后生成[0,1]區(qū)間的隨機(jī)數(shù)r，如果r滿足式（7），則選擇個(gè)體j。

其中：

2.3.2 自適應(yīng)交叉操作

在遺傳算法中，通過交叉操作既可以使種群多樣化，又可以引導(dǎo)搜索方向，促進(jìn)優(yōu)秀個(gè)體的產(chǎn)生。本文采取單變量交叉方式，即從種群中選出任意2 個(gè)個(gè)體（個(gè)體1、個(gè)體2）進(jìn)行交叉，當(dāng)產(chǎn)生的[0,1）區(qū)間的隨機(jī)數(shù)小于交叉概率Pc時(shí)，隨機(jī)選擇一個(gè)變量xi，使個(gè)體1 與個(gè)體2 中變量xi所對(duì)應(yīng)的基因位的值進(jìn)行交換。

交叉概率是影響遺傳算法性能的關(guān)鍵因素之一。盧明奇等[16]為了解決傳統(tǒng)的自適應(yīng)遺傳算法在進(jìn)化初期常存在停滯現(xiàn)象、進(jìn)化后期易陷入局部極值點(diǎn)的問題，提出了一種改進(jìn)的自適應(yīng)交叉概率。但該方法參數(shù)較多，增加了參數(shù)設(shè)置的難度。本文在考慮適應(yīng)度值及進(jìn)化代數(shù)的情況下，提出一種復(fù)合自適應(yīng)交叉概率，使交叉概率Pc隨個(gè)體適應(yīng)度值及進(jìn)化代數(shù)進(jìn)行調(diào)整，Pc的計(jì)算公式如下：

式中：fmax—當(dāng)前種群中的最大適應(yīng)度值；

f′—待交叉?zhèn)€體中較大的適應(yīng)度值；

—當(dāng)前種群的平均適應(yīng)度值；

T—總的進(jìn)化代數(shù)；

t—當(dāng)前進(jìn)化代數(shù)。

2.3.3 雙重變異操作

在傳統(tǒng)的遺傳算法中，變異操作本質(zhì)上是隨機(jī)搜索，并不能保證產(chǎn)生改進(jìn)的后代[17]，進(jìn)而常導(dǎo)致遺傳算法的收斂速度較慢。為此，本文借鑒模式搜索法中探測(cè)移動(dòng)的思想，提出探測(cè)變異操作，并與傳統(tǒng)的直接變異操作相結(jié)合，構(gòu)成雙重變異策略對(duì)遺傳算法進(jìn)行改進(jìn)，使遺傳算法的收斂速度及全局搜索能力得到提升。雙重變異操作包括探測(cè)變異操作與直接變異操作。其中，探測(cè)變異操作將提升算法的局部尋優(yōu)能力，直接變異操作將使算法具備全局尋優(yōu)能力。

為了敘述方便，將采用雙重變異策略的遺傳算法稱為雙重變異遺傳算法（DMGA），將僅采用探測(cè)變異策略的遺傳算法稱為探測(cè)變異遺傳算法（DEMGA），將僅采用直接變異策略的遺傳算法稱為直接變異遺傳算法（DIMGA）。

（1）探測(cè)變異操作

探測(cè)變異的具體操作是在隨機(jī)選擇1 個(gè)基因位后分5 種情況進(jìn)行。

情況1：若選擇的是符號(hào)位變異，則直接進(jìn)行變異探測(cè)，并計(jì)算適應(yīng)度值，判斷是否發(fā)生變異（圖2）。

圖2 探測(cè)變異操作情況1 示意圖Fig.2 Schematic diagram of the first case of detecting mutation operation

情況2：若選擇的基因位為變量xi的第1 個(gè)數(shù)字位，則直接進(jìn)行變異探測(cè)，并計(jì)算適應(yīng)度值，判斷是否發(fā)生變異。該變異探測(cè)分3 種子情況（圖3）。子情況1，若該位的值為0，則變異為[1,9]區(qū)間的隨機(jī)整數(shù)；子情況2，若該位的值為9，則變異為[0,8]區(qū)間的隨機(jī)整數(shù)；子情況3，若該位的值為[1,8]區(qū)間的整數(shù)m，則變異分兩個(gè)方向進(jìn)行探測(cè)，一個(gè)方向是變異為[m+1,9]區(qū)間的隨機(jī)整數(shù)，另一個(gè)方向是變異為[0,m?1]區(qū)間的隨機(jī)整數(shù)，然后將適應(yīng)度值較大者與變異前適應(yīng)度值進(jìn)行比較，判斷是否變異。

圖3 探測(cè)變異操作情況2 示意圖Fig.3 Schematic diagram of the second case of detecting mutation operation

情況3：若選擇的基因位位于變量xi第1 個(gè)數(shù)字位之后，且該位的值位于區(qū)間[1,8]，則直接進(jìn)行變異探測(cè)（同情況2 的子情況3），并計(jì)算適應(yīng)度值，判斷是否發(fā)生變異。

情況4：若選擇的基因位位于變量xi第1 個(gè)數(shù)字位之后，且該位的值為0，則該變異探測(cè)分2 種子情況（圖4）。子情況1，當(dāng)其前一位的值等于0 時(shí)，則該位變異為[1,9]區(qū)間的隨機(jī)整數(shù)；子情況2，當(dāng)其前一位的值大于0 時(shí)，則變異分兩個(gè)方向進(jìn)行探測(cè)，一個(gè)方向是該位變異為[1,9]區(qū)間的隨機(jī)整數(shù)m，另一個(gè)方向是該位值也變異為m，且前一位的值作減1 變異，然后將適應(yīng)度值較大者與變異前適應(yīng)度值進(jìn)行比較，判斷是否變異。

圖4 探測(cè)變異操作情況4 示意圖Fig.4 Schematic diagram of the fourth case of detecting mutation operation

情況5：若選擇的基因位位于變量xi第1 個(gè)數(shù)字位之后，且該位的值為9，則該變異探測(cè)分2 種子情況（圖5）。子情況1，當(dāng)其前一位的值等于9 時(shí)，則該位變異為[0,8]區(qū)間的隨機(jī)整數(shù)；子情況2，當(dāng)其前一位的值小于9 時(shí)，則變異分兩個(gè)方向進(jìn)行探測(cè)，一個(gè)方向是該位變異為[0,8]區(qū)間的隨機(jī)整數(shù)m，另一個(gè)方向是該位值也變異為m，且前一位的值作加1 變異，然后將適應(yīng)度值較大者與變異前適應(yīng)度值進(jìn)行比較，判斷是否變異。

圖5 探測(cè)變異操作情況5 示意圖Fig.5 Schematic diagram of the fifth case of detecting mutation operation

（2）直接變異操作

直接變異操作是指將待變異個(gè)體的每個(gè)變量xi所對(duì)應(yīng)的染色體等分為前后兩段，在每一段中隨機(jī)選擇1 個(gè)基因位進(jìn)行直接變異，變異后的適應(yīng)度值不控制變異的發(fā)生。直接變異能夠提供種群的多樣性，避免早熟的產(chǎn)生。

“昆南”陰平聲字“他”的唱調(diào)（《南西廂·佳期》【十二紅】“愛他兩意和諧”，687），其中的即為兩節(jié)型過腔。雖然這個(gè)過腔的音樂材料都相同，都來自于本唱調(diào)音階的級(jí)音，但由此組成的樂匯或句型可以分作兩節(jié)，為第一節(jié)級(jí)音性過腔，為第二節(jié)級(jí)音性過腔。這個(gè)過腔即是由同一種音樂材料組成的兩節(jié)型過腔。

（3）雙重變異操作

當(dāng)待變異個(gè)體適應(yīng)度值f接近當(dāng)前最佳個(gè)體時(shí)，希望待變異個(gè)體能夠具有局部尋優(yōu)能力，同時(shí)也具有跳出局部極值點(diǎn)，搜索全局最優(yōu)的能力；當(dāng)待變異個(gè)體適應(yīng)度值f距離當(dāng)前最佳個(gè)體較遠(yuǎn)時(shí)，希望個(gè)體向著適應(yīng)度值增加的方向變異。將探測(cè)變異與直接變異結(jié)合在一起形成雙重變異操作，便可以實(shí)現(xiàn)上述目標(biāo)。

雙重變異操作的實(shí)施策略：

（4）變異概率

變異操作是在一定的概率下發(fā)生的，當(dāng)產(chǎn)生的[0,1）區(qū)間的隨機(jī)數(shù)小于變異概率Pm時(shí)，將會(huì)發(fā)生變異。變異概率的大小決定著新個(gè)體產(chǎn)生的機(jī)會(huì)及計(jì)算量的大小。較小的變異概率，產(chǎn)生新個(gè)體的機(jī)會(huì)較?。惠^大的變異概率，遺傳算法的計(jì)算量較大，且趨于純粹的隨機(jī)搜索算法。

為了能達(dá)到好的尋優(yōu)效果，在進(jìn)化的早期期望能有較大的變異概率，以增加搜索的廣度；在進(jìn)化的晚期期望能有較小的變異概率，以提高進(jìn)化的速度。同時(shí)，也希望當(dāng)待變異個(gè)體適應(yīng)度值接近當(dāng)前最佳個(gè)體時(shí)，減小變異概率，使較好的個(gè)體得到保留；當(dāng)待變異個(gè)體適應(yīng)度值距離當(dāng)前最佳個(gè)體較遠(yuǎn)時(shí)，有較大的變異概率，使較差的個(gè)體容易被淘汰。

基于上述原因，本文在綜合考慮待變異個(gè)體的適應(yīng)度值f及進(jìn)化代數(shù)對(duì)變異概率影響的情況下，提出一種復(fù)合自適應(yīng)變異概率Pm，其計(jì)算公式如下：

其中，a為區(qū)間（0,1）內(nèi)的常數(shù)，a越大，變異概率Pm越大；k為區(qū)間[0,1]內(nèi)的常數(shù)，k越大，變異概率Pm也越大。反映了待變異個(gè)體的適應(yīng)度值距離當(dāng)前種群中的最大適應(yīng)度值的相對(duì)距離；當(dāng)一定時(shí)，Pm將隨t的增大而減小，即進(jìn)化早期變異概率較大，使算法有較大的搜索廣度；當(dāng)t一定時(shí)，越大，Pm就越小，即待變異個(gè)體適應(yīng)度值越接近當(dāng)前最佳個(gè)體，較好的個(gè)體就更容易得到保留。

2.4 雙重變異遺傳算法的計(jì)算步驟

步驟一：隨機(jī)生成實(shí)數(shù)編碼（圖1）的初始種群，初始種群中的個(gè)體數(shù)為N，且N為偶數(shù)。

步驟二：將個(gè)體染色體解碼為實(shí)數(shù)，根據(jù)式（6）計(jì)算初始種群中每個(gè)個(gè)體的適應(yīng)度值，并令當(dāng)前進(jìn)化代數(shù)t=0。

步驟三：根據(jù)式（8）計(jì)算選擇概率，并令j=0。

步驟四：如果j=0，將種群中適應(yīng)度值最大的一個(gè)個(gè)體直接選出，不進(jìn)行后續(xù)的交叉及變異操作，直接復(fù)制到下一代，然后根據(jù)式（7）采用轉(zhuǎn)盤式選擇方式選擇1 個(gè)個(gè)體，進(jìn)行步驟六；如果j>0，根據(jù)式（7）采用轉(zhuǎn)盤式選擇方式選擇2 個(gè)個(gè)體，進(jìn)行步驟五。

步驟五：將選出的2 個(gè)個(gè)體進(jìn)行單變量交叉操作，進(jìn)行步驟六。交叉操作方法：當(dāng)產(chǎn)生的[0,1）區(qū)間的隨機(jī)數(shù)小于交叉概率Pc時(shí)，隨機(jī)選擇一個(gè)變量xi，使個(gè)體1 中變量xi所對(duì)應(yīng)的基因位的值與個(gè)體2 的進(jìn)行交換。其中，交叉概率Pc根據(jù)式（9）計(jì)算得到。

步驟六：對(duì)個(gè)體進(jìn)行雙重變異操作，進(jìn)行步驟七。雙重變異操作方法：當(dāng)時(shí)，產(chǎn)生一個(gè)[0,1]的隨機(jī)數(shù) γ，若γ>0.5，采用探測(cè)變異；若γ ?0.5，采用直接變異。當(dāng)時(shí)，采用探測(cè)變異。探測(cè)變異操作、直接變異操作按2.3.3 節(jié)所述方法進(jìn)行。變異概率Pm按式（10）計(jì)算。

步驟七：如果（j+1）×2

步驟八：將新生成的種群個(gè)體的染色體解碼為實(shí)數(shù)，根據(jù)式（6）計(jì)算適應(yīng)度值。如果t==T，計(jì)算結(jié)束；如果t

說明：本文所述的直接變異遺傳算法、探測(cè)變異遺傳算法與雙重變異遺傳算法的計(jì)算步驟的不同之處主要體現(xiàn)在步驟六。若將步驟六改為：“對(duì)個(gè)體進(jìn)行直接變異操作”，則為直接變異遺傳算法的計(jì)算步驟；若將步驟六改為：“對(duì)個(gè)體進(jìn)行探測(cè)變異操作”，則為探測(cè)變異遺傳算法的計(jì)算步驟。

3 澳大利亞計(jì)算機(jī)應(yīng)用協(xié)會(huì)考核題分析

本文將雙重變異遺傳算法（DMGA）與簡(jiǎn)化Bishop法相結(jié)合編寫Python 程序，以1987年澳大利亞計(jì)算機(jī)應(yīng)用協(xié)會(huì)（ACADS）設(shè)計(jì)的考核題1（a）、考核題1（c）[18]為求解對(duì)象，對(duì)這兩道題分別進(jìn)行50 次獨(dú)立計(jì)算，并以ACADS 提供的裁判程序答案檢驗(yàn)算法的尋優(yōu)效果。然后，在分別僅進(jìn)行直接變異、探測(cè)變異的情況下，對(duì)上述兩題也進(jìn)行50 次獨(dú)立計(jì)算，對(duì)比分析DMGA 算法的性能。

3.1 考核題及其裁判程序答案

考核題1（a）為一均質(zhì)土質(zhì)邊坡，如圖6所示，坡高10 m，坡度為26.565°，黏聚力為3 kPa，內(nèi)摩擦角為19.6°，土體重度為20 kN/m3。要求確定邊坡的臨界滑動(dòng)面及最小安全系數(shù)。

圖6 考核題1（a）坡面示意圖（坐標(biāo)平移后）[18]Fig.6 Profile of slope in EX1（a）（after coordinates are modified）[18]

考核題1（c）為一非均質(zhì)土質(zhì)邊坡，如圖7所示，坡高10 m，坡度為26.565°，黏聚力、內(nèi)摩擦角及土體重度見表1。要求確定邊坡的臨界滑動(dòng)面及最小安全系數(shù)。

圖7 考核題1（c）剖面圖（坐標(biāo)平移后）[18]Fig.7 Profile of slope in EX1（c）（after coordinates are modified）[18]

表1 考核題1（c）的材料性質(zhì)[18]Table 1 Material characteristic in EX1（c）[18]

ACADS 給出的考核題1（a）邊坡最小安全系數(shù)的裁判程序答案有：0.990，0.991，1.000；考核題1（c）邊坡最小安全系數(shù)的裁判程序答案有：1.385，1.390，1.406。

3.2 計(jì)算參數(shù)

本文給定圓弧滑動(dòng)面半徑R、圓心橫坐標(biāo)xo、圓心縱坐標(biāo)yo的搜索區(qū)間分別為：R∈（0,99]，xo∈

[?99,99]，yo∈[0,99]。

計(jì)算所涉及的相關(guān)參數(shù)取值如下：

（a）滑體土條寬度按R/50 進(jìn)行取值，且單個(gè)條塊不跨越地形剖面的轉(zhuǎn)折點(diǎn)。當(dāng)含有地形轉(zhuǎn)折點(diǎn)、滑動(dòng)面剪出口及滑動(dòng)面后緣的條塊，其寬度不足R/50 時(shí)，上述相應(yīng)控制點(diǎn)即為條塊端點(diǎn)；

（b）變異概率計(jì)算公式（10）中參數(shù)a=0.5，k=0.6；

（c）實(shí)數(shù)編碼的串長(zhǎng)為18（每個(gè)變量的編碼長(zhǎng)度為6，且變量精確到3 位小數(shù)），如圖1所示；

（d）種群數(shù)N 為50；

（e）初次最大進(jìn)化代數(shù)T為200，當(dāng)最大進(jìn)化代數(shù)與最佳個(gè)體出現(xiàn)代數(shù)之差小于等于30 時(shí)，最大進(jìn)化代數(shù)增加30 代，最大進(jìn)化代數(shù)的最大值為500。

3.3 計(jì)算結(jié)果及過程分析

3.3.1 雙重變異遺傳算法計(jì)算結(jié)果

表2、表3 分別是基于雙重變異遺傳算法對(duì)考核題1（a）、1（c）的計(jì)算結(jié)果（限于篇幅，僅列出了部分計(jì)算結(jié)果）。

表2 考核題1（a）的計(jì)算結(jié)果（DMGA）Table 2 Calculated results of EX1（a）（DMGA）

表3 考核題1（c）的計(jì)算結(jié)果（DMGA）Table 3 Calculated results of EX1（c）（DMGA）

計(jì)算結(jié)果統(tǒng)計(jì)顯示，對(duì)于考核題1（a），50 次計(jì)算中，安全系數(shù)的最小值為0.985 2，最大值為0.993 7，平均值為0.985 7，與該考核題的裁判程序答案基本一致，說明對(duì)于均質(zhì)邊坡該算法可靠有效；對(duì)于考核題1（c），50 次計(jì)算中，安全系數(shù)的最小值為1.394 9，最大值為1.407 0，平均值為1.397 9，與該考核題推薦的裁判答案基本一致，說明對(duì)于非均質(zhì)邊坡該算法也可靠有效。

3.3.2 三種算法計(jì)算結(jié)果比較

基于雙重變異遺傳算法（DMGA）、直接變異遺傳算法（DIMGA）、探測(cè)變異遺傳算法（DEMGA），分別對(duì)兩道考核題計(jì)算50 次的統(tǒng)計(jì)結(jié)果（表4）顯示，DMGA算法搜索得到的安全系數(shù)平均值小于DIMGA 算法、DEMGA 算法的計(jì)算結(jié)果，說明DMGA 算法具有更強(qiáng)的全局搜索能力。并且DMGA 算法計(jì)算所得的標(biāo)準(zhǔn)差較小，說明該算法具有較好的魯棒性。

表4 考核題1（a）、1（c）的計(jì)算結(jié)果統(tǒng)計(jì)分析Table 4 Statistical analysis of computation results of EX1（a）and EX1（c）

3.3.3 三種算法計(jì)算過程對(duì)比

（1）單次計(jì)算的進(jìn)化過程對(duì)比

圖8 為采用DMGA 算法、DIMGA 算法、DEMGA算法對(duì)考核題1（a）、1（c）進(jìn)行的第1 次計(jì)算的進(jìn)化過程曲線；表5 為采用DMGA 算法、DIMGA 算法、DEMGA算法對(duì)考核題1（a）、1（c）進(jìn)行的第1 次計(jì)算的收斂過程對(duì)比表。圖8 顯示DEMGA 算法在計(jì)算考核題1（a）時(shí)，最大適應(yīng)度值曲線明顯偏低，表明陷入了局部極值點(diǎn)；DIMGA 算法計(jì)算的最大適應(yīng)度值曲線在進(jìn)化早期明顯偏低，且對(duì)考核題1（a）、1（c）分別在第116代和第151 代才收斂（表5），表明其收斂速度較慢；DMGA 算法分別在第66 代和第72 代時(shí)收斂，其計(jì)算得到的安全系數(shù)低于DIMGA 算法和DEMGA 算法，顯示出優(yōu)于DIMGA 算法和DEMGA 算法的尋優(yōu)能力。

表5 考核題1（a）、1（c）算法收斂過程對(duì)比表Table 5 Comparison of convergence processes of EX1（a）and EX1（c）

圖8 考核題1（a）、1（c）第1 次計(jì)算的進(jìn)化過程曲線Fig.8 Evolution curve of the first calculation of EX1（a）and EX1（c）

（2）50 次計(jì)算的平均進(jìn)化過程對(duì)比

為了進(jìn)一步對(duì)比三種算法進(jìn)化過程的統(tǒng)計(jì)變化規(guī)律，根據(jù)各算法50 次計(jì)算的歷代最佳個(gè)體的適應(yīng)度值計(jì)算出歷代最佳個(gè)體的平均適應(yīng)度值，繪制出前200 代的平均進(jìn)化過程曲線（圖9）。曲線顯示，與上述單次計(jì)算對(duì)比類似，DEMGA 算法的局部尋優(yōu)能力較強(qiáng)，在進(jìn)化的早期，能快速搜索到局部極值點(diǎn)附近，但其跳出局部極值點(diǎn)的能力弱，容易出現(xiàn)早熟現(xiàn)象；DIMGA 算法具有全局搜索的特點(diǎn)，局部尋優(yōu)能力較弱，收斂速度較慢；DMGA 算法將探測(cè)變異與直接變異結(jié)合，既具有較強(qiáng)的局部搜索能力，又具有較強(qiáng)的跳出局部極值點(diǎn)的能力，能夠在搜索的廣度與深度上達(dá)到較好平衡，因此其全局尋優(yōu)能力強(qiáng)。

圖9 考核題1（a）、1（c）的平均進(jìn)化過程曲線Fig.9 Average evolution curve of EX1（a）and EX1（c）

4 工程實(shí)例

某海堤邊坡[19?20] 為一非均質(zhì)土質(zhì)邊坡，其剖面及土層物理力學(xué)參數(shù)如圖10所示，現(xiàn)采用雙重變異遺傳算法結(jié)合簡(jiǎn)化Bishop 法搜索該邊坡的臨界滑動(dòng)面及最小安全系數(shù)。

圖10 邊坡剖面圖[19?20]Fig.10 Slope profile[19?20]

圓弧滑動(dòng)面半徑及圓心坐標(biāo)的搜索區(qū)間分別為：R∈（0,99]，xo∈[?99,99]，yo∈[0,99]。且需滿足約束條件[20]：6≤yo≤30。式（10）中的參數(shù)k=0.7，其他參數(shù)設(shè)置與3.2 節(jié)相同。

經(jīng)過1 次獨(dú)立計(jì)算，根據(jù)歷代的適應(yīng)度值（最大值、平均值）與進(jìn)化代數(shù)做出進(jìn)化過程曲線（圖11），該曲線顯示了歷代進(jìn)化的個(gè)體適應(yīng)度值的變化情況。在進(jìn)化的初期（第7 代之前），種群的最大適應(yīng)度值接近于0，說明進(jìn)化初期種群中的個(gè)體均位于非可行域內(nèi)。第7 代之后，遺傳進(jìn)化使種群中的一部分個(gè)體進(jìn)入可行域，種群的最大適應(yīng)度值及平均適應(yīng)度值總體上隨進(jìn)化代數(shù)的增加而增加。當(dāng)進(jìn)化到第83 代時(shí)，種群的最大適應(yīng)度值為0.370 89，隨后直到第200代進(jìn)化結(jié)束時(shí)均未更新，算法收斂。

圖11 進(jìn)化過程曲線Fig.11 Evolutionary process curve

計(jì)算結(jié)果顯示，邊坡的安全系數(shù)F=1.696 2，臨界滑動(dòng)面的圓心坐標(biāo)xo=6.859 m，yo=11.910 m，半徑R=14.910 m，臨界滑動(dòng)面左側(cè)與坡面交點(diǎn)的橫坐標(biāo)為xL=?2.11 m，臨界滑動(dòng)面右側(cè)與坡面交點(diǎn)的橫坐標(biāo)為xR=20.55 m。與文獻(xiàn)[20]采用簡(jiǎn)化Bishop 法計(jì)算的滑面位置（xL=?2.29 m,xR=20.17 m,yo=11.65 m）基本一致，且安全系數(shù)小于文獻(xiàn)[20]計(jì)算的最小值（Fmin=1.727），說明算法可靠有效，且尋優(yōu)效果更佳。

5 結(jié)論

（1）本文提出的雙重變異遺傳算法，能夠在搜索的廣度與深度上達(dá)到較好平衡，使待變異個(gè)體適應(yīng)度值距離當(dāng)前最佳個(gè)體較遠(yuǎn)時(shí)，能夠向著適應(yīng)度值增加的方向進(jìn)行搜索；待變異個(gè)體適應(yīng)度值接近當(dāng)前最佳個(gè)體時(shí)，算法既具有局部尋優(yōu)能力，也具有跳出局部極值點(diǎn)，搜索全局最優(yōu)的能力。

（2）與僅進(jìn)行直接變異或探測(cè)變異的遺傳算法相比，雙重變異遺傳算法具有更強(qiáng)的全局搜索能力及魯棒性，與簡(jiǎn)化Bishop 法結(jié)合，能夠有效地搜索出邊坡的圓弧形臨界滑動(dòng)面。