Bouc-Wen 模型參數識別的非線性自適應遺傳算法和試驗驗證

章紅梅，胡 帆，段元鋒

(浙江大學建筑工程學院，浙江，杭州 310058)

結構動力行為中滯回模型是表征結構動力響應行為的重要規則[1]。Bouc-Wen 模型是一種通用的非線性光滑滯回模型，可以表達結構構件滯回特征中的剛度退化、強度退化等多種力學特征，在表征和模擬結構復雜滯回行為上具有強大的適應性[2]。然而，Bouc-Wen 模型本身參數眾多，且其中的多數參數不具有明確的物理意義，這使得獲得恰當的參數組較為困難，Bouc-Wen 模型的應用因此也受到極大的阻礙。針對該問題，Ismail等[3]直接用解析的方法對模型進行參數求解，但是其中的大量假設，使得所得到的Bouc-Wen 模型失去了廣泛性。相比解析方法，利用數值方法從試驗中識別Bouc-Wen 模型參數具有更廣泛的適用性。

目前，從已知試驗數據中識別Bouc-Wen 模型參數的算法研究大致可以分為兩類：一類為傳統的數值算法，如Kalman 濾波器法[4-6]、最小二乘法[7]、梯度下降法[8]等，對于模型的參數識別效果不算優異；另一類是各種智能算法，如遺傳算法[9]、粒子群算法[10]、微分進化算法[11]、GSO 算法[12]等，該類算法雖然可以達到較好的效果，但是實現過程較復雜，還有巨大的改進空間。

其中應用較多的遺傳算法是1975 年由Holland[13]教授首次提出的，目前已廣泛應用于遺傳合成、VLSI 技術、土木工程和機器學習。劉永強等[14]采用遺傳算法對MR 阻尼器進行了參數識別，通過縮小參數取值范圍的方法取得了較好的擬合效果；李得民等[15]通過比較由遺傳算法識別MRE隔振器的Bouc-Wen 模型和遺傳算法優化對MRE力學行為仿真建模的BP 神經網絡兩者的結果，驗證了兩種模型的正確性；Negash 等[16]提出了一種針對環節優化的新型遺傳算法。然而，現階段應用遺傳算法對Bouc-Wen 模型進行參數識別仍然有許多不足，主要在于由局部收斂導致的擬合精度不足和由耗時較長導致的實用性下降。

本文提出了一種非線性自適應遺傳算法(NAGA)，通過考慮染色體庫策略、最優保存策略、改良的適應度函數和選擇算子以及改進的自適應交叉算子和變異算子，對變換軸壓比的鋼筋混凝土剪力墻低周反復加載試驗進行Bouc-Wen 參數識別并和標準遺傳算法進行對比，結果表明該方法顯著提高了識別的有效性。

1 遺傳算法識別模型

Bouc[17]提出了一種多變量、多用途的光滑滯回模型，Wen[18]在1976 年對其進行了改進，稱為Bouc-Wen 模型，并首先用于模擬單自由度系統的非線性滯回特性。如圖1 所示的單自由度系統的運動方程可以用式(1)表示。

圖1 單自由度Bouc-Wen 模型Fig. 1 Single-degree-of-freedom Bouc-Wen model

式中：u(t)為質量塊m相對于地面的位移；c為線性粘滯阻尼系數；F(t) 為 外部激勵荷載；z(t)為滯回位移； α為每個滯回環峰值處切線剛度與初始剛度之比，其范圍為0～1；k為初始剛度。

Bouc-Wen 模型采用式(2)來描述滯回位移。式中，n是控制滯回曲線屈服的尖銳程度參數；A、β、γ三個滯回參數決定了滯回環的基本形狀。

2 標準遺傳算法策略設置

標準遺傳算法(SGA)是由Goldberg[22]在研究控制優化應用中提出的，這為研究人員應用遺傳算法解決各類實際問題奠定了基礎。本文先按照SGA 流程依次進行策略設置，然后，針對其中出現的搜索方向不穩定、收斂效率低的問題進行改進，提出一種識別效果顯著改善的識別方法(NAGA)。主要算法策略介紹如下。

