基于粒子群-信賴域的金屬極化曲線擬合算法

孫 峰 孫 偉

(上海海事大學(xué)信息工程學(xué)院 上海 201306)

0 引 言

在由電極極化引起的電極反應(yīng)中，電極電位與電流密度之間的關(guān)系稱為極化曲線。動(dòng)電位極化曲線是腐蝕研究中有用的電化學(xué)測量技術(shù)，它可以提供關(guān)于腐蝕機(jī)理、動(dòng)力學(xué)參數(shù)和腐蝕速率的重要信息。其中腐蝕動(dòng)力學(xué)參數(shù)是機(jī)理模型的重要組成部分，它們的值可能會(huì)極大地影響建模結(jié)果的可靠性和準(zhǔn)確性，因此極化曲線的參數(shù)求解至關(guān)重要。在對金屬材料進(jìn)行仿真建模時(shí)，一個(gè)必須要考慮的問題就是對極化曲線的準(zhǔn)確描述，而極化曲線是非線性函數(shù)，故必須提高對極化曲線的擬合精度才能保證仿真模型的準(zhǔn)確性和可行性。

目前，針對如何解決極化曲線擬合的問題，國內(nèi)外學(xué)者已經(jīng)提出多種方法。曹楚南[4]利用高斯-牛頓曲線擬合法直接對腐蝕金屬電極的弱極化曲線擬合,從而解決這一問題，但是這種方法對參數(shù)初始值要求過高，初始值的好壞對結(jié)果有著極大的影響。周明華等[5]使用遺傳算法來求解曲線擬合問題，雖然無需設(shè)置參數(shù)初始值，但是計(jì)算耗時(shí)，并且擬合精度也相對較低。文獻(xiàn)[6-7]通過線性擬合極化曲線的塔菲爾區(qū)域來推導(dǎo)陽極和陰極塔菲爾斜率，并且通過使用塔菲爾外推法獲得腐蝕速率，但塔菲爾外推法僅適用于塔菲爾線性區(qū)域明顯且數(shù)據(jù)能夠獲得的反應(yīng)。李強(qiáng)[8]在文獻(xiàn)[6-7]基礎(chǔ)上，對于塔菲爾線性區(qū)域不明顯或不存在的非線性極化曲線，提出采用非線性最小二乘法來擬合曲線，常用的非線性最小二乘法為Levenberg-Marquardt算法與信賴域算法，雖然此方法解決了塔菲爾線性區(qū)域數(shù)據(jù)難以測量的問題，但是同樣面臨著對參數(shù)初始值要求過高的問題。

以上提及的極化曲線擬合方法存在擬合精度不高、對數(shù)據(jù)和參數(shù)初始值依賴度高、計(jì)算十分耗時(shí)和需要的塔菲爾線性區(qū)數(shù)據(jù)難以測量等問題。在李強(qiáng)的研究基礎(chǔ)上，本文受智能優(yōu)化算法的啟發(fā)，結(jié)合數(shù)值計(jì)算算法和智能優(yōu)化算法兩類算法的優(yōu)點(diǎn)，提出了一種基于粒子群-信賴域的金屬腐蝕極化曲線擬合算法，來解決對于塔菲爾線性區(qū)不明顯或不存在的非線性極化曲線。算法采用粒子群算法來優(yōu)化信賴域算法的參數(shù)，這樣可以避免人為選擇參數(shù)的盲目性，提高算法對極化曲線的擬合精度，并保證了每次擬合結(jié)果的魯棒性。

1 PSO-TR算法

1.1 算法思想

