基于Kriging模型及代理優化算法的老舊橋梁有限元模型修正

楊 鵬

(中國鐵建昆侖投資集團有限公司 四川成都 610000)

1 引言

當前橋梁有限元模型的建立多以橋梁設計圖紙為依據，在建模時引入了多種理論假設，并且存在邊界條件近似以及節點連接剛度模擬不恰當等問題，尤其當橋梁服役若干年之后，構件開裂等病害使得橋梁剛度下降，結構整體的力學性能發生變化。為了更好地模擬舊橋的力學性能，為日后的監測、保養以及加固提供參考依據，需要依靠模型修正技術來校正按照圖紙建立的理論有限元模型。

對于模型修正問題而言，可以將數值分析軟件看成一個輸入輸出系統，即作為輸出量的目標函數和約束函數可看成是待優化參數的函數。通過有限元軟件計算一次這個函數需要很長的時間，不利于很多啟發式算法的使用。如果可以通過建立一個模型繞過有限元計算，而直接擬合這兩個函數，則可以大大提高計算效率。Kriging模型作為一種無偏估計模型，具有出色的非線性擬合能力，在特定的隨機過程中，能夠給出最優線性無偏估計。在現有研究中，Kriging模型大量應用在地質、水文、氣象、航空航天以及飛行器設計領域[1]。

2 Kriging代理模型優化算法

2.1 Kriging代理模型

首先在優化空間通過抽樣方法，得到n個樣本點:

對這n個樣本點進行有限元分析，得出n個分析結果:

對這些結果進行插值，即:

只需要給出加權系數ω＝[ω1ω2…ωn]，便可得到任意參數對應的真實響應值。為此，Kriging模型先將預測值和輸入值之間的關系假設為:

式中，μ為常數項；z(x)是按N(0，σ2)分布的隨機變量。前者表示預測值的整體趨勢值，后者反映了預測值的局部波動偏差。

在優化空間的不同位置處，隨機變量存在一定的相關性，對于隨機過程z(x)，任意兩點隨著距離的接近，將趨于相同的函數值。Kriging模型將這種相關性表示為:

式中，θk為不同維度上相關性的衰減速率，決定了模型的擬合性能。

基于以上假設來尋找最優的加權系數，使得均方誤差MSE最小[2]:

可以得到未知點x處Kriging模型的預測值為:

采用極大似然估計的方法尋找θk，將問題轉化為尋找θk使式(8)的值最大[3]:

采用全局搜索算法尋找最優的θk，θk也隨著樣本庫的更新而更新來獲得更好的性能。

2.2 優化加點策略

建立好初始的代理模型之后，下一步就是選擇添加新的樣本點，以較少的樣本點最大程度地提高擬合精度，同時不斷添加樣本點直到優化收斂[4]。

目前國際上已經發展出多種基于代理優化方法的加點準則，有代理模型最優預測準則(MSP準則)、改善期望準測(EI準則)、均方誤差準則(MSE準則)等。

MSP準則尋找的樣本點是可能的最優解，可以快速提高最優解附近的精度；MSE準則尋找的樣本點是誤差期望值最大的點，可以快速提高全局精度。結合兩種準則的特點，通過計算均方誤差最大的點更新預測誤差最大的區域，同時計算更新代理模型預測的最優解。在迭代的前期，計算預測均方誤差最大的點可以快速改善模型的精度；在迭代次數足夠多時，可以認為模型擬合的精度已經達到要求，這時不再計算均方誤差最大的點而只計算代理模型預測的最優解進行加點，加快迭代收斂的速度[5]1716-1718。

2.3 尋找篩選更新點

首先選取好優化目標，針對每一個目標，分別建立一個代理模型。根據兩個加點準則，建立兩個多目標優化問題[6]。

本文采用NSGA-Ⅱ方法求解多目標優化問題，應用NSGA-Ⅱ算法求解式(9)和式(10)分別得到候選的更新樣本點集H1和H2，為了使Pareto解更加均勻，根據極大極小距離準則來選擇正式的更新點:

選擇的正式樣本與已有樣本的最小距離，在待更新樣本集中是最大的[7]，即:

式中，x*為選擇的樣本點。

更新樣本集S中，再尋找H1中極大極小距離對應的樣本點，這兩個樣本點就是篩選出的正式更新點。如果沒有求解式(10)，則從求解式(9)得到的樣本集中貫序選擇兩個點作為正式更新點[5]1720。

2.4 代理優化方法

多目標優化方法具體步驟如下[8]:

(1)首先確定好優化計算的迭代次數Tm，加點時計算預測方差最大點的次數Tmse和初始樣本數n，定義優化變量、約束條件、優化目標。