2.1 編碼策略

本文采用二進制編碼，使用基于字符集{0,1}構成的基因串表示Bouc-Wen 模型中的實際參數，下文中由MATLAB 程序所計算出的值均為二進制數轉化而來的實數值，不代表實際控制精度，同時令 β/|γ|≥1[19-20]。Bouc-Wen 模型中，[α0,k,γ,β,δν,δη,A0,δA]這8 個參數都是有一定范圍約束的，并且部分參數不具備明確的物理意義，具體取值也難以標準化，Ma 等[23]用理論推導的方式證明了這8 個參數是有冗余的，因此部分參數的范圍應該是動態的。如表1 所示，本文半定量半定性地選取了較大的解空間，因此，可以在廣泛的8 維參數空間中求解最優解。本文的染色體是這8 個參數的組合，根據每個參數的上下限和精度的要求，本文中染色體長度設置為121，各參數具體分布見圖2。

表1 Bouc-Wen 模型中各參數的上下限和精度Table 1 Limit and precision of Bouc-Wen model parameters

圖2 染色體上各參數對應基因的分布Fig. 2 Distribution of genes corresponding to parameters on chromosome

2.2 初始化種群

由于對于Bouc-Wen 模型的最優解數量及在解空間的分布不明確，為了有效地在實際參數空間中均勻選取樣本，本文利用MATLAB 自帶的隨機函數rand，在染色體基因位上隨機生成0 或1，形成一個初始種群。種群規模N的大小對遺傳算法的搜索質量有著非常大的影響，一般建議值N[24]為20～100?？紤]到整個參數空間比較大，達到2121的規模，為了避免算法發生“早熟”現象，本文的種群規模設定為500。

2.3 適應度函數

對于規模為N 的種群P={I1,I2,···,IN}，染色體Ii∈P的目標函數值E(Ii)可以通過式(9)計算，這個函數通過計算試驗數據和Bouc-Wen 模型擬合數據的均方根誤差(RMSE)值來評估該條染色體的優劣程度。適應度函數可以通過式(10)進行計算。

2.4 選擇算子

2.5 交叉算子和變異算子

本文采用均勻變異策略。作為能夠增加種群多樣性的有效手段，變異概率直接決定算法收斂速度。變異概率如果設置的太小，算法可能陷入早熟，如果過高的話，算法收斂速度過慢，難以生成一個穩定的優解。一般取pm=0.001 ～0.05[25]，本文采用Srinivas[25]推薦的0.005。

3 非線性自適應遺傳算法的提出

在上述SGA 流程的基礎上，本文提出染色體庫策略、最優保存策略，設計了一種雙選擇算子，引入自適應交叉變異策略和有界多點交叉方法來對其進行改進。方法介紹如下。

3.1 染色體庫策略

本文種群規模N=500，遺傳代數為G=200。在本文第4 節中4 片試驗剪力墻在3 次SGA 過程中的平均染色體重復率見圖3，由此可以看出重復出現的個體占30%左右，甚至還有重復出現14 次的個體，如果能夠記錄重復個體的適應度，整個優化過程將會節省大量算力。

圖3 3 次SGA 過程中各片剪力墻的染色體平均重復率Fig. 3 Chromosome average repetition rate of each shear wall in three SGA processes

本文采用牛頓-拉夫遜迭代求解Bouc-Wen 模型，要提高算法運算效率，就需要盡量減少耗時占總時長99%以上的數值求解函數的調用次數，由此提出染色體庫策略。染色體庫是將運算過程中出現過的染色體儲存到一個數組中，同時也將該染色體對應的RMSE(E)值儲存到同一個數組中，這樣形成了一個對于單個問題的染色體庫。在每次求解Bouc-Wen 方程之前，先查詢染色體庫，如果有相同染色體的結果，直接調用；如果沒有，則進行計算，同時將這些還未收錄的染色體連帶其RMSE(E)值收錄進染色體庫。若染色體數量達到染色體庫可收錄上限，則通過比較替換染色體庫中RMSE(E)值最大的，也就是適應度最差的染色體。這是由于隨著遺傳進化的推進，種群逐步收斂到適應度較高的個體上，那么這些適應度高的染色體出現地概率更大，如此染色體庫能夠得到更有效地調用。染色體庫大小可設計為種群規模的5 倍～10 倍，本文中由于N和G值較大，因此取10。