PSO-TR算法的基本思想：將數(shù)值計(jì)算算法和智能優(yōu)化算法相結(jié)合，綜合考量兩類算法的優(yōu)缺點(diǎn)。例如：粒子群算法的優(yōu)點(diǎn)是無需參數(shù)初始值，擬合精度高，缺點(diǎn)是結(jié)果不穩(wěn)定，每次結(jié)果不同，并且優(yōu)化耗時(shí)；而信賴域算法的優(yōu)點(diǎn)是結(jié)果穩(wěn)定，計(jì)算速度快，缺點(diǎn)是計(jì)算需要參數(shù)初始值，并且十分依賴參數(shù)初始值的好壞。發(fā)現(xiàn)兩者的優(yōu)缺點(diǎn)恰好可以互相彌補(bǔ)。故結(jié)合兩種算法，先使用粒子群算法來優(yōu)化參數(shù)初始值，再將得到的最優(yōu)參數(shù)初始值帶入信賴域算法中計(jì)算，這樣可以避免人為選擇參數(shù)的盲目性，提高算法對極化曲線的擬合精度，節(jié)約計(jì)算時(shí)間。PSO-TR算法主要分為兩步，第一步為PSO優(yōu)化TR參數(shù)過程，使用PSO算法求解極化曲線擬合問題，但是由于PSO算法求解此問題耗時(shí)和結(jié)果不穩(wěn)定，故而限制迭代次數(shù)，在短時(shí)間內(nèi)求得一個(gè)較優(yōu)的參數(shù)結(jié)果，并將其作為下一步的初始參數(shù)。第二步為TR求解參數(shù)過程，使用TR算法求解極化曲線擬合問題，由于算法依賴參數(shù)初始值，所以將第一步求得的優(yōu)化初始參數(shù)作為TR算法的參數(shù)初始值，以此來求解最終的最優(yōu)參數(shù)。

1.2 PSO優(yōu)化TR參數(shù)過程

PSO是Kennedy和Eberhart于1995年提出的一種全局優(yōu)化算法，是對于鳥群覓食過程中的遷徙和聚集的模擬[9]。粒子群算法作為智能進(jìn)化算法，其擁有搜索能力強(qiáng)、收斂速度快、設(shè)置參數(shù)少、程序?qū)崿F(xiàn)簡單和無梯度信息等優(yōu)點(diǎn)[10]，所有粒子都具有由目標(biāo)函數(shù)所確定的適應(yīng)值和確定其運(yùn)動(dòng)方向和運(yùn)動(dòng)距離的速度，該問題的解就是搜索空間中的一個(gè)粒子的位置。

本方法首先將極化曲線擬合問題看作是PSO目標(biāo)優(yōu)化問題，將每個(gè)參數(shù)的取值范圍作為粒子的一個(gè)搜索維度，將參數(shù)值定義為粒子的位置信息，因此可以看作是粒子在一個(gè)多維的空間內(nèi)尋找最優(yōu)值的過程。

PSO算法的主要思想如下[11-14]:在n維搜索空間中，有m個(gè)粒子，其中第i個(gè)粒子的位置是Xi=(xi1，xi2，…，xin)，其速度為Vi=(vi1，vi2，…，vin) ，i=1，2，…，m。第i個(gè)粒子搜索到的最優(yōu)位置(即個(gè)體最優(yōu)值Pbest)為Pi=(pi1，pi2，…，pin) ; 整個(gè)粒子群搜索到的最優(yōu)位置(全局最優(yōu)值Gbest)為Pg=(pg1，pg2，…，pgn)。在每次迭代期間，粒子通過對Pbest和Gbest兩個(gè)極值的比較來更新自己的狀態(tài)，每個(gè)粒子根據(jù)下式來更新自己的狀態(tài)：

(1)

(2)

PSO-TR算法的參數(shù)優(yōu)化部分是利用PSO算法獲得最優(yōu)初始參數(shù)，將每個(gè)參數(shù)的取值范圍作為粒子的一個(gè)搜索維度，將其定義為粒子的位置信息，因此可以看作是粒子在一個(gè)多維的空間內(nèi)尋找最優(yōu)值的過程。構(gòu)建樣本誤差平方和F作為PSO-TR算法的目標(biāo)函數(shù)，其表達(dá)式為：

(3)

1.3 信賴域算法求解參數(shù)過程