(2)運用拉丁超立方抽樣進行實驗設計，通過有限元軟件精確求解得到初始樣本集S，從初始樣本集S中篩選出Pareto最優解集A。

(3)基于樣本集S構建Kriging代理模型。

(4)如果迭代次數小于等于Tmse，則使用NSGAⅡ算法分別求解本文定義的式(9)和式(10)，分別得到候選樣本集H1和H2；如果加點個數大于Tmse，則只求解式(9)，得到候選樣本集H2。

(5)如果有H1，則先篩選H1中的點為樣本點，計算H1中每一個點與已有樣本點之間的距離的最小值，每個候選樣本點都有一個對應的最小值，從這些最小值中選取最大的那個對應的候選點加入樣本集S；再遵循這一準則從H2中選取一個點加入樣本集S和最優解集A。如果沒有計算H1，則從H2中選取兩次，得到兩個點加入樣本集S和最優解集A。最后維護最優解集A，只保留非支配個體。

(6)判斷終止條件，若迭代次數大于Tm，則輸出A作為最優解集合；否則，則轉到步驟(3)繼續迭代。

3 有限元模型

3.1 工程背景

背景工程上構為5×16 m預制鋼筋混凝土簡支T梁，下部結構為柱式墩、樁柱式臺、樁基礎，如圖1所示。設計荷載:汽車—20級，掛車—100級(1972標準)。

圖1 橋型示意

3.2 建立模型

利用有限元軟件對橋梁第一跨的上部結構建模，主梁、橋面鋪裝以及路緣石等使用實體結構單元模擬，鋼筋采用三維桿系單元模擬[9]，如圖2所示。

圖2 有限元模型(第一跨)

3.3 施加荷載

根據橋梁檢測報告，橋梁靜載試驗加載車采用290 kN雙后軸車進行加載，如圖3所示。

圖3 軸重及軸距示意

荷載布置如圖4所示。

圖4 荷載布置示意

4 模型修正

4.1 參數及目標選取

根據橋梁檢測報告，已知橋梁結構中存在損傷，實際結構與圖紙的理論模型的區別在于理論模型無法模擬出實際結構由于長期使用而產生的剛度變化，因此選擇4片梁的跨中以及1/4跨撓度(見圖5)和結構整體的一階自振頻率共9個量為修正目標[10]。根據檢測報告的損傷情況，選取3、4、5號梁(梁編號如圖4所示)跨中混凝土彈性模量e1以及密度ρ1，端部混凝土彈性模量e2以及密度ρ2，1、2、6、7號梁跨中混凝土彈性模量e3以及密度ρ3，端部混凝土彈性模量e4以及密度ρ4，橫隔板混凝土彈性模量e5，共計9個變量的系數作為修正參數，其變化區間設[0.2，2]，據此可以建立以下最優化問題:

圖5 跨中及端部混凝土位置示意

式中，yn(X)為數值模型計算的撓度以及頻率與檢測結果的殘差；X為修正系數。

考慮到混凝土中箍筋、架立鋼筋等的布置以及混凝土后期強度的增長，并結合現場觀察橫隔板等結構的損壞情況，各參數初始取值定為:

對有限元模型進行靜力分析以及特征值分析，得到中梁、邊梁跨中、1/4跨撓度以及結構一階自振頻率[11]，見表 1。

表1 有限元模型修正前結果與檢測結果對比

可以看出，修正前的有限元模型與橋梁的實際靜動力特性存在較大差異。

4.2 修正結果

在進行了76次迭代后，最終得到一個包括15個最優解的解集。為評估修正后模型的準確性，選取偏載工況下的位移檢測值與有限元模型計算值進行對比[12]，選取其中一組結果如表2、表3所示。

表2 修正后結果與檢測結果對比

表3 修正后偏載結果與檢測結果對比

結合有限元分析結果以及橋梁檢測報告，第一組修正參數相較于其他幾組更加符合橋梁實際狀況。

5 結論

(1)修正后的模型計算的最大相對誤差為3%，修正后的模型能夠較為精確地反映該橋的現有狀態，為后期該橋的加固計算提供了一個基準模型。

(2)改良后的Kriging代理優化方法可以通過較少的有限元計算次數得到較好的修正效果，且能夠忽略檢測數據中的奇異值，具有良好的魯棒性。

(3)相較于將多目標問題解決方法轉化為單目標問題進行求解，直接進行多目標優化能夠獲得多個可行的修正結果，結果之間具有均勻性、多樣性。能結合工程實際和進一步的檢測結果選擇較合理的參數，從而保證修正之后的結果不失真。