3.2 最優保存策略

Rudolph[26]通過齊次有限馬爾科夫鏈分析證明，由于SGA 的不可化歸性，無論如何初始化，目標函數、交叉變異算子如何選取，其都不會收斂到全局最優。因此，本文采用Negash 等[16]使用的最優保存策略，即在每一代中選取表現最好的染色體直接復制到下一代種群中，以此來保證每代中的最優個體不會在遺傳過程中遭到破壞。同時，Rudolph[26]證明了采用此策略后算法可以收斂到全局最優解。

3.3 適應度函數改進

圖4 SGA 進化過程中最差染色體RMSE(E)的常用對數Fig. 4 Common logarithm of RMSE(E) of worst chromosome during SGA evolution process

為了研究a值對變換后適應度值離散度的影響并進一步確定a的值，本文選擇剪力墻SW1-1在SGA 過程中第200 代染色體的適應度值作為分析對象。因為當a的值不同時，變換后的適應度的量級通常是不同的，所以不可能直接使用標準差來進行判斷。本文提出了離散度的概念來進行判斷，見式(14)。

SGA 進化過程后期的種群，染色體適應度值較為集中，因此 Υ較小。由表2 可以看出，a值從0.5 到1.0 時，種群的 Υ增加，而當a從1.0 到2.5 時，種群 Υ逐漸減小。當取1.0 時，分析對象的離散度是最大的，因此，本文中a取1.0。

表2 不同a 值對應的離散度Table 2 Measure of dispersion at different values of a

3.4 選擇算子改進

SGA 進化后期，由于染色體的適應度都非常接近，每個染色體被選擇的概率也相近，這會使得輪盤賭選擇過程變成隨機選擇過程，如此算法整體收斂性會大打折扣。因此本文提出一種雙選擇算子對此進行改進。在遺傳進化前中期，使用輪盤賭策略在種群中保持較小的選擇壓力，使得各種基因模式都能得到保留；而到進化后期，為加快算法收斂速度，提高選擇壓力，采用基于排序的確定性選擇算子。雙選擇算子的流程圖見圖5。引入遺傳代數g作為控制變量，當g≤0.6G時，采用輪盤賭算子；當g＞0.6G時，采用基于排序的確定性選擇算子，即根據適應度對每一代的染色體進行排序，將排在前25%的染色體復制兩次到下一代，排在中間50%的個體直接復制到下一代，排在末尾25%的個體直接淘汰。這樣后期通過確定性的篩選，種群會快速收斂。

圖5 雙選擇算子Fig. 5 Double-choice operator

3.5 交叉算子和變異算子改進

遺傳算法的全局搜索能力由交叉算子控制，而局部搜索能力由變異算子控制。對于SGA 來說，其交叉概率和變異概率都是恒定的，交叉概率和變異概率的合理設定一般需要通過大量的試驗，而且很難找到一個放在遺傳進化各個階段都合適的常數。同時固定的交叉概率和變異概率使得算法在遺傳進化后期，有一定可能性破壞表現較好的染色體，引起種群的“退化”現象。交叉概率和變異概率如果能夠依照每個染色體自身適應度值進行調整，使得較高適應度個體的基因組合能夠有較大概率保留，而較低適應度個體能夠通過大交叉概率和大變異概率突破現有的種群模式，這種方法使得整個算法能夠跳出局部收斂，將充分改善算法的運行效率和魯棒性。