信賴域算法(Trust Region，TR)是一種用于求解非線性優(yōu)化問題的數(shù)值計(jì)算方法。信賴域算法是一種迭代計(jì)算算法，從給定的初始解開始，不斷迭代，持續(xù)改進(jìn)，直到獲得令人滿意的近似最優(yōu)解。在每次迭代期間給出信賴域，其一般是當(dāng)前迭代點(diǎn)xk的小鄰域。然后解決該鄰域內(nèi)的子問題，獲得試探步長sk，然后使用評價(jià)函數(shù)來決定是否接受該試探步并確定迭代的信賴域。如果試探步長被接受，那么xk+1=xk+sk，否則，xk+1=xk。新的信賴域的大小取決于試探步的好壞，簡而言之，當(dāng)試探步比較好時(shí)，信賴域在下一步中擴(kuò)大或保持不變，否則信賴域減小[15-16]。

信賴域算法的基本思想：設(shè)當(dāng)前點(diǎn)xk的鄰域定義為：Ωk={x∈Rn|‖x-xk‖≤Δk}，其中，Δk是信賴域半徑。目標(biāo)函數(shù)近似為極值點(diǎn)附近的一個(gè)二次函數(shù)，故對于無約束優(yōu)化問題，由二次近似構(gòu)造的信賴域子問題如下：

(4)

式中：s=x-xk；gk為目標(biāo)函數(shù)f(x)在當(dāng)前迭代點(diǎn)xk處的梯度；Bk∈Rn×n，是f(x)在xk處Hesse陣▽2f(xk)或者其近似。

因此構(gòu)建樣本誤差平方和f(x)作為PSO-TR算法的目標(biāo)函數(shù)，其表達(dá)式為：

(5)

設(shè)sk為信賴域子問題的解。其中目標(biāo)函數(shù)f(x)在第k步的實(shí)際下降量，即真實(shí)下降量為：

Aredk=f(xk)-f(xk+sk)

(6)

二次模型函數(shù)q(k)(s)的下降量，即預(yù)測下降量為：

Predk=q(k)(0)-q(k)(sk)

(7)

定義比值rk為：

rk=Aredk/Predk

(8)

它用于測量二次模型和目標(biāo)函數(shù)的近似程度，rk越接近1，接近度越好。所以此數(shù)量也用來決定下一次迭代的信賴域半徑。

TR求解參數(shù)過程如下：

Step1將PSO優(yōu)化的參數(shù)初始值作為初始點(diǎn)x0，初始化信賴域半徑Δ0，并開始迭代。

Step2當(dāng)?shù)趉步時(shí)，計(jì)算梯度gk與Hesse陣Bk。

Step3求解信賴域模型，獲得試探步長sk，并計(jì)算比值rk。

Step5若rk≥0.75且‖sk‖=Δk，則表明這一步已經(jīng)到達(dá)了信賴域半徑的邊緣，同時(shí)步子有點(diǎn)小，可以嘗試擴(kuò)大信賴域半徑，令Δk+1=2Δk。

Step6若0.25

Step7若rk<0，則表明函數(shù)值向著上升而非下降的趨勢變化了(與最優(yōu)化的目標(biāo)相反)，表明這一步錯(cuò)得“離譜”了，此時(shí)不應(yīng)走到下一點(diǎn)，而應(yīng)“原地踏步”，即xk+1=xk，同時(shí)和上面rk≤0.25的情況一樣，縮小信賴域。反之，在rk>0的情況下，都可以到達(dá)下一點(diǎn)，即xk+1=xk+sk。

對于信賴域算法，參數(shù)初始值對算法的性能有極大影響，人為地選擇參數(shù)十分不合理，即使專家經(jīng)驗(yàn)也很難使結(jié)果最優(yōu)。因此PSO-TR算法先使用PSO算法優(yōu)化參數(shù)初始值，再將優(yōu)化的初始參數(shù)代入TR算法中計(jì)算，便可獲得最優(yōu)解。PSO-TR算法流程圖如圖1所示。

圖1 PSO-TR算法流程圖

2 實(shí) 驗(yàn)

2.1 實(shí)驗(yàn)數(shù)據(jù)集