Srinvivas 等[25]首先提出自適應遺傳算法，其交叉概率和變異概率只考慮了與適應度的關系，且是線性關系。本文在此基礎上引入遺傳代數和種群規模作為控制變量來調整交叉概率和變異概率。在遺傳初期，交叉變異概率維持在基本水平，到進化后期，抑制表現較好的個體交叉變異，加劇表現較差的個體交叉變異，如此在加快算法在當前最優解域收斂的同時，盡可能探索更多的基因模式。交叉概率和變異概率的具體形式見式(15)～式(16)[27]。

考慮到一點交叉有一定概率破壞較好的基因模式，不利于長距離模式的保留和重組。因此本文借鑒建筑模塊假說思想[13]，提出一種有界多點交叉方法，定向對每個參數所在編碼區域進行交叉操作，如此使得算法搜索具有一定方向性，詳細的操作過程見圖6。首先對于每個參數生成一個隨機數，若該隨機數小于交叉概率則該參數對應的基因片段要進行交叉，再根據各隨機數確定各交叉點。

圖6 有界多點交叉Fig. 6 Bounded multipoint crossover

4 RC 剪力墻低周反復加載試驗驗證

4.1 試驗簡介

作為一種重要的豎向承載和抗側力構件，鋼筋混凝土剪力墻被廣泛地應用在工業與民用建筑中。但是由于剪力墻在服役過程中受到彎、剪、扭的共同作用，和混凝土材料的復雜非線性行為等，使得剪力墻的滯回特性難以準確預測。本文中用于驗證的4 片變換軸壓比剪力墻[28]參數見表3。

表3 剪力墻構件參數Table 3 Shear wall component parameters

RC 剪力墻低周反復加載試驗見圖7。試驗采用擬靜力方式進行加載。正式試驗前，先在墻體頂部施加大小為預定荷載40%的豎向壓力，重復加載2 次～3 次，以消除試件內部組織不均勻性，然后將豎向壓力加至預定荷載并在試驗中保持不變。

圖7 RC 剪力墻低周反復加載試驗Fig. 7 Low-cycle cyclic loading test of RC shear wall

側向加載分兩個階段。第一階段采用力控制，先單調逐級加載至開裂，循環一次，以后各級荷載循環一次，直至屈服；第二階段采用位移控制，分別按屈服時頂點位移的倍數逐級加載，每級循環三次，直至構件的承載力下降到峰值承載力的80%左右為止。在整個加載的過程中，雖有加載方式的區別，但是無論是力還是位移控制加載，本文都是取用最后測量得到的不同荷載步上的一組力-位移值，大量力-位移值就構成了本文的識別數據組，在這樣的數據組空間中，本文采用了位移作為控制變量的目標函數式(9)，值得一提的是，該控制變量與試驗加載的控制方式的概念并不相同。SW1-1 的具體側向加載歷程如圖8所示。

圖8 SW1-1 的側向加載歷程Fig. 8 Specific lateral loading history of SW1-1

4.2 識別結果

本文半定量半定性地給8 個參數選擇了較大的解域，SGA 和NAGA 識別得到的多組8 參數未聚集到到所設定范圍的邊界，說明所選解域已包括全局最優解。NAGA 克服了SGA 易局部收斂的問題，可以不斷向全局最優解接近，從而提高識別效果。

4.2.1 參數識別結果

由于遺傳算法有一定的隨機性，本文對每片墻進行10 次識別，引入變異系數cv來衡量其結果的波動情況，具體表達式見式(17)。表4 為四片剪力墻平均RMSE(E)值μ及其變異系數，總體來說，NAGA 參數識別效果要明顯優于SGA，并且變異系數更小，表明算法收斂性更好。RMSE(E)值最小的SGA、NAGA 參數識別結果分別見表5、表6，并將基于該參數組合進行后續分析。

表4 4 片剪力墻參數識別波動情況Table 4 Parameter identification fluctuation of four shear walls

表5 SGA 參數識別結果Table 5 Parameter identification results based on SGA

表6 NAGA 參數識別結果Table 6 Parameter identification results based on NAGA

4.2.2 滯回曲線

由圖9 可知，所提出的非線性自適應遺傳算法識別效果更好，模擬滯回曲線的飽和度和試驗滯回曲線更接近。從4 片剪力墻的滯回曲線可以看到，加載過程中沒有出現明顯的捏攏現象，說明第二節中對于Bouc-Wen 模型的簡化是合理的。

圖9 基于SGA 和NAGA 的模擬滯回曲線比較Fig. 9 Comparison of simulated hysteretic curves based on SGA and NAGA

決定系數是用來評價模型的擬合程度的，一般來說決定系數越大，說明模擬的效果越好，其值不大于1，具體可由式(18)計算。由表7 可知改進之后算法識別的效果更好。

表7 SGA 和NAGA 的決定系數Table 7 Coefficient of determination of SGA and NAGA

4.2.3 力與位移響應過程

從圖10 可以看到，對于SGA 由于識別出的初始剛度過大使得初期計算出的側向力較大，而剛度退化和強度退化過大使得在加載中后期側向力計算值小于試驗值；而對于NAGA 來說，雖然由于識別出的初始剛度偏小使得初期模擬曲線的側向力小于試驗力，但當數據量足夠多時識別和擬合效果優異。

圖10 基于SGA 和NAGA 的力與位移響應過程比較Fig. 10 Comparison of force and displacement response process based on SGA and NAGA

4.2.4 算法效率分析

由圖11 可知，引入最優保存策略后，種群的最優個體穩定地朝著使RMSE(E)更小的方向進化，最優個體不會出現波動甚至“退化”現象。NAGA 在100 代左右就收斂到當前最優解，收斂速度較為快速，也表明G的取值是合理的。

圖11 SGA 與NAGA 進化過程中最優染色體的RMSE(E)變化對比Fig. 11 Comparison between change of RMSE(E) of best chromosome during evolution of NAGA and SGA

各片剪力墻最優的RMSE(E)及其識別的時間見表8。由此可見，在識別效果方面，NAGA 最優染色體的RMSE(E)比SGA 的值小了40%；在識別速度方面，由于染色體庫策略的引入使得數值求解函數的調用次數減少，同時算法收斂性的增強也使得數值求解過程中可以更快地滿足各荷載步滯回位移的精度要求，如此NAGA 的速度相較SGA 有了一定的優勢。

表8 SGA 和NAGA 實用性比較Table 8 Practicality comparison between SGA with NAGA

5 結論

本文在標準遺傳算法的基礎上，通過引入染色體庫策略、最優保存策略，改良適應度函數和選擇算子，并在現有自適應交叉算子和變異算子的基礎上進一步改良，提出了一種非線性自適應遺傳算法。通過4 組變換軸壓比的RC 剪力墻低周反復加載試驗數據驗證了該方法的有效性，與已有的SGA 方法相比，本方法在Bouc-Wen 模型參數識別上具有明顯的優越性。本文主要結論如下：

(1) 本文所提出的非線性自適應遺傳算法(NAGA)對于8 參數Bouc-Wen 模型參數的識別效果良好，通過所識別的參數能夠較好地模擬剪力墻這種結構構件在低周反復加載中的滯回行為；

(2) 與試驗數據相比，在本文識別表征RC 剪力墻滯回行為的Bouc-Wen 模型參數的擬合效果顯著優于現有的標準遺傳算法(SGA)；

(3) 研究中染色體庫策略的引入和算法收斂性的增強使得參數識別速度有了一定提升，算法實用性得到有效增強。

由于遺傳算法的隨機性，要想通過改變其超參數，用分析滯回曲線、RMSE(E)或時間等定量指標來分析是困難的，目前只能做到定性分析來得到一組合適的遺傳算法參數。本文所述方法的提出可用于混合模擬試驗的模型參數識別中，通過試驗子結構的滯回性能參數實時識別，并及時反饋給計算子結構進行整體性能模擬分析，也同時適用于大型復雜結構在服役期間根據所識別得到的關鍵結構構件的滯回特性預測整體結構的抗災能力。