本文中使用的實(shí)驗(yàn)數(shù)據(jù)集為某型汽車車門部件在常溫下的電化學(xué)腐蝕實(shí)驗(yàn)數(shù)據(jù)，分別為“點(diǎn)焊鋼板”、“Q235鋼”和“鍍鋅鋼板”三種材料在PH值為7的0.1 mol/L的NaCl溶液環(huán)境中通過電化學(xué)腐蝕裝置得到的真實(shí)數(shù)據(jù)集。

2.2 實(shí)驗(yàn)內(nèi)容

本文采用文獻(xiàn)[8]所提供的金屬在大氣環(huán)境下的電化學(xué)腐蝕機(jī)理和理論模型。為了計(jì)算理論電流密度，使用基于混合電位理論的疊加模型方法，總的電流密度通過反應(yīng)中涉及的所有電流密度的總和來計(jì)算。對于含氧水溶液中的金屬腐蝕，通過式(9)-式(11)計(jì)算每個(gè)反應(yīng)電流，并且可以通過式(12)計(jì)算總的電流密度i。

(9)

(10)

(11)

i=iMe-iO2-iH2O

(12)

式中:iMe、iO2和iH2O分別為金屬氧化反應(yīng)、氧還原反應(yīng)和水還原反應(yīng)的電流密度;bMe、bO2和bH2O分別為金屬氧化反應(yīng)、氧還原反應(yīng)和水還原反應(yīng)的塔菲爾斜率;iL，O2是溶解氧的限制擴(kuò)散電流密度; 超電勢η=E-Erev，Erev是可逆的潛力。iMe、iO2、iH2O、bMe、bO2、bH2O和iL，O2便是所需要求得的7個(gè)腐蝕動(dòng)力學(xué)參數(shù)。

2.3 實(shí)驗(yàn)過程

本文提出的基于粒子群-信賴域的金屬極化曲線擬合算法的流程如下：

Step1輸入所需實(shí)驗(yàn)樣本數(shù)據(jù)。

Step2確定待擬合的極化曲線表達(dá)式。依據(jù)金屬腐蝕的極化行為特征確定極化曲線表達(dá)式。

Step3定義粒子位置信息。將每個(gè)參數(shù)的取值范圍作為一個(gè)維度，因此可以看作是粒子在一個(gè)多維的空間內(nèi)尋找全局最優(yōu)值的過程。

Step4定義目標(biāo)函數(shù)表達(dá)式。構(gòu)建樣本誤差平方和F作為PSO-TR算法的目標(biāo)函數(shù)，其表達(dá)式為：

(13)

Step5設(shè)置PSO-TR算法的初始參數(shù)，如最大迭代次數(shù)k、種群規(guī)模m等參數(shù)。

Step6執(zhí)行PSO-TR算法，獲得最優(yōu)參數(shù)解。

2.4 實(shí)驗(yàn)結(jié)果及分析

本文從PSO-TR算法的有效性和效率性兩個(gè)方面進(jìn)行了實(shí)驗(yàn)驗(yàn)證。

(1) PSO-TR算法在不同金屬材料下的極化曲線擬合實(shí)驗(yàn)及分析 為驗(yàn)證本文提出方法的有效性，針對以上三種金屬材料的實(shí)驗(yàn)數(shù)據(jù)集，采用PSO-TR算法進(jìn)行實(shí)驗(yàn)驗(yàn)證，并將各自的擬合曲線與實(shí)驗(yàn)數(shù)據(jù)進(jìn)行對比。其中7個(gè)參數(shù)的取值上下限為([1e-4, 30, 1e-12, 30, 0, 1e-8, 30], [1e-2, 500, 1e-5, 500, 1, 1e-3, 500]) ，三種材料的PSO-TR算法擬合曲線與實(shí)驗(yàn)數(shù)據(jù)對比如圖2-圖4所示。顯然，從圖2-圖4中也可以看出PSO-TR算法擬合的極化曲線和實(shí)際測量的極化曲線偏差較小，可較好地?cái)M合實(shí)驗(yàn)數(shù)據(jù)。

圖2 點(diǎn)焊鋼板的PSO-TR擬合曲線與實(shí)驗(yàn)數(shù)據(jù)圖

圖3 Q235鋼的PSO-TR擬合曲線與實(shí)驗(yàn)數(shù)據(jù)圖

圖4 鍍鋅鋼板的PSO-TR擬合曲線與實(shí)驗(yàn)數(shù)據(jù)圖

將總的極化曲線進(jìn)行分解，便可得到單個(gè)反應(yīng)對極化曲線的貢獻(xiàn)，如圖5-圖7所示，可知金屬的氧化反應(yīng)和氧的還原反應(yīng)對總極化曲線有很大貢獻(xiàn)，且水的還原反應(yīng)對總極化曲線貢獻(xiàn)很小。

圖5 點(diǎn)焊鋼板中單個(gè)反應(yīng)對極化曲線的貢獻(xiàn)圖

圖6 Q235鋼中單個(gè)反應(yīng)對極化曲線的貢獻(xiàn)圖

圖7 鍍鋅鋼板中單個(gè)反應(yīng)對極化曲線的貢獻(xiàn)圖

(2) PSO-TR算法對比實(shí)驗(yàn)及分析 為驗(yàn)證本文提出方法的效率，針對以上三組實(shí)驗(yàn)數(shù)據(jù)集，分別使用Levenberg-Marquardt(LM)算法、信賴域算法、PSO算法和本文提出的PSO-TR算法進(jìn)行實(shí)驗(yàn)對比，其中7個(gè)參數(shù)的取值上下限為([1e-4, 30, 1e-12, 30, 0, 1e-8, 30], [1e-2, 500, 1e-5, 500, 1, 1e-3, 500]) ，同時(shí)Levenberg-Marquardt(LM)算法和信賴域算法所需要的參數(shù)初始值為[3.1e-3, 128, 3e-6, 243, 3.3e-2, 5.05e-5, 126]，實(shí)驗(yàn)得到的擬合精度結(jié)果如表1所示。

表1 不同數(shù)據(jù)集擬合精度

從表1中可以看出，參數(shù)初始值對于LM算法與TR算法的性能有著極大影響，同一組參數(shù)初始值在三組數(shù)據(jù)的測試中，使用LM算法和TR算法時(shí)，結(jié)果時(shí)好時(shí)壞，并且三種金屬材料的實(shí)驗(yàn)結(jié)果中PSO-TR算法的擬合精度比LM算法和TR算法的擬合精度都要高。同時(shí)可以看出PSO算法的擬合精度與PSO-TR算法的擬合精度十分接近，但這是在PSO算法迭代時(shí)間較長的情況下才有的結(jié)果，并且每次的結(jié)果都不相同，十分不穩(wěn)定。顯然，本文提出的PSO-TR算法相比于TR、LM和PSO三種算法，優(yōu)點(diǎn)較為明顯，無需設(shè)置參數(shù)初始值，結(jié)果穩(wěn)定，擬合精度高，迭代次數(shù)少、速度較快。

3 結(jié) 語

本文在針對強(qiáng)極化曲線擬合問題時(shí)，提出了一種基于粒子群-信賴域的金屬腐蝕極化曲線擬合算法。PSO-TR算法采用粒子群算法來優(yōu)化信賴域算法的參數(shù)，先用粒子群算法求的最優(yōu)初始參數(shù)，再將求得的最優(yōu)初始參數(shù)帶入信賴域算法中求得最終結(jié)果，這樣可以避免人為選擇參數(shù)的盲目性，提高算法對極化曲線的擬合精度。實(shí)驗(yàn)證明該方法可以有效地解決傳統(tǒng)的非線性曲線擬合方法中的缺陷。首先無需設(shè)置參數(shù)初始值，這樣避免了傳統(tǒng)方法中大量調(diào)參的過程或者對相關(guān)專家經(jīng)驗(yàn)的依賴；其次這種方法的曲線擬合精度較高，不會(huì)陷入局部最優(yōu)，并且每次的擬合結(jié)果十分穩(wěn)定；最后解決了傳統(tǒng)方法中所需的塔菲爾線性區(qū)數(shù)據(jù)難以測量的問題，這樣在接下來的有限元仿真建模中會(huì)有效地提高建模結(jié)果的可靠性和準(zhǔn)確性。